不定方程求解c語言編寫自己的求解器

如果不寫自己的求解器，那麼學這個軟體將毫無意義。而寫自己的求解器，常見的，無外乎幾種目的：

這些需求裡面，其實除了第一條是需要改動求解器框子的，剩下的無非是在一般計算的中間插入一些計算步驟罷了。所以下面寫的時候我先從最簡單的求解器（簡單方程、沒有源項、常係數）開始寫起，然後單獨處理上述若干種不同情況。

ps 還有一種高階情況需要修改求解器，那就是不同的離散格式。這個我不太懂，感覺比較高階，以後再慢慢考慮吧。。 pps 這個筆記侷限在不可壓粘性流體裡面。可壓縮和尤拉方程我基本不懂。。。

所謂順序，就是乙個求解器從頭到尾要幹的事情。初學者對這個往往沒概念導致看求解器**的時候雲山霧罩的。這裡寫清楚一點：

#include "createmesh.h"

#include "createfields.h"

但實際上沒那麼easy，因為具體要讀取什麼變數和係數每個求解器不太一樣，有的時候所以要自己修改一些內容。

peqn.solve();

或者

fvvectormatrix ueqn(...);
solve(ueqn == -fvc::grad(p));

，但是of畢竟是個通用計算框架，裡面有大量的考慮鬆弛、非正交修正的語句。目前這個筆記只涉及最簡單的矩形結構化網格，多餘的因素盡量扣除。 + 最後就是輸出了，呼叫語句：

runtime.write();

即可。當然可以從乙個空檔案從第乙個字元開始敲，但是of提供了很多生成摸板的工具。19年初去學軟體的時候就有所耳聞，但是如何建立乙個空的求解器還是今天看了乙個教程才知道的，那就是使用命令：

在新建的空求解器中加入一句話：

#include "createmesh.h"

即可。這句話應該是將求解器路徑下的/constant/polymesh下的網格檔案讀入程式，程式設計乙個物件。

#include "createfields.h"

info<< "reading transportpropertiesn" << endl;
iodictionary transportproperties
( ioobject
("transportproperties",
runtime.constant(),
mesh,
ioobject::must_read_if_modified,
ioobject::no_write));

info<< "reading diffusivity dtn" << endl;
();

ok，接下來搞初場。我們的方程只有乙個變數，也就只有乙個初場。這個初場放在./0路徑下面，就叫做h。讀取的方式跟上面差不多。上面我們讀取係數的時候建立了乙個係數屬於的物件，那麼現在我們也要建立乙個變數屬於的物件。此處是h，是個標量所以我們要建立乙個標量場volscalarfield。在of裡，寫成：

info<< "reading field tn" << endl;
volscalarfield h
( ioobject
("t",
runtime.timename(),
mesh,
ioobject::must_read,
ioobject::auto_write
),mesh
);

關於為什麼這裡有兩個mesh，以及物件序號產生器的寫詳細說明，參閱資料

openfoam中的註冊機制是什麼意思？www.cfd-china.com

openfoam輸入輸出中的註冊機制 | 大專欄www.dazhuanlan.com 和

openfoam 中的物件註冊機制marinecfd.xyz

。從我個人的感覺來說，這部分屬於軟體架構和資料結構的內容，涉及的是如何把軟體寫的高效，並不涉及太多cfd的內容，如果僅僅是在of的框架下編寫自己的求解器或許不需要太過在意。。。。

接下來的任務是建立乙個離散的方程組。方程組在of中一般以fvvectormatrix或者fvscalarmatrix的類的物件存在。前者是用來儲存速度這類向量，壓力、溫度、vof相函式等標量儲存在後者。此處我們建立乙個fvscalarmatrix的物件heqn：

fvscalarmatrix heqn
( fvm::ddt(h)
);

可以看見這裡的寫法很類似與真實的數學寫法。常見的微分算符，其實就是瞬態微分，散度和laplace，分別寫作ddt,div,laplacian。這三者又分別屬於fvm和fvc兩個命名空間，fvm是隱式，fvc是顯式，取決與需要的離散格式。當然還有乙個要注意的地方是時間迴圈。在of裡沒有真正的穩態計算（？）所有的問題都是瞬態的，只是在穩態問題中可以放大時間步讓方程迅速收斂到穩態而已（此處存疑）。所以在寫穩態求解器的時候一樣要寫上瞬態項。既然有瞬態項就必須有步進策略，而在of中步進一般由******演算法的框架來完成，其實這種簡單求解器完全不需要******演算法，此處只是為了通用性而呼叫******而已。具體的語句是，在此前，就使用mesh生成乙個叫做******的物件，屬於******control類，語句是：

******control ******(mesh);

接下來，在上述離散步驟的外面需要套乙個迴圈：

while (******.loop(runtime))

構建好方程組之後，只需要乙個簡單的solve就能搞定求解：

heqn.solve();

runtime.write();

整個求解器就算初具雛形了。看一下整體的樣子：

#include "fvcfd.h"// 必備標頭檔案
#include "fvoptions.h"// 與源項有關，目前不涉及可以刪掉
#include "******control.h"// 時間推進有關
int main()// 主函式
runtime.write();//輸出
}info<< "endn" << endl;//提示計算結束
return 0;
}

接下來要做的是把源**程式設計乙個可執行程式。所有的solver都是可執行程式。所以需要編譯。編譯通過wmake即可，之前的make資料夾就是為這個準備的。

實際上示例的方程還是有點麻煩，因為

這樣子有2種處理方法：

fvscalarmatrix thetaeqn
(fvm::ddt(theta)
);h = -1*pow(-1600000+341532/(theta-0.075),0.252525);//代入求解h

沒想到竟然支援這麼狂野的表示式，竟然編譯通過了。。。等後面網格畫好跑一下試試。

wmake

就vans了。這部分有個很強的

可以參考。

fvoptions

不定方程求解c語言編寫自己的求解器

C語言方程求解

求解方程語言 c

c語言求解非線性方程組 C語言求解線性方程組

不定方程求解c語言 編寫自己的求解器

C語言 方程求解

求解方程 語言 c

c語言求解非線性方程組 C語言求解線性方程組

相關推薦

不定方程求解c語言編寫自己的求解器

C語言方程求解

求解方程語言 c