メニュー

技術情報 Techinicalinfo

  1. ホーム
  2. 技術情報
  3. 【技術情報】有限要素法入門
  4. 7.5 有限要素法のプログラム(その7)

【技術情報】有限要素法入門

7.5 有限要素法のプログラム(その7)


ここで有限要素法のサンプルプログラムを作成して実行することにします。
まず次のメインプログラムをホルダー Source の中にファイル名 femMain.cpp として作成してください。


#include "../include/fem.h"

int main()
{
    //有限要素法のメインプログラム

    //解析クラスの作成
    Analysis* analysis = new Analysis();

    int fieldType = 1;
    double eps = 1.0e-6;

    //入力ファイル
    ifstream fp_in; fp_in.open("input.csv", ios::in);
    if(!fp_in)
    {
        cout << "入力ファイル input.csv がありません!" << endl;
        return 1;
    }
    //出力ファイル
    ofstream fp_out; fp_out.open("output.lpr", ios::out);
    //チェック用ファイル
    ofstream fp_chk; fp_chk.open("output.chk", ios::out);

    //解析データの読み込み
    analysis->read(&fp_in); 

    bool lnflg = analysis->isLinearanalysis(); 

    //節点スカラーFEMの作成
    NodalScalarfem* nodalScalarfem = new NodalScalarfem(analysis,lnflg,fieldType,eps,&fp_chk);

    nodalScalarfem->setAddress();                //アドレスの作成

    nodalScalarfem->setIndexmatrix();            //インデックス行列の作成

    int nstep = nodalScalarfem->getNstep();      //時間テーブルの作成

    for (int istep=1; istep<=nstep; istep++)
    {
        nodalScalarfem->setTime(istep);          //時間の設定

        nodalScalarfem->exec();                  //実行

        nodalScalarfem->output(&fp_out);         //結果出力

        nodalScalarfem->renew();                 //次のステップの設定
    }

    return 0;
}

次にこのプログラムを実行する環境を作ります。MinGW の場合について述べます。
コマンドプロンプトを立ち上げたときのカレントディレクトリーに次のバッチファイルを作成してください。
ここでは、D ドライブのホルダー Gcc の下の MinGW64 にシステムが入っている場合を例にしていますので、実際に WinGW をインストールしたホルダーを指定してください。
また、ここで作ったプログラムのホルダー FEM が D ドライブの直下に作られていることを想定してますのでそうでない場合はその部分をかきかえてください。


@ECHO OFF
setlocal
rem ; 次の行をMinGWをコピーしたフォルダに書き換えてください。
SET MINGW_DIR=D:\Gcc\\MinGW64

SET PATH=%MINGW_DIR%\bin;%MINGW_DIR%\mingw32\lib;%MINGW_DIR%\lib\gcc\mingw32\4.4.0;%PATH%;
SET CLANG=C-SJIS
SET LANG=ja_JP
D:
cd \FEM\compile
%COMSPEC%
endlocal

このバッチファイルを実行するとカレントディレクトリーは FEM の下のディレクトリー Compile に移ります。
ここに次のバッチファイルを fem.bat という名で作成してください。


G++ -Wall -fopenmp -c ..\Source\FEM.cpp -o ..\Object\FEM.obj
G++ -Wall -fopenmp -c ..\Source\solver.cpp -o ..\Object\solver.obj
G++ -Wall -fopenmp -c ..\Source\femMain.cpp -o ..\Object\femMain.obj
G++ -Wall -fopenmp -o ..\Go\FEM.exe ..\Object\*.obj

コマンドプロンプトから fem と打ち込むと実行用モジュール frm.exe が FEM の下のホルダー Go に作成されます。
これで準備が出来ました。
次にデータですが、実行用モジュールと同じディレクトリー Go にファイル名 input.csv として次の用に作成してください。


2,
0,1,1,1,16,9,
1,
誘電体
1.0,1.0,
0,0,
1,
0.0,
1,2,
300.0,0.0,0.0,
0,
1,1,0,0.0,0.0,0.0,
2,1,0,0.1,0.0,0.0,
3,1,0,0.2,0.0,0.0,
4,1,0,0.3,0.0,0.0,
5,1,0,0.0,0.1,0.0,
6,0,0,0.1,0.1,0.0,
7,0,0,0.2,0.1,0.0,
8,0,0,0.3,0.1,0.0,
9,1,0,0.0,0.2,0.0,
10,0,0,0.1,0.2,0.0,
11,0,0,0.2,0.2,0.0,
12,0,0,0.3,0.2,0.0,
13,1,0,0.0,0.3,0.0,
14,0,0,0.1,0.3,0.0,
15,0,0,0.2,0.3,0.0,
16,0,0,0.3,0.3,0.0,
1,24,1,1,1,2,6,5,
2,24,1,1,2,3,7,6,
3,24,1,1,3,4,8,7,
4,24,1,1,5,6,10,9,
5,24,1,1,6,7,11,10,
6,24,1,1,7,8,12,11,
7,24,1,1,9,10,14,13,
8,24,1,1,10,11,15,14,
9,24,1,1,11,12,16,15,

このデータは2次元の一辺が0.3mの誘電体を縦横3分割にした有限要素法用のデータで、節点数が16,要素数が9のモデルです。誘電体の誘電率は1で各要素に電荷密度300を与えています。
実行はディレクトリー Go に移り fem と打ち込むと開始します。

実行が終わるとこのディレクトリに二つのファイル output.lpr と output.chk が作成されます。
output.lpr に計算結果として節点の電位と要素の電場ベクトルが出力されています。
output.chk はチェックファイルです。

求まったポテンシャルの値を以下に示します。
計算結果を以下に示します。
\begin{equation}
\begin{split}
&\phi_{13}=0 \hspace{5mm} \phi_{14}=4.81585 \hspace{6mm} \phi_{15}=7.35435 \hspace{5mm} \phi_{16}=8.13784 \\
&\phi_{09}=0 \hspace{5mm} \phi_{10}=4.41488 \hspace{6mm} \phi_{11}=6.67133 \hspace{5mm} \phi_{12}=7.35435 \\
&\phi_{05}=0 \hspace{5mm} \phi_{06}=3.06264 \hspace{6mm} \phi_{07}=4.41488 \hspace{5mm} \phi_{08}=4.81585 \\
&\phi_{01}=0 \hspace{5mm} \phi_{02}=0 \hspace{21mm} \phi_{03}=0 \hspace{21mm} \phi_{04}=0
\end{split} \tag*{$(7.5-1)$}
\end{equation}
この結果は1章で計算した(1.4-48)式と一致しています。