■LU Decomposition(三角行列分解法)の適用■
3次元に入る前に2次元のNSEQ8DD.FORにLU Decompositionを組み込むことから話をはじめます。 LU Decompositionで説明した手法を流体ソフトに適用する場合、時間ステップ(Δt)を固定しなければなりません。 なぜなら、[A]=(1/Δt)[GEOM] + [K]ですからΔtを途中で変える度にCALL LLDECOMPVTが必要になりLU Decomposition の効果が発揮されません。。

LU Decompositionを用いNSEQ8DD.FORを改造すると、計算速度は数段速くなります。ただし、 CALL LLDECOMPVTでちょっと時間がかかります。今後は、LLDECOMPVTの効率化を行う必要がありますね。 作成したプログラム名は、以下になっています。

おさらいですが、ここのサイトで紹介している流体のそふとでは、時間方向の積分に Crank-Nicolson methodを使っていました。この方法は精度の向上には役立ちますが プログラムが複雑になるという欠点がありました。上の表にあるように、 陰解法(Implicit)を時間軸積分に用いています。3次元のプログラムも陰解法を用いてあります。
また、3次元に向けて入力データのフォーマットも変更しました。上で示した入力データ作成プログラム では新しフォーマットでデータを作成します。 おさらいですがNSEQ8DD.FOR(上表のNSEQ8DD4.FORのこと)では、時間方向の積分にCrank-Nicolson methodを使っています。 この方法は精度の向上には役立ちますがプログラムが複雑になるという欠点がありました。 そのために、t時で[K]{u}の計算を行わなければなりません。 つまり、[K]またはそれに順ずるメモリーが余分に必用になります。 また、プログラムも複雑になっているので、今回3次元プログラムを作成する前準備として、 もっとシンプルな陰解法(Implicit)のNSEQ8DD.FORを作成しました。 プログラム名は上に示してあります。 入力データの作成には、SET4NS12.FOR、SET9NS12.FORを使っています。 計算結果をエクセルで表示できるプログラムも用意してあります。

■3次元流体プログラム■
3次元の流体プログラムは、2次元陰解法プログラムと3次元弾性解析プログラムをベースに作成しました。 基本的には、カットアンドペーストと少々の変更のみで3次元流体プログラムが出来上がりました。 プログラムには8節点ヘキサ要素用と27節点ヘキサ要素用の2つがあります。両者とも連立方程式の解法には、 三角行列分解が使われています。プログラム名は以下です。