「純粋」と「応用数学者」の関係は、信頼と理解に基づいています。 「純粋な数学者」は「応用数学者」を信頼せず、「応用数学者」は純粋な数学者を理解していません。
しばらく前、私は、浸漬境界法が偉大で強力なもので記述されるようなアクセス可能な資料を見つけることができなかったという事実に直面しました。 要するに、これは計算流体力学の手法であり、形状と動力学が非常に複雑なオブジェクトの周りの流れを計算することができます。 したがって、このトピックに関するロシア語の出版物は非常に不十分でした。 「それは問題ではありません。外国人の同僚の作品を読みます」と私は思いました。 しかし、ここでは小さなキャッチが待っていました-この方法で利用可能なすべての資料と出版物は非常に理論的でした。 したがって、私と同じ不幸な人々（そして賢明な人からのアドバイスをある程度期待している人々）のために、この方法を簡単に説明し、それを実装する最も簡単な方法を提供することにしました。
オリジナルのPeskinメソッドのみを検討します。 私の意見では、ダミーセル（ゴーストセル）の方法、切り捨てられたセル（カットセル）の方法、およびその他すべての変更よりも簡単です。
したがって、ウィキペディアは、「没入境界法」は液体と流線形オブジェクト間の相互作用をモデル化するアプローチであることを示しています（生体系の細い繊維をモデル化するためによく使用されるため、「繊維」 ） そして、これはいまいましいです=）
このアプローチの主な特徴は、流体の流れを計算するため（オイラー座標で）と水没境界のパラメーターを計算するために2つの別個のグリッドが導入されることです（ラグランジュ座標で同じ繊維繊維）。 これにより、単純なグリッドと高速の計算方法を使用して流体を計算し、シミュレーションの複雑さ全体を没入境界との相互作用の段階に移すことができます。
純粋に理論的な観点から、このようなトリックは、Navier-Stokes方程式に液体上の浸漬境界の作用の力（より正確には力の密度）を追加することによって行われます。
この力は、水没境界の計算された点の位置に基づいて何らかの方法で計算され、流れを変更します。 完全な方程式系は次のようになります。
たとえば、次の形式の水没境界上の境界条件（接着条件、別名「滑り止め条件」）：
この連立方程式
-エリア内の流体の流れを記述する変数があります 
座標付き 
、そして 
-領域内の没入境界の状態を記述する変数 
座標付き 
。
実際、流体と境界の相互作用は最後の2つの方程式に隠れています。
実装の観点からこの不名誉をすべて考慮すると、各タイムステップのアルゴリズムはおよそ次のようになります。
- 没入境界の各点で、現在の状態に基づいて電圧の強度を計算します
- 得られた力を何らかのアルゴリズムに従って液体の点に分配します
- 受けた力を考慮して、流体の流れのパラメーターを計算します
- 流体のポイントから水没境界のノードまでの速度を補間します（境界条件により、境界の速度は流体の速度と一致する必要があります）
- 速度に応じて境界ノードの位置を変更します
ご覧のとおり、アルゴリズムには2つの不確実な点があります。電圧の強さを決定する方法、および流体の流れのパラメーターを計算する方法については何も言われていません。 これはこのアプローチの深刻なプラスです-これらの方法はほとんどすべての=になりますが、実装に近づくために、少なくとも境界の電圧を計算する方法について話し合う必要があります（2番目の方法は、それだけであれば、任意に選択できます）十分に正確でした）。
実際、タフガイは、エネルギー汎関数を使用して境界で形成される力を見つけ、引張力と曲げ力を考慮します。 しかし、残念なことに、私自身はまだこの瞬間を完全に理解していません。 したがって、カンニングして、精度は低くても計算には十分なペナルティメソッドを使用できます。 その本質は単純です-境界線が元の位置から離れるほど、出現する力は大きくなります。 その結果、偏差を計算し、特定の剛性係数を乗算するだけで、それを測定できます。
アルゴリズムの2番目の微妙な点は、連立方程式に現れるデルタ関数です。 数値版が必要です。良い意味で、これが曲全体です。 離散化のオプションは巨大なクラウドです。 ただし、この場合、これの使用を勧める同僚の経験に頼ることができます（3次元の場合も同様）。
これで、1000個の単語の代わりに、結果のアルゴリズムをコードの形で説明できます。
ノードでの電圧の計算は次のようになります。
for(int n = 0; n < boundary->nodes_count; ++n) { Node *node = &boundary->nodes[n]; // get_area , node->x_force = -boundary->stiffness * (node->x - node->x_ref) * boundary->get_area(); node->y_force = -boundary->stiffness * (node->y - node->y_ref) * boundary->get_area(); node->z_force = -boundary->stiffness * (node->z - node->z_ref) * boundary->get_area(); }
液体の点にかかる力の分布：
for(int n = 0; n < boundary->nodes_count; ++n) { // int x_int = index(boundary->nodes[n].x, COORD_X); int y_int = index(boundary->nodes[n].y, COORD_Y); int z_int = index(boundary->nodes[n].z, COORD_Z); // .. - , // , // ( ) for(int i = x_int-1; i <= x_int + 2; ++i) { for(int j = y_int-1; j <= y_int + 2; ++j) { for(int k = z_int-1; k <= z_int + 2; ++k) { // long double dist_x = fabs(boundary->nodes[n].x - coord(i, COORD_X)); long double dist_y = fabs(boundary->nodes[n].y - coord(j, COORD_Y)); long double dist_z = fabs(boundary->nodes[n].z - coord(k, COORD_Z)); long double weight_x = delta(dist_x, COORD_X); long double weight_y = delta(dist_y, COORD_Y); long double weight_z = delta(dist_z, COORD_Z); // force_X[i][j][k] += (boundary->nodes[n].x_force * weight_x * weight_y * weight_z) * boundary->get_area(); force_Y[i][j][k] += (boundary->nodes[n].y_force * weight_x * weight_y * weight_z) * boundary->get_area(); force_Z[i][j][k] += (boundary->nodes[n].z_force * weight_x * weight_y * weight_z) * boundary->get_area(); } } }
適切な方法で、フローを計算します。 次に、液体のノードに速度を補間します。
for(int n = 0; n < boundary->nodes_count; ++n) { // , int x_int = index(boundary->nodes[n].x, COORD_X); int y_int = index(boundary->nodes[n].y, COORD_Y); int z_int = index(boundary->nodes[n].z, COORD_Z); // for(int i = x_int-1; i <= x_int + 2; ++i) { for(int j = y_int-1; j <= y_int + 2; ++j) { for(int k = z_int-1; k <= z_int + 2; ++k) { long double dist_x = fabs(boundary->nodes[n].x - coord(i, COORD_X)); long double dist_y = fabs(boundary->nodes[n].y - coord(j, COORD_Y)); long double dist_z = fabs(boundary->nodes[n].z - coord(k, COORD_Z)); long double weight_x = delta(dist_x, COORD_X); long double weight_y = delta(dist_y, COORD_Y); long double weight_z = delta(dist_z, COORD_Z); boundary->nodes[n].x_vel += (velocity_U[i][j][k] * weight_x * weight_y * weight_z) * CB(Hx[i]); boundary->nodes[n].y_vel += (velocity_V[i][j][k] * weight_x * weight_y * weight_z) * CB(Hy[j]); boundary->nodes[n].z_vel += (velocity_W[i][j][k] * weight_x * weight_y * weight_z) * CB(Hz[k]); } } } }
浸された境界のノードの座標を更新する：
for(int n = 0; n < boundary->nodes_count; ++n) { boundary->nodes[n].x += boundary->nodes[n].x_vel; boundary->nodes[n].y += boundary->nodes[n].y_vel; boundary->nodes[n].z += boundary->nodes[n].z_vel; }
実際には、それだけです。 このアプローチでは、境界は柔軟であると想定されていますが、十分に高い剛性（実際には非常に大きく、実験で考慮する必要があります）を使用して、固体境界の動作をシミュレートできます。
現時点では、最小残差を不完全に近似する方法を使用して流体の流れを計算し、一連の定常計算によって初期の非定常問題を近似していることに注意してください。 もちろん、これは最良のオプションではありません。そうするべきではありませんが、それでも動作します。
すべて、私はこのテキストを慈悲に捧げます。 合理化に関するコメントと提案はすべて受け入れられます。 仕事はまだ理想からはほど遠い。
手始めに、十分に高い粘度で（およびかなり粗いグリッド上で）球の周りにいくつかの写真があります。
