Portable, Extensible Toolkit for Scientific Computation/PETScの応用事例

はじめに

本章では、PETScを実際の科学技術計算問題に適用した事例を紹介します。具体的な問題を通じて、PETScの実践的な使い方を学びます。取り上げる事例は以下の通りです。

1. 熱伝導問題
2. 流体力学問題
3. 構造解析問題
4. 最適化問題

熱伝導問題

問題設定

熱伝導問題は、物体内の温度分布を求める問題です。以下の偏微分方程式（熱方程式）を解きます。 $\frac{\partial u}{\partial t}=\alpha {\mathrm{\nabla }}^{2}u$ ここで、$u$ は温度、$\alpha$ は熱拡散率です。

PETScによる解法

PETScの時間ステッピングソルバー（TS）を使用して、熱方程式を解きます。

例: 熱伝導問題のコード

#include "petsc.h"

PetscErrorCode FormRHSFunction(TS ts, PetscReal t, Vec U, Vec F, void *ctx) {
    // 右辺関数の実装
    Mat A;
    TSGetMatrix(ts, &A);
    MatMult(A, U, F);
    return 0;
}

int main(int argc, char **argv) {
    PetscInitialize(&argc, &argv, NULL, NULL);

    TS ts;
    Vec U;
    Mat A;

    // 初期化と設定
    TSCreate(PETSC_COMM_WORLD, &ts);
    VecCreate(PETSC_COMM_WORLD, &U);
    MatCreate(PETSC_COMM_WORLD, &A);

    // 行列とベクトルの設定
    // ... (省略)

    // 時間ステッピングの設定
    TSSetProblemType(ts, TS_LINEAR);
    TSSetRHSFunction(ts, NULL, FormRHSFunction, NULL);
    TSSetFromOptions(ts);

    // ソルバーの実行
    TSSolve(ts, U);

    // 終了処理
    TSDestroy(&ts);
    VecDestroy(&U);
    MatDestroy(&A);
    PetscFinalize();
    return 0;
}

流体力学問題

問題設定

流体力学問題は、流体の速度場や圧力場を求める問題です。ナビエ-ストークス方程式を解きます。 $\frac{\partial 𝐮}{\partial t}+(𝐮\cdot \mathrm{\nabla })𝐮=-\mathrm{\nabla }p+\nu {\mathrm{\nabla }}^2𝐮$ ここで、$𝐮$ は速度場、$p$ は圧力、$\nu$ は動粘度です。

PETScによる解法

PETScの非線形ソルバー（SNES）を使用して、ナビエ-ストークス方程式を解きます。

例: 流体力学問題のコード

#include "petsc.h"

PetscErrorCode FormFunction(SNES snes, Vec U, Vec F, void *ctx) {
    // 非線形関数の実装
    Mat J;
    SNESGetJacobian(snes, &J, NULL, NULL, NULL);
    MatMult(J, U, F);
    return 0;
}

int main(int argc, char **argv) {
    PetscInitialize(&argc, &argv, NULL, NULL);

    SNES snes;
    Vec U, F;
    Mat J;

    // 初期化と設定
    SNESCreate(PETSC_COMM_WORLD, &snes);
    VecCreate(PETSC_COMM_WORLD, &U);
    VecCreate(PETSC_COMM_WORLD, &F);
    MatCreate(PETSC_COMM_WORLD, &J);

    // 行列とベクトルの設定
    // ... (省略)

    // 非線形ソルバーの設定
    SNESSetFunction(snes, F, FormFunction, NULL);
    SNESSetFromOptions(snes);

    // ソルバーの実行
    SNESSolve(snes, NULL, U);

    // 終了処理
    SNESDestroy(&snes);
    VecDestroy(&U);
    VecDestroy(&F);
    MatDestroy(&J);
    PetscFinalize();
    return 0;
}

構造解析問題

問題設定

構造解析問題は、物体の変形や応力を求める問題です。以下の弾性方程式を解きます。 $\mathrm{\nabla }\cdot \sigma =𝐟$ ここで、$\sigma$ は応力テンソル、$𝐟$ は外力です。

PETScによる解法

PETScの線形ソルバー（KSP）を使用して、弾性方程式を解きます。

例: 構造解析問題のコード

#include "petsc.h"

int main(int argc, char **argv) {
    PetscInitialize(&argc, &argv, NULL, NULL);

    KSP ksp;
    Vec U, F;
    Mat A;

    // 初期化と設定
    KSPCreate(PETSC_COMM_WORLD, &ksp);
    VecCreate(PETSC_COMM_WORLD, &U);
    VecCreate(PETSC_COMM_WORLD, &F);
    MatCreate(PETSC_COMM_WORLD, &A);

    // 行列とベクトルの設定
    // ... (省略)

    // 線形ソルバーの設定
    KSPSetOperators(ksp, A, A);
    KSPSetFromOptions(ksp);

    // ソルバーの実行
    KSPSolve(ksp, F, U);

    // 終了処理
    KSPDestroy(&ksp);
    VecDestroy(&U);
    VecDestroy(&F);
    MatDestroy(&A);
    PetscFinalize();
    return 0;
}

最適化問題

問題設定

最適化問題は、与えられた制約条件下で目的関数を最小化または最大化する問題です。以下の最適化問題を解きます。 $\underset{𝐱}{\min }f(𝐱)$ ここで、$f(𝐱)$ は目的関数です。

PETScによる解法

PETScの最適化ツール（TAO）を使用して、最適化問題を解きます。

例: 最適化問題のコード

#include "petsc.h"
#include "petsctao.h"

PetscErrorCode FormObjective(Tao tao, Vec X, PetscReal *f, void *ctx) {
    // 目的関数の実装
    Vec g;
    VecDuplicate(X, &g);
    VecNorm(X, NORM_2, f);
    VecDestroy(&g);
    return 0;
}

int main(int argc, char **argv) {
    PetscInitialize(&argc, &argv, NULL, NULL);

    Tao tao;
    Vec X;

    // 初期化と設定
    TaoCreate(PETSC_COMM_WORLD, &tao);
    VecCreate(PETSC_COMM_WORLD, &X);

    // ベクトルの設定
    // ... (省略)

    // 最適化ソルバーの設定
    TaoSetObjective(tao, FormObjective, NULL);
    TaoSetFromOptions(tao);

    // ソルバーの実行
    TaoSolve(tao);

    // 終了処理
    TaoDestroy(&tao);
    VecDestroy(&X);
    PetscFinalize();
    return 0;
}

練習問題

1. 熱伝導問題: 2次元領域での熱伝導問題を解くプログラムを作成し、温度分布を可視化してください。
2. 流体力学問題: 簡単な流れ（例: ポアズイユ流れ）をシミュレートするプログラムを作成してください。
3. 構造解析問題: 梁の変形を解析するプログラムを作成し、変形状態を可視化してください。
4. 最適化問題: 二次関数の最小化問題を解くプログラムを作成し、最適解を表示してください。

まとめ

本章では、PETScを実際の科学技術計算問題に適用した事例を紹介しました。熱伝導問題、流体力学問題、構造解析問題、最適化問題を通じて、PETScの実践的な使い方を学びました。これらの事例を参考に、自身の問題にPETScを適用してみてください。