有限要素プログラミング

はじめに

NDSolveは,有限要素法を使って偏微分方程式を解くための,高レベルでワンステップのインターフェースを提供する.しかし,解くプロセスのステップを詳しく制御したいという場合もある.NDSolve`FEM`パッケージは,解くプロセスの各ステップのおいて広範な制御を与える,低レベルインターフェースを提供する.

有限要素関数を使うためには,パッケージをロードしなければならない.

NDSolve`FEM`パッケージの低レベル関数は,さまざまな用途に使用できる.

NDSolve内の有限要素データ

NDSolveフレームワークには,NDSolve`ProcessEquationsNDSolve`IterateNDSolve`ProcessSolutionsという低レベル関数と,すべてのタイプの微分方程式に対して解法プロセスへの何らかのアクセスを提供するデータオブジェクトNDSolve`StateDataが含まれる.これらの関数の一般的な用法は,「コンポーネントとデータ構造」に説明されている.

残りのNDSolveと一貫して,有限要素法は,そのデータにNDSolve`StateDataを通してアクセスすることができる.次のセクションでは,NDSolve`StateDataオブジェクトから有限要素データにアクセスする方法について見てゆく.有限要素データはNDSolve`FEM`FiniteElementDataオブジェクトに保存され,"FiniteElementData"の特性を通してアクセスされる.

NDSolve`StateDataオブジェクトを設定する.

この問題のように,問題に時間依存がない場合には,NDSolve`StateDataオブジェクトの表示形式は,このことを"SteadyState"を表示することによって示す.非定常方程式については,積分が行われる時間区分が示される.

有限要素データは,状態オブジェクトから抽出される.

有限要素データを得る.

解法プロセス中に生成されたデータオブジェクトはすべて,例えばFiniteElementDataのように,そのコンストラクタ関数として文書化されている.このチュートリアルでは,これらの関数の使用方法について説明する.

FiniteElementDataオブジェクトの表示形式は,離散化された偏微分方程式の自由度を示す.自由度は通常,メッシュ内のノード数に従属変数の数を掛けたものである.自由度は,組み立てられる行列の行と列の数を示す.

定常状態問題については,NDSolve`Iterateを呼び出すことで,LinearSolveを使って方程式系の解を求める.

系の解を計算する.

解はその後,解データオブジェクト内の状態オブジェクトに保存される.

解の値を抽出する.

解は,×1の行列として表される. は自由度である.

InterpolatingFunctionとして解を抽出し,結果をプロットする.

有限要素オプションをNDSolveに渡す

次のセクションでは,NDSolveの内部の有限要素関数について述べる.これらの関数にはオプションを含むものもある.内部の有限要素関数を使用せずに,これらのオプションをNDSolveに渡すことは常に可能である.これらは,NDSolveMethodオプションで与えられる.

NDSolve[,Method{"PDEDiscretization"{"FiniteElement", femopts}}]定常状態問題について,オプション femopts とともに有限要素法を使う
NDSolve[,Method{"PDEDiscretization"{"MethodOfLines", "SpatialDiscretization"->{"FiniteElement", femopts}}}]非定常問題の空間離散化について,オプション femopts とともに有限要素法を使う

有限要素法のオプションを指定する.

有限要素オプションを指定して,NDSolveを呼び出す.

同じことが非定常の偏微分方程式でも可能である.

NDSolveで有限要素オプションを指定して,非定常偏微分方程式を解く.

ToBoundaryMeshToElementMeshのオプションは,上で示すように,"MeshOptions"を指定することによって,NDSolveに渡すことができる.さらに,InitializePDECoefficientsPDESolveのオプションおよびそのサブオプションは,それぞれ"InitializePDECoefficientsOptions""PDESolveOptions"を通してNDSolveに渡すことができる.また,積分次数と補間次数はオプションを使って指定できる.内部的な有限要素法関数のオプションの用途と目的については,それぞれの参照ページおよびオプションについてのチュートリアル,有限要素のためのNDSolveオプションに記載されている.

ワークフローの概要

次のセクションでは,NDSolveが有限要素モデルを解く方法についてステップごとに示す.NDSolveは基本的にこれらのセクションで見られるように,同じ関数を使う.さらに,NDSolve`FEM`FiniteElementDataオブジェクトのコンテンツが下で紹介する関数を通して作成される.

解は,次の3つの段階で求められる.

偏微分方程式問題の設定

偏微分方程式を有限要素法のような数値計算法によって解きやすい状態にするには,3つの成分が必要である.

定常偏微分方程式

初期化段階

初期化段階では,NDSolve`ProcessEquationsの実行中に作成されるすべてのデータオブジェクトが作成される.このセクションは,そのプロセスについて示す.このセクションを飛ばして,離散化段階に進み,NDSolve`ProcessEquationsが作成する初期化されたデータ構造を利用することも可能である.これで,例えば新しい数値離散化法に,NDSolve`ProcessEquationsを方程式プリプロセッサとして使用することができる.

変数データと解データ

