【OpenFOAMとOpenMoelica連携（３）】chtMultiRegionFoamによるマルチリージョンソルバ
モデルの作成
モデル形状はFreeCADで行いました。
流体側の領域(fluid）と固体側の領域（solid）を別々で作成してstlファイルとして出力します。
出力したstlファイルを元にメッシュ作成を行いました。
主に使用したコマンドを載せておきます。
スケール変換 surfaceConvert -scale 0.001 model/solid/solidwalls.stl constant/triSurface/solidwalls_m.stl surfaceConvert -scale 0.001 model/fluid/fluid.stl constant/triSurface/fluid_m.stl 特徴量抽出 surfaceFeatureExtract -dict system/surfaceFeatureExtractDict.solid surfaceFeatureExtract -dict system/surfaceFeatureExtractDict.fluid メッシュ作成 blockMesh snappyHexMesh -overwrite 領域分割 splitMeshRegions -cellZones -overwrite
chtMultiRegionのマルチリージョンソルバは作り方なども癖があるので、なかなか理解ができないかもしれません。
温度の境界条件としては以下のように室内は固体壁で囲まれており、その固体壁のまわりのモデル化は行っていないため、固体壁と外気は熱伝達率$10[\text{W}/\text{m}^2 \text{K}]$としています。
計算の流れ
計算の流れは前回の記事と同じです。
【OpenFOAMとOpenMoelica連携（２）】buoyantPimpleFoamソルバ
こんにちは（@t_kun_kamakiri）
前回の記事ではOpenModelicaのコマンド実行の方法について紹介しました。
...
ただし今回はマルチリージョンソルバなのでフォルダ構成が少し異なります。
0, constant, systemのそれぞれのフォルダにfluid領域とsolid領域を作成する必要があります。fluid領域とかsolid領域などの名前は任意でつけることができます。
好みに応じてAir（流体領域）やwalls（固体領域）などにしても良いです。
計算毎時に実行する独自プログラム
system/controlDictは計算毎時に独自のプログラムを実行するcodedという機能があります。
設定ファイルの全文を載せておきます。
system/controlDict
application chtMultiRegionFoam;
startFrom latestTime;
startTime 0.001;
stopAt endTime;
endTime 300;//21600;// 6 hours
deltaT 0.001;
writeControl adjustable;
writeInterval 5;//600;// 10 minutes
purgeWrite 0;
writeFormat ascii;//binary;//
writePrecision 12;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
maxCo 20.0;//0.9;
maxDi 10.0;
adjustTimeStep yes;
functions
{
#include "outletTempAverageValue"
continuityError1
{
type continuityError;
libs
(
fieldFunctionObjects
); // Mandatory entries (unmodifiable)// Optional entries (runtime modifiable)
phi phi; // Optional (inherited) entries
writePrecision 8;
writeToFile yes;
useUserTime yes;
region region0;
enabled yes;
log yes;
timeStart 0;
// timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeControl timeStep;
writeInterval 1;
}
residuals
{
type solverInfo;
libs
(
"libutilityFunctionObjects.so"
);
fields
(
".*"
);
}
fieldMinMax1
{
type fieldMinMax;
libs
(
fieldFunctionObjects
); // Mandatory entries (unmodifiable)// Mandatory entries (runtime modifiable)
mode magnitude;
fields
(
U
p
T
); // Optional entries (runtime modifiable)
location yes;
writePrecision 8;
writeToFile yes;
useUserTime yes;
region fluid;
enabled yes;
log yes;
timeStart 0;
// timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeControl timeStep;
writeInterval 1;
}
fieldMinMax2
{
type fieldMinMax;
libs
(
fieldFunctionObjects
); // Mandatory entries (unmodifiable)// Mandatory entries (runtime modifiable)
mode magnitude;
fields
(
U
p
T
); // Optional entries (runtime modifiable)
location yes;
writePrecision 8;
writeToFile yes;
useUserTime yes;
region solid;
enabled yes;
log yes;
timeStart 0;
// timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeControl timeStep;
writeInterval 1;
}
yPlus1
{
type yPlus;
libs
(
fieldFunctionObjects
); // Mandatory entries (unmodifiable)// Optional (inherited) entries
writePrecision 8;
writeToFile yes;
useUserTime yes;
region region0;
enabled yes;
log yes;
timeStart 0;
timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeControl timeStep;
writeInterval 100;
}
YMin
{
type surfaceFieldValue;
libs
(
"libfieldFunctionObjects.so"
);
writeControl timeStep;
log yes;
writeFields no;
regionType patch;
name ZMax;
operation sum;
fields
(
phi
); // Optional (inherited) entries
writePrecision 8;
writeToFile yes;
useUserTime yes;
region air;
enabled yes;
timeStart 0;
timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeInterval 100;
}
YMax
{
$YMin;
name YMax;
}
//https://doc.openfoam.com/2306/tools/post-processing/function-objects/field/wallHeatFlux/
wallHeatFlux1
{
// Mandatory entries (unmodifiable)
type wallHeatFlux;
libs (fieldFunctionObjects);
// Optional entries (runtime modifiable)
//patches (<patch1> ... <patchN>); // (wall1 "(wall2|wall3)");
//qr qr;
// Optional (inherited) entries
// writePrecision 8;
// writeToFile true;
// useUserTime true;
region solid;
// enabled true;
// log true;
// timeStart 0;
// timeEnd 8000;
executeControl timeStep;
executeInterval 1;
writeControl writeTime;//timeStep;
writeInterval 1;
}
大事なのは#include "outletTempAverageValue";の部分です。
ここに独自のプログラムを書いておき、includeで読み込んでいます（直接記述すると煩雑になるためファイルを分けています）
特定の境界面の温度、特定の場所の温度を取得
outletTempAverageValueの全文を載せておきます。
system/outletTempAverageValue
outletTempAverageValue1
{
type coded;
name codeTempFieldValue;
libs ("libutilityFunctionObjects.so");
// 重要: どのリージョンで動かすか
region fluid;
// 重要: いつ実行するか
enabled yes;
executeControl timeStep;
executeInterval 1;
codeExecute
#{
Info << "===== codeTempFieldValue (coded) start =====" << nl;
// ---- ここから元の処理 ----
const volVectorField& C = mesh().C();
const scalarField& V = mesh().V();
word outletBoundaryName = "f-outlet";
const volVectorField& U = mesh().lookupObject<volVectorField>("U");
const volScalarField& T = mesh().lookupObject<volScalarField>("T");
// targetNameのIDの取得
const label outletID = mesh().boundaryMesh().findPatchID(outletBoundaryName);
// targetNameのIDから速度場の取得
const fvPatchVectorField& outletUList = U.boundaryField()[outletID];
// targetNameのIDから温度場の取得
const fvPatchScalarField& outletTList = T.boundaryField()[outletID];
// targetNameのIDから面積ベクトルの取得
const scalarField& outletmagSf = mesh().magSf().boundaryField()[outletID];
// プローブ点の温度を取得
scalar Tprobe = GREAT; // 初期化（並列用）
scalar xs = 1.0;
scalar ys = 2.5;
scalar zs = 2.5;
Info << " === probe Temperature === " << endl;
// probeT = T.internalField()[mesh().findCell(C.findNearest(point(xs, ys, zs)))];
Tprobe = T[ mesh().findCell(point(xs, ys, zs)) ];
Info << mesh().findCell(point(xs, ys, zs)) << " Tprobe = " << Tprobe << endl;
Info << " === probe Temperature end === " << endl;
// ================================== Start WeightAreaAverageTemp==================================
// Open the file to write the results
std::ofstream outletWeightAreaTempFile;
outletWeightAreaTempFile.open("constant/outletWeightAreaTempFile");
if (!outletWeightAreaTempFile.is_open())
{
FatalErrorIn("codeExecute") << "Cannot open outletWeightAreaTempFile file for writing!" << exit(FatalError);
}
scalar Toutlet;
// 流出温度の計算
Toutlet = (gSum(mag(outletUList)*outletmagSf*outletTList) )
/ (gSum(mag(outletUList)*outletmagSf) );
// Write the header part of the file
outletWeightAreaTempFile << "Tout," << Toutlet << ";" << std::endl;
outletWeightAreaTempFile << "Tprobe," << Tprobe << ";" << std::endl;
outletWeightAreaTempFile << "Time," << mesh().time().timeName() << ";" << std::endl;
// Close the file
outletWeightAreaTempFile.close();
// ================================== End WeightAreaAverageTemp==================================
// ================================== Python同期処理 Start ==================================
Info << "===== Execute Python updateMOS_and_Run.py =====" << nl;
// ロックファイル作成
const fileName lockFile("constant/updateMOS.lock");
{
std::ofstream lock(lockFile.c_str());
if (lock.is_open())
{
lock << "LOCKED" << std::endl;
lock.close();
Info << "Lock file created: " << lockFile << endl;
}
else
{
FatalErrorIn("codeExecute") << "Cannot create lock file!" << exit(FatalError);
}
}
// Pythonスクリプトの実行（同期モードでも非同期でもOK、ここでは同期で実行）
const fileName pyScript("updateMOS_and_Run.py");
int ret = std::system(("python3 " + pyScript).c_str());
if (ret != 0)
{
Info << "Warning: Python script exited with nonzero status: " << ret << endl;
}
// ロックファイルが削除されるまで待機
while (Foam::isFile(lockFile))
{
Info << "Waiting for Python script to finish..." << endl;
Foam::sleep(0.01); // 0.01秒スリープ
}
Info << "Python script finished. Lock removed." << endl;
// ================================== Python同期処理 End ================================
#};
}
いろいろ書いていますが、やっているのは以下の事です。
- 流量で重みづけしたoutlet境界面の温度$\overline{T} = \frac{\sum_i \rho v_iA_i C_p T_i}{\sum_i \rho v_iA_i C_p }= \frac{\sum_i \rho Q_i C_p T_i}{\sum_i \rho Q_i C_p }$を求める。
- 座標$(1, 2.5, 2.5)$の温度を所得
- それを
constant/outletWeightAreaTempFileファイルに記述 - updateMOS_and_Run.pyのPythonファイルを実行
ただし、Pythonの実行が完全に終わるまで待機させる必要がある。事前にconstant/updateMOS.lockを作成しておき、Pythonファイルの実行が終わるとconstant/updateMOS.lockが削除される（この記述はupdateMOS_and_Run.pyに書かれている）ことでPythonの実行終了のフラグとしている。
途中作成されたconstant/outletWeightAreaTempFileは以下のような記述になる。
constant/outletWeightAreaTempFile
Tout,293.153; Tprobe,297.223; Time,1.21601;
このToutとTprobとTimeはOpenModelicaの実行時に必要な変数です。
Follow @t_kun_kamakiri
