C++

【オブジェクト指向C++】1次元拡散方程式の数値計算

こんにちは(@t_kun_kamakiri

本記事ではC++を使った1次元拡散方程式を題材にオブジェクト指向でのコーディングを解説します。
コードの解説は以下の記事に詳しく書いています。

本記事の内容

1次元の拡散方程式をC++で数値計算

1次元拡散方程式
\begin{align*}\frac{\partial u(x,t)}{\partial t}=\nu\frac{\partial^2 u(x,t)}{\partial x^2}\tag{1}\end{align*}
数値計算にはC++を使います。

カマキリ

C++の勉強中です

C++は3か月前に勉強を始めたばかりなのでコーディングスタイルの良し悪しに関してはこれからもっと磨いていかないといけないと思っています。今回はアウトプットの場として勉強したことをメモとして残しておきます。

【環境】
Microsoft Visual C++ (Visual Stido 2022)

解く方程式

1次元拡散方程式
\begin{align*}\frac{\partial u(x,t)}{\partial t}=\nu\frac{\partial^2 u(x,t)}{\partial x^2}\tag{1}\end{align*}
以下のように離散化します。

コードは前回の記事のを利用します。
コードの解説は以下の記事に詳しく書いています。

プログラム保存先

https://github.com/kamakiri1225/1DNS

コードの修正

まずは2階微分の中心差分による離散化$\frac{\partial^2 u}{\partial x^2}=\frac{u_{i+1}-2u_{i}+u_{i-1}}{dx^2}$をDiff1d.hに定義します。

Diff1d.h

$dx^2$を計算するために冒頭で「#include <cmath>」をインクルードし、「pow((vec[i].x1d – vec[i – 1].x1d), 2.0);」で2乗を計算しています。

Solution1d.hのSolution1d関数の引数3つにします。

同様にSolution1D.cppのSolution1d関数の引数も3つにします。

次にmain.cppの、

を以下に変更します。

先ほどSolution1d.hとSolution1d.cppでSolution1d関数の引数を3つにしたので、main.cppでも以下のように粘性係数viscを定義し、「Solution1d sol1d(ntimestep, dt, visc);」でインスタンス化するときの引数を3つにします。

拡散方程式による初期状態からの拡散を見るためにタイムステップ数を大きくします。

出力するファイル名を以下に変更します。

修正は以上です。

全体のコード

前回の記事から変更したのはmain.cppとDiff1d.hだけです。

main.cpp

creatFields.hを作ってUの設定は別のヘッダーファイルで定義するようにしました。

creatFields.h

さらにプログラム作成途中でU[i].vauleの値を確認したい場合に出力するプログラムをヘッダーファイルに書きました。

main.cppで

によりインクルードしています。

PrintVector.h

main.cppで

の2行でよりでU[i].vauleの値を確認できます。

Mesh.h

Mesh.cpp

Fields.h

Fields.cpp

Diff1d.h

Solution1D.h

Solution1D.cpp

FileWriter.h

FileWriter.cpp

出力結果の確認

計算を実行し、計算が終了すると以下の2つのファイルができています。

  • 1ddiffusion_initial.dat:初期状態
  • 1ddiffusion_final.dat:最終状態

今回はさくっと確認するためExcelで重ねてみました。

  • オレンジ:初期状態
  • 薄青:最終状態(10000ステップ目)

参考書

1週間でC++の基礎が学べる本

1週間でC++の基礎が学べる本

亀田 健司
1,430円(03/28 11:00時点)
発売日: 2021/03/12
Amazonの情報を掲載しています
独習C++ 新版

独習C++ 新版

高橋 航平
3,346円(03/28 11:00時点)
発売日: 2019/11/11
Amazonの情報を掲載しています
江添亮のC++入門

江添亮のC++入門

江添 亮
5,500円(03/28 19:52時点)
Amazonの情報を掲載しています

COMMENT

メールアドレスが公開されることはありません。 が付いている欄は必須項目です