ofsp2025-06
参考:
- https://github.com/UnnamedMoose/BasicOpenFOAMProgrammingTutorials
- https://www.topcfd.cn/simulation/solve/openfoam/openfoam-program/
- https://www.tfd.chalmers.se/~hani/kurser/OS_CFD/
- https://github.com/ParticulateFlow/OSCCAR-doc/blob/master/openFoamUserManual_PFM.pdf
- https://www.youtube.com/watch?v=KB9HhggUi_E&ab_channel=UCLOpenFOAMWorkshop
- http://dyfluid.com/#
- https://ss1.xrea.com/penguinitis.g1.xrea.com/study/OpenFOAM/index.html
感谢原作者们的无私引路和宝贵工作。
OpenFOAM 中,当某一些类实现某个特定的功能,可以把它们开发成库。
对“类”“库”链接有疑问的,请参考 OpenFOAM开发编程基础00 基本实现和开发 | 𝓐𝓮𝓻𝓸𝓼𝓪𝓷𝓭 (aerosand.cn)
我们将上篇讨论中的压力场、速度场的伪计算开发成独立的库,然后在应用中调用它,当然开发库也可以被任何应用调用。
建立本文的项目文件夹并进入
// terminal
ofsp
mkdir 06_devLib
cd 06_devLib
应用准备
// terminal
foamNewApp 06_01_devLib
cd 06_01_devLib
cp -r $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity debug_case
cp -r debug_case/0 debug_case/0_orig
在应用文件夹中新建开发库
// terminal
foamNewApp computeVelocityPressure
脚本和说明
略。
文件结构
整个应用的文件结构如下
|- 06_01_devLib/
|- computeVelocityPressure/
|- Make/
|- files
|- options
|- computeVelocityPressure.H
|- computeVelocityPressure.C
|- debug_case/
|- 0/
|- constant/
|- system/
|- Make/
|- files
|- options
|- createFields.H
|- 06_01_devLib.C
|- caserun
|- caseclean
|- README.md
开发库
头文件
头文件 /userApp/computeVelocityPressure/computeVelocityPressure.H
,内容如下
#include "fvCFD.H"
scalar computeR(const fvMesh& mesh, volScalarField& r, dimensionedVector x0);
volScalarField computePressure(const fvMesh& mesh, scalar t, dimensionedVector x0, scalar scale,scalar f);
void computeVelocity(const fvMesh& mesh, volVectorField& U, word pName = "p");
// 第三个参数给了默认值
- 各个函数的原型
- 为什么开发库的头文件要包含
fvCFD.H
?读者可以尝试注释掉这行,重新编译,并阅读报错信息。这样可以加深对fvCFD.H
文件内容的理解。
源代码
源文件 /userApp/computeVelocityPressure/computeVelocityPressure.C
,我们采用更“OpenFOAM”风格的实现方式,内容如下
注意和 OpenFOAM开发编程基础05 场的基本操作 | 𝓐𝓮𝓻𝓸𝓼𝓪𝓷𝓭 (aerosand.cn) 比较各个函数的实现
#include "computeVelocityPressure.H"
scalar computeR(const fvMesh& mesh, volScalarField& r, dimensionedVector x0)
{
r = mag(mesh.C() - x0); // 场 r 保存各个单元网格中心点向量到参考向量的距离
return max(r).value(); // 返回各个距离的最大值
}
volScalarField computePressure // 计算压力场,注意返回类型
(
const fvMesh& mesh, // 传入 mesh 的引用,引用传递效率高
scalar t,
dimensionedVector x0, // 不仅是 vector,而且是有单位制的 vector
scalar scale, // 纯 scalar
scalar f // 纯 scalar
)
{
volScalarField r(mag(mesh.C()-x0)/scale);
// 上一篇讨论中的 r 只是 scalar 类型,这里是 volScalarField 类型
// 注意计算中变量类型的统一
volScalarField rR(1.0/(r+dimensionedScalar("tmp",dimLength,1e-12)));
// 赋予 1e-12 单位,以通过OpenFOAM对计算的单位检查
// 注意计算中结果单位的统一
return Foam::sin(2.0*Foam::constant::mathematical::pi*f*t)*rR*dimensionedScalar("tmp",pow(dimLength,3)/pow(dimTime,2),1.0);
// 单位可以按数学方式计算
// 同样为了保证结果单位一致,额外乘以含单位临时量 tmp
// 该函数返回的计算结果的单位,等于,主函数要接收该值的变量的单位
}
void computeVelocity(const fvMesh& mesh, volVectorField& U, word pName)
// 注意函数的实现不能再写上形参的默认值
{
const volScalarField& pField = mesh.lookupObject<volScalarField>(pName);
// 按 pName 字段查找场量,并赋值给 pField
U = fvc::grad(pField) * dimensionedScalar("tmp",dimTime,1.0);
// 计算 U,保证单位统一
}
也许会有同学对比数据类型有困惑,比如
scalar
类型和volScalarField
类型。简单理解就是,scalar
只是一个值,而volScalarField
是一组数据,包含了计算域内所有离散点的scalar
值,相当于是一个矩阵,也就是表示了一个“场”。
库 Make
/userApp/computeVelocityPressure/Make/files
内容如下
computeVelocityPressure.C
LIB = $(FOAM_USER_LIBBIN)/libcomputeVelocityPressure
- 注意
LIB
关键字的变化
/userApp/computeVelocityPressure/Make/options
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools
- 该库,强调是这个库,在编译的时候,不会使用到其他更多的库,所以
options
中有这些基础库就够了,不需要修改
编译
终端使用命令
// terminal
wmake computerVelocityPressure/
// wclean computerVelocityPressure/
- 这个库也是独立的,可以被其他任何应用调用,只要写对调用路径
createFields.H
创建场 /userApp/createFields.H
,内容上和上一篇讨论的场文件相同
Info<< "Reading transportProperties\n" << endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
dimensionedScalar nu("nu",dimViscosity,transportProperties);
Info<< "Reading field p\n" << endl;
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field T\n" << endl;
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField zeroScalarField
(
IOobject
(
"zeroScalarField",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("zeroScalarField",dimless,Zero)
);
volVectorField zeroVectorField
(
IOobject
(
"zeroVectorField",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedVector("zeroVectorField",dimless,Zero)
);
volVectorField gradT
(
IOobject
(
"gradT",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::grad(T)
);
Info<< "Reading/calculating face flux field phi\n" << endl;
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
fvc::flux(U)
// fvc::interpolate(U)*mesh.Sf();
);
主源码
主源码 /userApp/06_01_devLib.C
,内容如下
#include "fvCFD.H"
#include "computeVelocityPressure.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
dimensionedVector originVector("x0",dimLength,vector(0.05,0.05,0.005));
// 构造 dimensionedVector 类的对象 originVector,并赋单位和初始值
scalar f(1.0);
volScalarField r // 只是计算的过程量,所以没有放进 createFields.H 中
(
IOobject
(
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("r",dimLength,0.0)
);
const scalar rFarCell = computeR(mesh,r,originVector);
// 计算得到距离的最大值,并传给 rFarCell 变量保存
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop()) // 时间推进
{
Info<< "Time = " << runTime.timeName() << nl << endl;
p = computePressure(mesh,runTime.time().value(),originVector,rFarCell,f);
// 计算 p 场
computeVelocity(mesh,U);
// 计算 U 场
phi = fvc::flux(U);
// 计算 phi 场
runTime.write(); // 写出场
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< nl;
runTime.printExecutionTime(Info);
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //
应用 Make
/userApp/Make/files
06_01_devLib.C
EXE = $(FOAM_USER_APPBIN)/06_01_devLib
/userApp/Make/options
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-IcomputeVelocityPressure/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-L$(FOAM_USER_LIBBIN) -lcomputeVelocityPressure
- 注意开发库的处理,既要“包含”,也要“链接”。
编译运行
wmake
./caserun
运行结果如下
Create time
Create mesh for time = 0
Reading transportProperties
Reading field p
Reading field T
Reading field U
Reading/calculating face flux field phi
Starting time loop
...
ExecutionTime = 0.01 s ClockTime = 0 s
End
我们也可以通过 paraview
打开查看计算结果。
// terminal
paraFoam -case debug_case
小结
研究中涉及到的通用方法可以抽象成开发库,但是初学者不能追求把方法都写成开发库,具体问题应当具体对待。通过本篇讨论,最主要的还是熟悉代码架构和代码语法。
此外,编程中注意单位统一,注意数据类型统一。