NDSolveが解を求めるには,従属変数(u)と独立変数({x,y}),および領域(Ω)を指定しなければならない.このことは,有限要素関数でも同じである.変数データと解データは,別々の成分に対応する特定のデータ構造を作成するリストに保存される.これらのデータリストは,NDSolve`VariableDataNDSolve`SolutionDataを使って生成することもできる.データリストの要素は,NDSolve`SetSolutionDataComponentを使って挿入または変更することができる.

従属と独立の変数名を持つ変数データを作成する.

変数データに加わるのが解のデータ構造である.解のデータ構造は,基本的に,変数データが数値的に生まれ変わったものである.例えば,解データの"Space"成分の数値バージョンは,ElementMeshも含むNumericalRegionである.「要素メッシュの作成」に領域Ωからのメッシュの作成についての詳細が記載されている.

NumericalRegionを指定する.
"Space"成分の集合を持つ解データを作成する.
InitializePDECoefficients

1つの従属変数 で次の偏微分方程式を考える.

偏微分方程式はで定義される.係数 はスカラーである. は,長さ のベクトルであり,× 行列である.InitializePDECoefficientsへの入力係数は,上の偏微分方程式とモデル化される偏微分方程式を対応させることによって求めることができる.

モデル偏微分方程式の係数.

明示的に指定されていない係数はすべてゼロであると想定される.

偏微分方程式の係数を初期化する.

初期化された係数は,PDECoefficientDataのデータ構造に保存される.表示形式は,初期化された係数の系の大きさ,つまり従属変数の数(1),および演算子の大きさ,つまり空間次元(2)を示す.

初期化された係数の系の大きさと空間次元を抽出する.
生の係数を抽出する.

初期化注に,係数は定常係数,非定常係数等,さまざまな係数のカテゴリに分類される.このことについては詳しく後で説明する.

InitializeBoundaryConditions

偏微分方程式の初期化と同じように,境界条件も初期化されなければならない.

境界条件を初期化する.

境界条件の係数は,偏微分方程式の係数と同じようなカテゴリに分類される.

InitializePDEMethodData

後に続くステップで必要となる,有限要素法データを初期化するには,InitializePDEMethodDataを使う.現在,このフレームワークで使用できる離散化法は,有限要素法だけである.したがって,デフォルトで,InitializePDEMethodDataFEMMethodDataオブジェクトを生成する.初期化中にメソッドに特定の設定が行われる.例えば,メソッドの補間と積分の次数が他の諸々と一緒に設定され,ElementMeshNumericalRegionから生成される.

変数データと解データで有限要素データを初期化する.

初期化された偏微分方程式のメソッドデータは,FEMMethodDataオブジェクトに保存される.FEMMethodDataの表示形式は,自由度,各従属変数の補間次数,および積分次数を示す.

自由度,補間次数,積分次数のクエリを行う.

「要素メッシュの生成」のチュートリアルに説明されているように,メッシュ生成オブションはすべて指定することができる.

有限要素データを変数データと解データで,メッシュ生成のオプションを使って初期化する.
自由度,補間次数,積分次数のクエリを行う.

メソッドの初期化中にElementMeshオブジェクトが生成され,NumericalRegionオブジェクトの中に保存される.

NumericalRegionからElementMeshオブジェクトを抽出する.

方程式系の自由度は,複数の従属変数が与えられている場合には,要素メッシュのノード,従属変数の数,従属変数の補間次数を組み合せたものから求められる.従属変数が1つだけの場合には,自由度は,ElementMeshオブジェクトのノードに対応する.

1つの従属変数については,自由度は,メッシュ内の座標数に対応する.

"IntegrationOrder"は,有限要素の演算子を積分するのに使われる確度の次数である.

積分次数を増やす.

一次要素メッシュについては,デフォルトの積分次数は,二次の正確な積分点を使用するが,二次要素メッシュについては,デフォルトでは四次の正確な積分点を使用する.サポートされる最大の積分次数は5である.

積分次数を減らすと,解が速く解かれるが,この場合には,あまり正確ではない解になる可能性がある.積分次数を増やすと,有限要素モデルを解くのにより長い時間がかかる.

NDSolve`ProcessEquationsを通して作成された状態オブジェクトが同じFEMMethodDataオブジェクトを生成することを調べる.
初期化された有限要素データをNDSolve`StateDataから抽出する

有限要素法のデータを手作業で初期化する代りに,偏微分方程式,境界条件,NDSolve`ProcessEquationsをこのために利用することができる.

NDSolve`StateDataから有限要素データを抽出する.
作成されたオブジェクトが同じであることを確かめる.

離散化段階

離散化段階では,連続偏微分方程式と境界条件は,離散変数で近似される.このプロセスは,偏微分方程式を解くプロセスの中心となるものである.

DiscretizePDE

初期化された偏微分方程式の係数は,離散形式に変換される.

離散化された偏微分方程式を計算する.

DiscretizedPDEDataの表示形式は,系の行列の一次元を示す.

系の行列をすべて抽出する.
剛性行列を可視化する.
系の行列を別々に抽出する.

行列の組立て過程において,有限要素にアクセスすることは可能である.有限要素の形式は,行列の組立て過程が効率的でなるものであることに注意する.

行列の組立て中に計算された有限要素を抽出する.
DiscretizeBoundaryCondition

DiscretizeBoundaryConditionsは,与えられた境界条件の離散化バージョンを計算する.後のステップで,離散化された境界条件が実際の系の行列に配備される.

初期化された境界条件を離散化する.

DirichletConditionNeumannValueの境界指定は,違った形で表示される.NeumannValueからの寄与は,系の負荷行列と系の剛性行列に対してである.

一般化されたノイマン(Neumann)境界値からの成分.

離散化された境界条件の"StiffnessMatrix"成分は,NeumannValueが従属変数に依存する場合には非零の寄与を持ち,それ以外の場合にはゼロである.言い換えれば,ノイマン値が一般化されたノイマン値である場合には,"StiffnessMatrix"成分への寄与がある.

DirichletCondition表示は,位置を含む行列を通して与えられる.これらの位置は,条件にどの自由度が適用されるべきかを指定する.返される2つ目の行列は,適用される値を含む.

DirichletConditionの境界条件からの成分.

解の段階

DeployBoundaryConditions

離散化された境界条件が効果を持つためには,系の行列に配備されなければならない.メモリを節約するために,これはその場で行われ,新しい行列は生成されない.

境界条件をその場に配備する.

DeployBoundaryConditionsは属性HoldFirstを持つので,系の行列はコピーをまず作成するということをせずに,その場で変更されることがある.このことはしかし,DeployBoundaryConditionsへの入力が,実際の行列としてではなく値として行列を持つシンボルで与えられなければならないことを意味する.変更されていない系の行列が後で必要となる場合には,DiscretizedPDEDataから再度抽出することができる.言い換えれば,変更されていない系の行列がもはや不要であり,コンピュータのメモリ消費が問題となるような場合には,DiscretizedPDEDataを例えばClearAllを使ってクリアする必要がある.

境界条件の配備中に,NeumannValueからの寄与が剛性行列と負荷行列に加えられる.その後,系の行列は,DirichletConditionが満たされるように変更される.

境界条件の配備の後,系の行列が変更される.系の行列は通常より少ない要素を含む.

一般に,境界条件の配備の後でのみ,方程式の系を解くことができるということに注意する.

方程式を解く

方程式の系は,LinearSolveを使って解くことが可能である.

方程式の系を解く.

後処理

後処理のステップとして,補間関数をLinearSolveから得た解から作成することができる.この解のデータは求めた解でアップデートされ,ProcessPDESolutionsがその後InterpolatingFunctionを構築する.

解で解のデータ sd をアップデートする.
InterpolatingFunctionオブジェクトを作成する.

ProcessPDESolutionsは,複数の従属変数についても使えるという点が利点である.この場合には,ProcessPDESolutionsを使う代りに,ElementMeshInterpolationを使って補間関数を直接作成することができる.ここで第1引数はメッシュ,第2引数はメッシュノードにおける値である.

ノードでの値で補間関数を作成する.
補間関数を可視化する.

通常,補間関数はその領域外でクエリが行われた場合には,extrapolate外挿する.この動作を変更することは可能である.

外挿を使って領域外で評価する.
補間関数がIndeterminateを返すように設定することで,補間関数の自動外挿をオフにする.

"DefaultExtrapolationHandler"のシステムオプションを設定すると,デフォルトの動作を変更することができる.

補間関数からElementMeshを抽出することも可能である.

補間関数からElementMeshを抽出する.

非線形の偏微分方程式

以下のセクションでは,非線形の定常偏微分方程式の解について話し合う.これは複数のレベルで行う.まず,トップレベルの例について紹介し,これをどのように使うのか,どのように実装するのかについてさらに深く見ていくことにする.

トップレベルの例

非線形の偏微分方程式を解く過程を説明するために,特定の境界を持つ領域上での最小曲面を求める必要がある.この最小曲面は,次の非線形方程式を解くことによって求めることができる.

係数は,従属変数 を含み,特に従属変数 の勾配を含むので,非線形である.非線形の偏微分方程式を公式化するためには,方程式1に詳しく示されている

の任意の係数は, および のような一次偏微分係数の両方に従属できる.しかしこれらは, のような一次より高次の微分係数に従属することはできない.

一般に,解法として有限要素法を宇買う非線形方程式は,NDSolveを普通に使う場合のように設定することができる. NDSolveは,その引数をホールドしないので,入力方程式は,NDSolveに到達する前に評価される.必要以上に早く評価されることを避けるために,偏微分方程式はInactiveの形式で指定することができる.Inactiveでラップされていない非線形方程式は,好ましくない不要に高次の導関数を生み出すことがあり,そのような導関数はNDSolveで解くことはできない.特に,この場合のように,従属変数の導関数が拡散係数 内に存在する場合には解けない.そのような場合には,偏微分方程式を非アクティブ化しなければならない.偏微分方程式の非アクティブ化は,偏微分方程式を指定する最も一般的な方法であり,どのメソッドを使っていいのか分からない場合には,このメソッドを選ぶべきである.

非線形方程式を設定する.
領域とその領域の境界条件を設定する.

FindRootによく似て,NDSolveは定常非線形偏微分方程式を有限要素法で解く場合には,それぞれの従属変数に対して,初期シード値が必要である.初期シード値は,InitialSeedingで指定することができる.

初期シード値を指定して非線形方程式を解く.

初期シード値が指定されていない場合には,ゼロであると見なされる.

解を可視化する.

ここで方程式が非アクティブな形式ではない場合には,NDSolveはこの偏微分方程式を解くことができない.

アクティブな形式の非線形方程式を解くと,以下のメッセージが発せられる.

非アクティブな偏微分方程式についての詳細は,「形式的な偏微分方程式」についてのセクションを参照されたい.

以下のセクションでは,上のような非線形方程式が与えられたときにNDSolveが内部で行う解の一般的な過程を紹介する.線形方程式の場合と同じ段階が非線形方程式の場合にも使われる.

初期化段階

非線形方程式の初期化段階は,線形の場合とほとんど同じである.非線形の場合に唯一異なるのは,FindRootの場合と同じように,初期シード値が指定されなければならない点である.これは,解データのコンポーネント"DependentVariables"を初期シードに設定することで行う.

数値領域と変数データを設定する.
初期シード値をとして解データを設定する.

この場合,初期シード値は,数値でも独立変数である に基づく任意式でもよい.かなり非線形な問題の場合には,あまり非線形ではない問題を解いて,結果として得られるInterpolatingFunctionをより高度な非線形問題の初期シード値として使うという方法がある.初期シード値が指定されていない場合には,が初期シード値として使われる.

メソッドデータ,境界条件,偏微分方程式の係数を指定する.

離散化と解の段階

非線形方程式を解くことは,インタラクティブな過程である.つまり,最下層レベルでは,インタラクティブな過程にFindRootが使われる.まず最初のステップは,偏微分方程式係数を線形化することである.この過程の詳細は下で説明する.その後初期シード値を使って,線形化された偏微分方程式を評価する.これでLinearSolveで解ける線形方程式系が作成される.この過程の結果は,偏微分方程式に代入し直され,アップデートされた,線形化された偏微分方程式は解の新しい近似値で離散化される.この過程は,理想的には解が収束するまで繰り返される.

関数PDESolveFindRootの呼出しを調整し,線形化と繰返しの離散化を行い,解を含む解データオブジェクトを返すか,アルゴリズムが収束に失敗した場合にはPDESolveが求めることができた解の最良の近似値を返すかする.

PDESolveを使って非線形偏微分方程式を解く.

求まった解は,返された解データオブジェクトに保存される.

後処理

InterpolatingFunctionオブジェクトは,解データから作成できる.

ProcessPDESolutionsを使ってInterpolatingFunctionオブジェクトを作成する.
解を可視化する.

線形化過程

FindRootはデフォルトでニュートン・ラフソン(NewtonRaphson)法を使って非線形方程式を解く.このメソッドは, 想定される解,初期シード値で始め,方程式の線形化された形式でこれを改善しようとする.

線形化過程をよりよく理解するためには,スカラー関数のためのニュートン法の微分を詳しく見ると役に立つ. の非線形関数とし, を以下が成り立つ未知の解とする.

初期シード値 を取り,それを未知の解 から引くと,結果として剰余 が得られる.

方程式を並べ替えて以下の形にする.

の周りでテイラー展開を行い,以下を得る.

より高次の項を無視すると,以下が得られる.

の微分がゼロではないと想定すると,以下の反復が得られる.

.

について一般化すると,以下が得られる.

2つの代数方程式の例を見てみよう.

.

非線形関数 を設定する.
関数行列式 を計算する.
関数を設定して関数行列式を評価する.
初期シード を指定する.
初期シード からの最初のニュートンステップを計算し,その結果を に保存する.
からの第2ニュートンステップを計算し,その結果を に保存する.

解が求まるまで,同じ過程を続ける.次のステップとして,反復過程をベクトルで与えられるステップを計算する関数で囲む.2つの連続するステップのノルムが閾値より大きい間に,新しい次のステップを計算する.

ニュートンステップを計算する関数を設定する.
非線形問題を解く.
計算した結果が希望する確度で非線形方程式を解くことを確かめる.

与えられた例で1点気を付けなければならないことは, と関数行列式 の両方が一定の部分 に従属する部分があることである.有限要素法を使って非線形偏微分方程式を解く場合には,反復過程の前に一定部分を一度計算し,それらを各ステップで計算される非一定部分に加えるとよい.

方程式1で得た係数形式の非線形偏微分方程式の初期化過程も同じステップで行われる.方程式2から始めて, を設定して時間独立係数形式を得る.

使いやすくするために,方程式を以下のように並べ替える.

係数の名前を

および

に変える.これで時間独立係数形式は以下のように書き換えられる.

線形化を行うためには,上のスカラーの場合と同じように行う. および の関数であり,以下が得られる未知の解である.

任意のシード値 を取り,これを未知の解 から引くと,以下のように剰余 が得られる.

方程式を並べ替えて以下の形にする.

の周りでテイラー級数の展開を行い,以下を得る.

より高次の項を無視すると,以下が得られる.

および に対して, および についての微分が計算される.

線形化された偏微分方程式の実際の解は,方程式3で与えられた係数形式を,任意数の従属変数の以下の係数と一緒に利用する.

ここで は整数である.から (従属変数の数)までである. 個のスカラー方程式に対して が得られる.

ニュートンの求根法とそれをどのように偏微分方程式に適用するかについての説明の後は,方程式

をより深いレベルで解くことを続けて行い,「離散化と解の段階」 のセクションで紹介したPDESolveへの呼出しがどのように使えるかを示す.最初のステップでは,線形化された係数が構築される.

従属変数と独立変数を取り出す.
従属変数と独立変数の数に対して変数を設定する.
従属変数のリストを構築する.
従属変数の微分のリストを構築する.
任意の非線形性も含めて,PDECoefficientDataオブジェクトから係数を取り出す.
方程式4から を構築する.
方程式5から を構築する.

方程式6についての偏微分係数 は,Dを使って計算することができる.

線形化された偏微分方程式の拡散係数 を調べる.

方程式7についての偏微分係数 は,Dを使って計算することができる.

線形化された偏微分方程式の保守的な対流係数 を調べる.

方程式 8 についての偏微分係数 は,Dを使って計算することができる.

線形化された偏微分方程式の対流係数 を調べる.

方程式9 についての偏微分係数 は,Dを使って計算することができる.

線形化された偏微分方程式の反応係数 を調べる.

これらの係数は,ユーティリティ関数LinearizePDECoefficientsを使って,便利なことに有限要素のコンテキストから求めるとこができる.

線形化された偏微分方程式の係数を調べる.
線形化された拡散係数を調べる.
線形化された負荷微分係数を調べる.

線形化された偏微分方程式を解く

一旦偏微分方程式が線形化されると,係数は分割される.これは,求根アルゴリズムの異なる部分で係数が必要となるからである.係数形式

は以下のようにグループ化し直される.

左辺のすべての項は,剛性行列にまとめられる.剛性行列は偏微分方程式の関数行列式である.剛性行列は,上のスカラーの例からの の項に類似するものである. の接線であるので,剛性行列も接線剛性行列と呼ばれる.偏微分方程式の右辺の項は, に対応する.PDECoefficientDataの分割は,関数SplitPDECoefficientsで行う.

線形化された偏微分方程式の係数を分割する.

係数 はこれで変数linLoadPDEC内にあるが,残りはlinStiffnessPDEC内にある.減衰行列と質量行列の係数は,ここで扱っている問題については存在しないので,linDampingPDECおよびlinMassPDECは空である.

メッシュの初期シード値を評価する.
線形化された偏微分方程式の定常成分を離散化する.
定常境界条件を離散化する.
右辺の関数を指定する.

ここで与えられた関数femRHSは,ここで取り扱っている例には存在しないが,一般化されたノイマン値も取り扱う.実際,ここで与えられた関数は,実際の実装に大変近いものである.

関数femRHSでは負荷ベクトルのみが利用されることに注意する.これは,SplitPDECoefficientsコマンドで生成されるlinLoadPDECには,負荷係数 と負荷微分係数 のみが含まれるからである.

一旦nonlinearおよびnonlinearLoadの変数が必要なくなると,これらはNullに設定され,これらが指すメモリが解放される.

関数を指定して関数行列式を計算する.

ここでも,関数femJacobianは実際の実装に大変近いので,もはや必要がない変数はNullに設定される.

関数行列式の構築には,線形化された剛性成分が使われる.これには,拡散係数 ,保守対流係数 ,対流係数 ,反応係数 が含まれる.これらはすべてSplitPDECoefficientsコマンドで生成されたlinStiffnessPDEC PDECoefficientDataに保存される.

実装では,線形である係数はすべて一度だけ集められることに注目されたい.関数femRHSおよびfemJacobianのどちらでも,非線形係数のみが計算され,先に計算された線形成分に加えられる.これで各反復において変化しない線形係数を無駄に計算せずに済む.最後の非線形と線形の境界条件が配備される.

境界条件を初期シード値に配備する.

最後に,線形化された偏微分方程式の解を求める.このためには,アフィン共変ニュートン法が使われる.

FindRootを使って非線形偏微分方程式を解く.
解データをアップデートする.

InterpolatingFunctionオブジェクトが解データから作成できる.

ProcessPDESolutionsを使ってInterpolatingFunctionオブジェクトを作成する.
解を可視化する.

前に説明した手順は,PeriodicBoundaryConditionが存在する場合にはうまく作動しないという点において制約されている.

非定常偏微分方程式

非定常偏微分方程式は,時間依存の偏微分方程式である.最も簡単な非定常偏微分方程式は,偏微分方程式の一部として,依存変数の時間導関数を持つ.この場合は,離散化された時間導関数を表す追加の系の行列が作成される.その後,系の行列の時間積分を,例えばNDSolveを使って行うことができる.それ以外の点では,ワークフローは定常偏微分方程式の場合とほぼ同じである.より複雑な例には,依存する係数および/あるいは境界条件が含まれる.後ろのセクションでそれらを時間積分する方法について述べる.

定常係数と定常境界条件を持つ非定常偏微分方程式 はじめに

最初の例として,1Dの熱方程式を見てみよう.空間領域が0から1で,時間領域が0から1である,次のモデル偏微分方程式を考えてみる.

境界条件と初期条件はどこでも0である.

直接NDSolveを使って偏微分方程式の解を求める.

この場合,メソッドオプションは,空間変数 が有限要素法で離散化されるべきであり,時間積分が線分法で行われることを指定する.このアプローチは,半離散化と呼ばれる.反対に,完全な離散化は,時間変数 と空間変数の両方を離散化する.

数値解をプロットする.

偏微分方程式の解を求める過程は,定常の場合に似ている.その大きな理由として,偏微分方程式の係数も境界条件も時間に依存しないことが挙げられる.非定常系の解を求めるために,偏微分方程式を半離散化する方法が取られる.偏微分方程式は空間()内で,有限要素法を使って離散化される.そうすると,その有限要素法の離散化の結果は,時間変数()における常微分方程式系である.NDSolveを使って,その常微分方程式系が時間積分される.

まず,変数データが作成され,加えられる.

従属変数と独立変数の名前を持つ変数データを作成する.

次に,解のデータを初期時間とメッシュで設定する.

NumericalRegionを指定する.
"Space""Time"の成分集合を持つ解のデータを作成する.

偏微分方程式係数,境界条件,メソッドデータが初期化される.

偏微分方程式係数を初期化する.
境界条件を初期化する.
変数と解のデータで有限要素データを初期化する.
NumericalRegionからElementMeshを抽出する.

次に離散化を行う.非定常系が考慮されているのだが,偏微分方程式の係数も境界条件も時間に依存しないので,離散化は定常成分のみを持つ.

離散化された偏微分方程式を計算する.
初期化された境界条件を離散化する.
系の行列をすべて抽出する.
境界条件を配備する.

次に初期条件が設定され,NDSolveは離散化された方程式系の時間積分を行う.

境界条件に基づいて,初期条件を設定する.
NDSolveで方程式系の時間積分を行う.

下の方のセクションで,関数行列式を指定することによって時間積分の過程を速める方法について説明する.

NDSolveが解の過程のすべての面を制御していた上のNDSolveへの呼出しから自動的に選択されたオプションにマッチするように,NDSolveオプションを選ぶ.

時間 が与えられると,補間関数を構築して記憶する関数を設定する.
自動と手動の解の違いを可視化する.

定常係数と定常境界条件を持つ非定常偏微分方程式

以下のモデル偏微分方程式を考える.

境界条件は,では であり,では である.マーカーについての詳細は,「要素メッシュの生成」のチュートリアルを参照されたい.

偏微分方程式は,一階の時間導関数を持つが,偏微分方程式の係数 と境界条件は時間依存ではない.偏微分方程式を解く1つの方法は,まず偏微分方程式と境界条件を離散化してから,結果の常微分方程式を時間積分器としてNDSolveで時間積分する方法である.この過程は,NDSolveが内部で非定常偏微分方程式を解く方法であり,これについては次のセクションで見る.

この例では,マーカー付きの境界メッシュが事前定義されている.形状は,Scanning Electrochemical Microscopy (SECM)プローブをモデル化する.このプローブは,7つの電極が円状のグラス構造に埋め込まれているので,七極管である.

事前定義された境界メッシュを使う.
BoundaryMeshRegionを境界ElementMeshに変換する.

BoundaryMeshRegionオブジェクトとElementMeshオブジェクトの間の関係についての詳細は, 「要素メッシュの生成」チュートリアルを参照のこと.

事前定義された境界メッシュ領域では,整数マーカーが境界の表面と点に設定されている.これらのマーカーは,境界条件が解のプロセス中に適用される位置として可視化し,使うことができる.

点マーカーをグループに従って可視化する.

この境界条件はディリクレ条件だけを含むので,点の要素とマーカーだけが重要である.ディリクレ条件は常に点要素に作用する.境界メッシュも実装された境界要素についてマーカーを持つ.これらは,ノイマンタイプの境界値を指定するのに必要であり,これらの値は常に境界要素に作用する.

最適化されたノード接続性を持つ完全な一次メッシュを作成する.

厳密に言うと,ノードの並び替えはオプションである.ノードの並び替えメソッドのオプションがTrueに設定されている場合には,系の行列が持つバンド幅は最小化される.低いバンド幅は,反復ソルバに有益である.ノードの並び替えは1度行うだけでよく,ソルバが各時間ステップに使われるので,非定常問題にノードの並び替えを使うことが役に立つ場合が多い.

変数データを設定する.

この場合,ElementMeshオブジェクトがNumericalRegionオブジェクトを作成するために直接使われたかもしれない ImplicitRegionの指定は存在しない.

解のデータを設定する.

解のデータを追加する."Time"成分は,初期時間t0であり,積分が始まる."DependentVariables""TemporalDerivatives"の成分を設定すると,初期条件と初期条件の導関数が設定される.

偏微分方程式の係数を初期化する.
偏微分方程式の境界条件を初期化する.
メソッドデータを初期化する.
系を離散化する.
系の行列をすべて抽出し,境界条件を配備する.

最初の非定常行列である減衰行列が追加されたことに注意する.これは後で時間積分を行う際に必要となる.

剛性行列を可視化する.

剛性行列の構造に注意する."NodeReordering"->Falseをメッシュ生成の過程に使うと,あまり構造化されていない系の行列が返される.

非定常の解を計算する前に,定常の解を調べる.

定常の解を求める.
または,反復ソルバを使って定常の解を求める.
曲面プロットとして,解を可視化する.
領域外でクエリを行い,それを警告しない場合に, Indeterminateを返す補間関数を構築する.
x 方向にContourPlotを作成する.
モデル内に3D等高線プロットとして解を可視化する.

偏微分方程式の離散化が常微分方程式の系で起る.常微分方程式の時間積分は,NDSolveで行われる.これには,境界条件にマッチする初期条件を設定する必要がある.

境界条件に基づく初期条件を設定する.
初期条件を可視化する.

NDSolveの時間積分の速度を上げるために,オプションを使って関数行列式の疎なパターンを指定することが可能である.系が全体としてNDSolveを通して解かれる場合には,これは不要である.ここでは,NDSolveがベクトル変数の構造を読み解こうとはしないため,必要である. 関数行列式の疎なパターンは,減衰行列のパターンを通して使用可能である.

減衰行列と剛性行列に基づいて関数行列式の疎なパターンを設定する.

疎性は,境界条件が配備されない場合にも使える.

NDSolveを時間積分器として,非定常偏微分方程式の解を求める.

前のバージョンのNDSolveでは,関数行列式が"Jacobian{Automatic,Sparsesparsity}"と指定されていた."sparsity" は,減衰行列の疎なパターンであった.バージョン12.2以降では,剛性行列と減衰行列の両方に対して,疎なパターンの指定が優先されるようになった.古い関数行列式の設定も引き続き使用可能であるが,時間積分がゆっくりになるため,剛性行列と減衰行列の両方に対して関数行列式を指定することに切り替えることを強く推奨する.

EvaluationMonitorを使うと,解の過程が少しゆっくりになるが,解の過程をライブで管理する事ができる.

解は,非定常の系がどのように上で計算された定常状態に到達するかを示す.

さまざまな時間ステップで定常状態への到達を可視化する.

時間についての単一の補間関数と空間変数を構築することができる.このため,非定常補間が設定された時点のステップが抽出されなければならない.

補間関数が構築された時点を抽出する.

次に,補間の値を設定する.それぞれの時間のステップについて,ノードにおける値のリストが必要である.

それぞれの時間のステップにおけるソリューションの値を抽出する.

実際の補間関数を構築する.

単一の補間関数を構築する.
すべての時間のステップについて空間における点で単一の補間関数を調べる.

定常係数と定常境界条件を持つ非定常偏微分方程式のモデル次数の簡約

系の行列へのアクセスが重要である非定常の場合,モデル次数の簡約を考慮する.モデル次数の簡約が問題とするのは,非定常方程式系において、もとの特性のほとんどを保ちつつ,もとのものより小さい同等の方程式系が存在するかどうかということである.

形式の の自由度で,離散化された方程式系を考慮する:

減衰行列 および剛性行列 × の大きさで,負荷行列は ×1の大きさである.今度は, の自由度のずっと小さい2つ目の方程式系を設定する:

このとき,である.

2つ目の自由度のため,簡約された系 は,もとの系 よりもずっと小さく,時間積分はかなり効率的である.以下のアルゴリズムは, の間でマップする射影ベクトル を設定する:

ここで の大きさは の大きさは ,そして は近似の相違である.

射影ベクトル を簡単に説明すると,非零の開始ベクトル と行列 を与えられた場合に,アルゴリズムは以下の式のようになる正規直交 を求める:

モデル次数の簡約の射影コード.

余談だが,プラットフォームによっては,Method"Pardiso"LinearSolveのオプションとして使うことができるものもある.このソルバは非常に効率的であり,非正則行列に近い行列もうまく処理することができる.

手続きは,簡約された系の次元にもとの系の次元を掛けた射影ベクトル を返す.

たった の自由度を持つ近似の方程式系を求める.
簡約された系の行列を抽出し,その次元を調べる.
簡約された系の時間積分を効率的に行う.
完全な系と簡約された系の ノルムを比べる.

非定常係数を持つ非定常偏微分方程式

今度は,1つあるいは複数の係数が時間依存である非定常偏微分方程式について見てみよう.時間依存負荷 を持つ次のモデル偏微分方程式を考えてみる.

これは,上の例と同じ定常境界条件を持つ.この例では,以下の負荷関数を使う.

変数と解のデータも上の例と同じである.有限要素法のメソッドデータも再利用できる.

偏微分方程式の初期化は,追加の成分である負荷係数を持つようになった.これは時間変数 に依存する.

偏微分方程式の係数を初期化する.

初期化された係数は,指定の係数を保存する.

保存された係数を抽出する.

効率的な時間積分を行うために,係数は定常係数と非定常係数に分類される.定常係数は上の例と同じように,1度だけ評価される.時間依存成分を持つ係数は,時間ステップから時間ステップまで再評価され,定常離散化に追加されなければならない.

偏微分方程式の係数のInitializationは,すべての係数の分類を内部で行う.この分類は,チェックして利用することができる.それぞれの分類(例えば,定常係数)は,位置を含む2つのリストで表される.最初のリストにはスカラー係数が含まれ,2つ目のリストにはベクトルと行列の係数が含まれる.これをさらに詳しく見るために,一般的な偏微分方程式の形式を考えてみよう.

係数 はスカラー, はベクトル,× 行列である.返された2つのリストのうち最初のものは,スカラー係数(),2つ目はベクトルと行列の係数()を表す.

初期化された非定常偏微分方程式の係数の位置を抽出する.

モデルの偏微分方程式には,非定常成分が1つ含まれる.

すべての係数から非定常スカラー係数を抽出する.

この例には,非定常の行列成分は含まれていない.

すべての係数から非定常の行列係数を抽出する.

"Stationary"係数は,指定されず,ゼロであるものも含めて,定常であるすべての係数を表す.

定常と非定常の成分に分類することによって,偏微分方程式の定常と非定常の離散化を別々に効率的に計算することが可能になる.

偏微分方程式の定常成分を離散化し,保存する.

上の例では,文字列の引数"Stationary"はデフォルトであるが,必須ではない.離散化された定常偏微分方程式の負荷ベクトルは空である.負荷ベクトルには非定常の成分しか含まれていない.

系の非定常成分を特定に時点で離散化する.

離散化された非定常偏微分方程式の係数の負荷ベクトルが追加される.時点 は,数値3で置き換えられている.その他の系の行列には非定常係数が含まれていないので,加えられていないことに注意する.

偏微分方程式の係数を分類するのは,定常係数はいったん離散化し,保存してから,時間積分のループで再利用することができるからである.つまり,時間変数 の依存を示した係数だけを各ステップにおいて離散化すればいいということである.これによって,計算にかかる多くの労力を節約することができる.重要な幾何情報が事前に計算され,FEMMethodDataオブジェクトに保存されているので,それぞれの非定常成分の計算時間も最小に抑えることができる.

離散化された偏微分方程式の境界条件を再利用する.
系の行列を計算し,定常境界条件を再配備する関数.

右辺の計算には,すでに計算された偏微分方程式の定常係数が利用される.この例では,元帥行列は非定常成分を持たず,因数分解はキャッシュできる.減衰行列が時間依存の成分を持つ場合には,それぞれの時間ステップで,減衰行列がアップデートされ,LinearSolveステップがキャッシュされないので,上の右辺の計算は変更されなければならない.すると,因数分解が各時間ステップにおいて行われなければならない.

非定常境界条件がないので,初期条件は上と同じように設定される.

境界条件に基づいて初期条件を設定する.

疎なパターンは,また減衰行列から導き出して設定される.この場合は,境界条件は減衰行列に配備されていない.

剛性行列と減衰行列に基いて,関数行列式の疎なパターンを設定する.
NDSolveを時間積分器として,非定常偏微分方程式の解を求める.

可視化は,上と同じように行われる.

さまざまな時間ステップにおいて定常状態への到達を可視化する.

非線形の非定常係数を持つ非定常偏微分方程式

このセクションでは,1つあるいは複数の係数が非線形である非定常偏微分方程式について見てみる.非線形の減衰行列 と非線形の拡散係数 を持つ以下のモデル偏微分方程式を考えてみよう.

非線形の熱方程式が特定の例として使われる.温度 ,密度 ,温度依存の熱容量 ,熱依存の伝導性 で以下の式を使う.

ここでは,物事を単純に保ち,有限要素法のプログラミング面に出来るだけ焦点を置くために,位相変化は一次元領域でモデル化され,単純性のためにパラメータが選ばれている.

位相変化の物理特性についての詳細は,位相変化がNDSolveレベルでどのようにモデル化できるかを説明する熱伝導の偏微分方程式モデリングのチュートリアルも参照されたい.

上の偏微分方程式を使って,位相変化をモデル化する.1つの位相から次の位相への変化をモデル化するのによく使われるのは,平滑な推移関数を使う方法である.

平滑化されたステップ関数を作成する.
位相1から位相2への,平滑化されたステップ関数の推移を可視化する.
非線形の係数とパラメータを設定する.
1D領域を指定する.

左の境界には一定流入量のエネルギーが存在し,右の境界には固定の温度が設定されている.時間の経過とともに,ものの初期位相は領域の左から右の端に変化する.

境界条件を設定する.

初期化段階は,上のセクションに示されるのと同じ方法で行われる.

数値領域,変数,解データを設定する.
偏微分方程式の係数,境界条件,メソッドデータを初期化する.

非線形係数 は直接使えることに注意する.

NDSolveの方程式を設定するために,今度も上のセクションとほぼ同じ方法を取る.ここでは,上のセクションで と呼ばれたものを拡張し,同時に非線形係数を計算して非線形の減衰行列を含む.これにより,すべての偏微分方程式係数が1つの関数(ここでは微分代数方程式ソルバの剰余を計算するため,discretizePDEResidualと呼ぶことにする)に含まれるようになる.

上と同じように,定常成分を計算することから始める.

偏微分方程式の係数と境界条件の定常成分を計算する.

離散化関数は,時間 ,従属変数 ,従属変数 の導関数を引数として得る.次に,時間従属偏微分方程式の係数と境界条件が計算され,定常成分に加えられる.その後,非線形成分が計算され,これも加えられる.そうすると,関数は離散化されたバージョンを返す.

関数を作成して,各時間ステップで離散化された偏微分方程式を計算する.
値がどこでもに設定されるように初期条件を評価する.

疎なパターンを計算するために,ヘルパー関数を利用して,特定の時点における減衰行列と剛性行列を構築する.積分時間全体を通して有効である減衰行列と剛性行列の疎なパターンを作成する必要がある.したがって,任意の時点において構築された行列の疎なパターンを抽出し,その時点において,可能な時間従属係数が確実にアクティブであるようにする.しかし,ここでは,剛性行列と減衰行列のパターンだけが問題であるので,使われた厳密な非零の時間は重要ではない.

特定の時点において,剛性行列と減衰行列を構築するヘルパ関数.
時間10^-6における疎なパターンを計算する.
方程式を解く.
時間ステップからInterpolatingFunctionを作成する.
解を可視化する.

熱方程式と適用可能な境界条件についての詳細は,「熱移動」のチュートリアルおよび「熱移動の偏微分方程式と境界条件」のガイドページに記載されている.

積分係数を持つ非定常偏微分方程式

非線形の非定常係数を持つ非定常偏微分方程式のセクションで提示されたメソッドは拡張して,非定常の積分偏微分方程式を解くことができる.各時間ステップでElementMeshInterpolationを構築し,NIntegrateを利用して現行の補間関数上で積分するというアイディアである.このアプローチは,積分だけに限らない.実際,数字を与えるどのような演算でも補間関数に使える.

例として,以下の非定常積分微分方程式を考える.

方程式を拡張する.

積分方程式と微分方程式の領域はどちらも0から5までである.

1D領域を指定する.
領域の左端と右端についてDirichletConditionを指定する.

初期化段階は,上のセクションで示される方法と同じように行われる.

数値領域,変数データ,解のデータを設定する.

積分係数以外の係数をすべて初期化する.

偏微分方程式係数,境界条件,メソッドデータを初期化する.

上と同様,定常成分を計算することから始める.

偏微分方程式係数と境界条件の定常成分を計算する.

discretizePDEResidualと呼ばれる離散化関数で,線形成分をロードし,非定常で非線形の成分を加えてから,現行の補間関数を構築する.次に,構築された関数上で数値積分を行い,積分を含む,見付からない偏微分方程式係数を構築する.見付からない係数を離散化し,その結果を系の行列に加える.そうすると,関数は離散化されたバージョンを返す.

各時間ステップでElementMeshInterpolationを作成するために,メッシュをNumericalRegionから抽出する.

ElementMeshを抽出する.
各時間ステップで離散化された偏微分方程式を計算する関数を作成する.

厳密には,この場合,積分係数は線形であるので,積分係数の非定常で非線形な成分の離散化は不要である.一般化のためにコードに非定常で非線形な積分係数の計算も含む.

値がすべてにおいてに設定されるように初期条件を評価する.

上と同じ方法で疎なパターンを計算する.

特定の時点における剛性行列と減衰行列を構築するヘルパ関数.
時間10^-6における疎なパターンを計算する.
方程式を解く.
時間ステップからInterpolatingFunctionを作成する.
解を可視化する.

結合偏微分方程式

このセクションでは,結合偏微分方程式がどのように解かれるかについて説明する.基本的には,ワークフローは単一の偏微分方程式と同じである.違いを下に述べる.

負荷による梁の変形

例として,構造力学の問題として,梁の変形について考える.変形は,2つの従属変数における結合偏微分方程式としてモデル化できる.

結合偏微分方程式としての平面応力の形成.

平面応力の演算子では,YはYoungの係数であり,ν はポアソン(Poisson)の率である.これらは物質の特性である.

あるいは,関数SolidMechanicsPDEComponentを使って,偏微分方程式の成分を生成することもできる.SolidMechanicsPDEComponentの関数ページにもテキスト形式での方程式が含まれている.

変形される梁の長さは5単位で,高さは1単位である.

明示的な境界ボックスとメッシュの可視化で,メッシュ領域を作成する.

境界条件は以下の通りである.左辺では,梁は xy の両方向に固定されている.右辺の境界では,1単位の下向きの圧力が適用される.NeumannValueは,負の y 方向に梁にかかる圧力を指定することによって境界負荷を強制する.

領域の左右で境界条件を設定する.

初期化は,上のセクションでの場合と同じ成分を持つ.まず,変数と解のデータが設定される.

変数データに2つの従属変数 ,および独立変数名を加える.

結合偏微分方程式については,2つ以上の従属変数が与えられなければならない.

解データの"Space"成分を設定する.

メソッドデータの初期化は,単一の従属変数の場合と同じようになされる.

有限要素データを変数と解のデータで初期化する.

FEMMethodDataの表示形式は,自由度,それぞれの従属変数の補間次数,そして積分次数を示す.

方程式系の自由度は,要素メッシュのノード,従属変数の数,および従属変数の補間次数を組み合せた結果である.

自由度は,従属変数の数とメッシュ内のノードに関係している.

結合偏微分方程式の場合は,係数がもう少し関係してくる.結合偏微分方程式は,基本的に,単一の偏微分方程式が集まったものである.しかし,それぞれの成分がそれぞれの従属変数に対して存在する.次の2つの従属変数 の偏微分方程式を考えてみよう.

結合偏微分方程式はで定義される.係数 は,それぞれの偏微分方程式でスカラーである.各偏微分方程式において, は長さ のベクトルであり,× 行列である.InitializePDECoefficientsの入力係数は,上の結合偏微分方程式とモデル化される結合偏微分方程式との対応によって求められる.

モデル偏微分方程式の係数.

結合偏微分方程式の場合には,各偏微分方程式の各従属変数について拡散係数行列がある.

偏微分方程式の係数を初期化する.

PDECoefficientDataの表示形式は,初期化された係数について,系の大きさ,つまり従属変数の数(2個)と演算子の大きさ,つまり空間次元(2)を示す.

次に,境界条件を設定する.従属変数のそれぞれには,境界条件のリストが含まれる.これが境界条件を偏微分方程式の系の特定部分に関連付ける.境界条件の最初の集合は,結合系の最初の偏微分方程式に,境界条件の2つ目の集合は結合系の2つ目の偏微分方程式に,という風に関連付けられる.これによって,境界条件を一意的に偏微分方程式の部分に関連付けることが可能である.

境界条件を初期化する.

初期化された偏微分方程式の係数は,離散形式に変換される.

離散化された偏微分方程式を計算する.

DiscretizedPDEDataの表示形式は,一次元の系の行列を示す.

系の行列を抽出する.
剛性行列を可視化する.

剛性行列の可視化は,もとになる偏微分方程式の係数が2つの従属変数の系であることを示す.自由度は,2つの従属変数間で分けられる.

2つの従属変数間でどのようにインシデントが分けられるかを計算する.
初期化された境界条件を離散化する.
境界条件を配備する.

方程式系はこれでLinearSolveを使って解くことができるようになる.

方程式系を解く.

それぞれの従属変数の補間関数を作成することができる.解の形式LinearSolveは, への両方の解を含む1つの解ベクトルであるので,解ベクトルは分けられなければならない.

ノードに値を持つ補間関数を作成する.
補間関数を可視化する.
変形した構造を赤で可視化する.

固体力学と適用可能な境界条件についての詳細は,「固体力学」のモノグラフおよび「固体力学の偏微分方程式と境界条件」のガイドページに記載されている.

揺れる梁非定常の結合偏微分方程式

非定常の結合偏微分方程式には,これ以上はほとんど何も必要ではない.以下の例では,梁を偏向させ,レイリー(Rayleigh)の減衰で動きを弱める方法を示す.レイリーの減衰では,減衰行列は,質量行列と剛性行列の線形の組合せである.2つの従属変数を持つ非定常の結合偏微分方程式についての系の行列の構造は,次のように表される.

ここで,減衰系の行列は,質量行列と剛性行列の組合せである.

非定常のモデル偏微分方程式の追加の係数.

非定常係数も質量行列の係数集合を持つ.

非定常の偏微分方程式の係数を初期化する.
離散化された偏微分方程式を計算する.

DiscretizedPDEDataの表示形式は,一次元の系の行列を示す.

系の行列を抽出する.

減衰行列の疎な配列は空である.これは,質量行列の係数が設定され,ここで使用されるデフォルトの減衰行列の係数がゼロであるからである.

質量行列を可視化する.

質量行列は,ブロック対角が追加され,対角外のブロックがゼロである,偏微分方程式系数の構造を示す.

レイリーの減衰は,質量行列と剛性行列の線形の組合せとして計算される.

レイリーの減衰行列を計算する.
境界条件を配備する.

初期条件と初期条件の導関数は,ゼロの設定される.これは,ディリクレ境界条件がすべてゼロであり,梁が静から始まるので妥当である.

初期条件を設定する.
関数行列式について減衰行列の疎なパターンを設定する.
関数行列式について剛性行列の疎なパターンを設定する.

次のステップは,NDSolveを使って行う実際の時間積分である.シミュレーションのために,単一の記号的変数 を導入する.記号的変数 は, のベクトルを組み合せたものを表し,NDSolveには記号 が初期条件と系の行列の大きさの両方から, の長さを組み合せたベクトルを表さなければならないことが分かっている.言い換えれば,それぞれの自由度について変数を導入する必要はない.

NDSolveを時間積分器として,非定常偏微分方程式の解を求める.

それぞれのタイムスタンプで同じ量の変形を表示するには,"ScalingFactor"を1に設定する.

揺れる梁を可視化する.

固体力学と適用可能な境界条件についての詳細は,「固体力学」のモノグラフおよび「固体力学の偏微分方程式と境界条件」のガイドページに記載されている.

揺れる動的に負荷された梁

このセクションでは,時間依存のボディ負荷 と時間依存の境界条件を加えることによって,上の例をさらに1ステップ発展させる.

そのためには,変数と解のデータは時間変数集合を持つ必要がある.

変数データの時間変数を設定する.
解のデータの最初の時間を設定する.
メソッドデータを初期化する.

新しい偏微分方程式モデルとして,まずボディ負荷に働く項

を加えて,正弦波運動での x 方向への梁の膨張と収縮をモデル化する.

非定常モデルの偏微分方程式のための追加の時間依存のボディ負荷.

次に,NeumannValueを通して圧力を指定することによって,y 方向に適用される時間依存の境界の圧力を加える.

時間依存の境界負荷.

偏微分方程式の係数と境界条件を初期化する.

非定常偏微分方程式の係数を初期化する.
境界条件を初期化する.

離散化された定常偏微分方程式の係数と境界条件を事前に計算する.

離散化された定常偏微分方程式を計算する.
離散化された定常境界条件を計算する.

方程式を少し並べ替える.

現行時間 の値を得るためのヘルパ関数を書く.これらは1つのベクトル および に結合される.現行時間は解のデータで更新される必要がある.定常の系の行列が離散化された偏微分方程式から抽出される.これらはそれぞれの時間ステップで変化しないので,再度計算されることはないということを覚えておくことが重要である.次に,偏微分方程式の非定常部分と境界条件が離散化される.非定常と定常の系の行列が足し合わされ,レイリー減衰が計算される.定常と非定常の境界条件はどちらも配備される.関数は離散化された を返す.

非定常の系の行列を計算し,境界条件を配備する関数.

初期条件と初期条件の導関数は零に設定される.ディリクレ条件はすべて零であり張りは残りから始まるため,この設定は妥当である.

境界条件に基づく初期条件を設定する.

疎なパターンの指定には,偏微分方程式のの定常離散化を使う.定常コンポーネントは負荷ベクトルにのみ存在するので,これは妥当である.言い換えれば,質量行列と剛性行列に時間依存のコンポーネントがない場合にのみ以下は正しい.

関数行列式について減衰行列の疎なパターンを設定する.
関数行列式について剛性行列の疎なパターンを設定する.
NDSolveを時間積分器として,非定常偏微分方程式のの解を求める.
動的に負荷された揺れる梁を可視化する.

固体力学と適用可能な境界条件についての詳細は,「固体力学」のモノグラフおよび「固体力学の偏微分方程式と境界条件」のガイドページに記載されている.

大規模な有限要素法分析

NDSolveで実装される有限要素法は,スピードの点で最適化されている.この最適化により失われるものもある.数値計算法では,スピードとメモリの消費が両立することがないことはよく知られている.つまり,実装が速ければ速いほど,手元の問題を解くのにより多くのメモリを必要とする.しかしながら,実装された有限要素法をさまざまな方法を通して調整し,偏微分方程式の解を求めるのに追加の時間を費やすことによって,より少ない量のメモリを使うようにすることは可能である.「有限要素法で偏微分方程式を解く」にはいくつかの方法が提案されている.以下のセクションでは,有限要素解法で必要なメモリ量を減らすためのさまざまなその他のオプションについて述べる.これらの提案には,一般的なものから特定のオプションまでいろいろある.

クリーンスタートでカーネルを始める.

一般にWolfram言語は,先に計算された結果への参照を保存する.この仕組みにより,将来の計算には不要であるデータが累積されることになってしまう.履歴の長さを制御する仕組みを設定し,前の出力を保存しないようにすることが可能である.

履歴の長さを0に設定して,前の出力を保存しないようにする.

有限要素の離散化における一般的な手順は以下の通りである.