OpenFOAM的程序开发初步.docx
《OpenFOAM的程序开发初步.docx》由会员分享，可在线阅读，更多相关《OpenFOAM的程序开发初步.docx（16页珍藏版）》请在冰豆网上搜索。
OpenFOAM的程序开发初步
OpenFOAM的程序开发初步
一.OpenFOAM应用的类型：
使用OpenFOAM进行CAE模拟的，大致可分为三种类型：
1）直接利用OpenFOAM的标准的求解器进行模拟，把OpenFOAM替代商业软件来使用，OpenFOAM已基本具有这样的功能和人气，与Fuent，Star-CD等相比较，OpenFOAM显然具有更高的求解效率和灵活性。
2）用户自定义求解器，即利用OpenFOAM的基本类库，如finiteVolume，OpenFOAM库来按照自己的求解流程来编写针对某类应用的求解器。
用户需要开发的求解器就是类似于在OpenFOAM的applications中所看到的标准求解器icoFOAM，simpleFOAM等。
显然这一需求是非常大的，从OpenFOAM问世以来，已有很多用户定义了自己的求解器。
这类需求的特点是，并不需要特别关心，离散和求解的最底层的知识，如时间项离散，空间项离散等，关注的重点是求解的步骤或者流程。
在编程中，通常是顶层的求解流程的开发，在多数情况下可以不编译OpenFOAM的finiteVolume和OpenFOAM库。
这种顶层的求解器的开发，是我们以前常常忽略的，或者是以前没有能力做到的。
需要指出的是，商业软件中的所谓udf，usersubroutine和这是不可相比的。
3）用户自己定义离散方法等。
对于研究离散格式、代数求解器等人来说，更关注时间项ddt，扩散项Laplacian，对流项div是如何离散的，能否有更高效更高精度的离散方法，这需要修改finiteVolume库和OpenFOAM库中对应的代码。
尤其是对流项，尽管OpenFOAM已经提供了基于NVD和TVD的模板和40多种有名的高阶高精度格式，但可以预见，这仍然是不够的，毕竟对流项的离散仍然是目前CFD的重点研究方向。
可以肯定的是，目前有很多人关注类型2的应用，毕竟将OpenFOAM当成Fluent或Star-CCM来使用，并不见得方便。
但是将OpenFOAM作为类库来构建自己的求解器，这是其它软件无法实现的
OpenFOAM的程序开发初步
二．OpenFOAM程序开发的基本知识
2.1OpenFOAM的基本术语
重要的环境变量:
$WM_PROJECT_USER_DIR――OpenFOAM的用户目录
$FOAM_TUTORIALS------OpenFOAM的算例目录
$FOAM_SRC------OpenFOAM库的源程序目录
$FOAM_APP------OpenFOAM的求解器目录
$FOAM_APPBIN-------OpenFOAM的求解器执行文件目录
$FOAM_RUN------用户的算例目录
重要的shell：
run＝cdto$FOAM_RUN
src=cdto$FOAM_SRC
app=cdto$FOAM_APP
util=cdto$FOAM_APP/utilities
sol=cdto$FOAM_APP/solvers
tut=cdto$FOAM_TUTORIALS
求解器的基本文件结构
appName包含求解器源代码的目录
+appName.C求解器主程序
＋CreateFields.H场变量的声明和初始化
+Make/编译指令
+files编译需要的源程序文件和生成的目标文件
+options编译选项,如链接库等
appName/appName.C是求解器的主程序
appName/createFields.H声明变量，并从文件中读入初值，如p，物性。
appName/Make/files所有源程序的名称，一个文件一行，最后一行是目标代码的名称和存放位置，EXE=$（FOAM_USER_APPBIN）/appName
appName/Make/options设定查找头文件和库的路径，EXE_INCS,和需要链接的库EXE_LIBS
算例的基本文件结构
case/算例目录
+0/包含初始和边界条件
+constant/包含初次读入后，不随时间变化的数据
+polyMesh/包含多面体网格数据
+transportProperties/包含物性数据
+system/包含计算控制和离散格式设定
+controlDict包含计算控制，如时间步长等
+fvSchemes包含离散格式设定
+fvSolutions包含代数求解器或SIMPLE，PISO算法设定
具体而言
case/0每个需求解的变量需要一个文件设定其初始边界条件
case/constant/polyMesh网格数据，如ownerneighbourpointsfacesboundary
case/system/transportProperties物性数据
case/system/controlDict设定起始终止时间，时间步长，输出控制
case/system/fvSchemes设定程序用到的每个微分算子的离散格式
case/system/fvSolution为每个变量选择代数方程求解器/收敛精度及PISO等算法设定
三．OpenFOAM程序开发的理论知识
作求解开发，必须能写出需要求解的控制方程及其定解条件，并且对于如何求解方程或方程组的步骤已经明确。
这些流体力学、传热学以及相关的理论是必需的，所谓连续介质力学中的数学模型，控制方程和定解条件就是表示它的语言。
在这里是不可能说清楚的，这要看个人的功底了。
四.OpenFOAM程序开发的最简单的例子
下面采用OpenFOAM来开发一个用户自己的求解器。
主要是利用OpenFOAM的标准求解器icoFoam，用户不需要写任何代码，只为为了熟悉OpenFOAM程序开发的环境和步骤。
步骤：
1）将icoFoam目录拷贝到新的目录
可采用下面的Linux的命令实现：
到OpenFOAM的incompressible目录
cdapplications/incompressible
cp–ricoFoammyicoFoam
以上只是复制目录icoFoam到新的位置，并且新目录名为myicoFoam
cdmyicoFoam
进入新的目录，查看一下，可以看到里面的文件和icoFoam中是否一样
2）原文件改名，并且删除依赖文件
将icoFoam.C改名myicoFoam.C
mvicoFoam.CmyicoFoam.C
删除依赖文件
rmicoFoam.dep
3）修改编译文件files和options
进入Make目录，打开files文件
将
icoFoam.C源程序文件名
EXE=$（FOAM_APPBIN）/icoFoam可执行文件名
修改为
myicoFoam.C源程序文件名
EXE=$（FOAM_APPBIN）/myicoFoam可执行文件名
此例中options不需修改，可以打开看看
EXE_INC=\头文件包含
-I$（LIB_SRC）/finiteVolume/lnInclude
EXE_LIBS=\链接库
-lfiniteVolume
4）删除原来的obj文件
rm–rflinuxGccDPOpt
cd..
5）编译
wmake
6）检验一下
到tutorial目录，检验一下
myicoFoam.cavity
六.OpenFOAM程序开发――例子一：
在icoFoam中加入温度场求解
准备：
能量控制方程：
dT/dt+div（den*U*T）=div（agradT）
在壁面上给定值条件。
需要解决的问题：
a）如何创建标量场，T
b）如何创建物性，a
c）如何定义温度方程，并求解
d）如何在算例中设定T和a
e）如何设定T的离散格式
f）如何设定T的求解器的收敛标准等
步骤：
1）创建程序需要的新物性和新变量场
打开myicoFoam.C可以看到，程序开始运行时调用CreateFields.H，创建变量场。
打开CreateFields.H，可以看到程序首先从transportProperties文件中读入物性，
Info<<"ReadingtransportProperties\n"<IOdictionarytransportProperties
（
IOobject
（
"transportProperties",从字典文件transportProperties读入
runTime.constant（）,//transportProperties文件位于目录runTime.constant（）中
mesh,网格对象
IOobject:
:
MUST_READ,
IOobject:
:
NO_WRITE
）
）;创建了Iodictionary类型对象transportProperties
dimensionedScalarnu//首先读入粘性系数
（
transportProperties.lookup（"nu"）
）;创建有量纲标量nu，nu通过从字典transportProperties查找”nu”来赋值
可以加上新方程需要的物性
dimensionedScalarDT//首先读入热扩散率
（
transportProperties.lookup（"DT"）
）;创建有量纲标量DT，DT通过从字典transportProperties查找”DT”来赋值
此外还要从createFields中读入p,U场，我们要加入的新的变量场为温度场T，最快的加入温度场的方法是拷贝p场的代码，修改为
Info<<"ReadingfieldT\n"<volScalarFieldT
（
IOobject
（
"T",
runTime.timeName（）,
mesh,
IOobject:
:
MUST_READ,
IOobject:
:
AUTO_WRITE
）,
mesh
）;
这样，创建了新的vol标量场T，从文件T中读入。
对于T的创建具体解释如下：
a）创建了标量场T
b）T通过读（IOobject:
:
MUST_READ）在runTime.timeName（）目录下名称为“T”的文件创建，在开始计算时，runTime.timeName（）是contorlDict中设定的startTime值决定的。
c）T将自动写入（IOobject:
:
AUTO_WRITE）计算结果到runTime.timeName（）目录中，runTime.timeName（）随迭代是变化的，写入控制由contorlDict中设定。
d）T是定义在mesh对象上的，这意味着T在内部cell上有值internalField，在边界上还需要边界条件，这与polyMesh/boundary中要一致。
2）在求解器中加入新的求解方程
下一步回到myicoFoam.C加入新的微分方程，由于温度场依赖于速度场，可放在PISO循环后面。
#include"continuityErrs.H"
U-=rUA*fvc:
:
grad（p）;
U.correctBoundaryConditions（）;
//Addthetemperatureequation
fvScalarMatrixTeqn温度是标量方程
（
fvm:
:
ddt（T）
+fvm:
:
div（phi,T）要用到界面流量
-fvm:
:
laplacian（DT,T）扩散项
）;
TEqn.solve（）;求解
3）编译
wmake
4）在算例中加入新方程的初始和边界条件
4.1拷贝一个cavity算例到mycavity
4.2修改transportProperties字典文件，设定DT
cdconstant
修改transportProperties文件，前面已提到DT要从该字典文件读入。
设定DT=0.002m2/s
DTDT[02-10000]0.002;
4.3修改T文件，设定初始值和边界
cd0进入0目录
拷贝一个T文件
cppT
修改T文件为
FoamFile
{
version2.0;
formatascii;
classvolScalarField;
objectT;
}
//*************************************//
dimensions[0001000];
internalFielduniform300;初始内部点为300℃
movingWall
{
typefixedValue;
valueuniform350.;边界为350℃
}
fixedWalls
{
typefixedValue;
valueuniform300.;边界为300℃
5）修改离散格式和代数求解器求解控制文件
A进入system目录
由于温度方程有非稳态项，对流项，扩散项，分别要在ddt，div，laplacian中设置
打开fvSchemes文件，添加
divSchemes
{
defaultnone;
div（phi,U）Gaussupwind;
div（phi,T）Gaussupwind;
}
laplacianSchemes
{
defaultnone;
laplacian（nu,U）Gausslinearcorrected;
laplacian（DT,T）Gausslinearcorrected;
laplacian（（1|A（U））,p）Gausslinearcorrected;
}
在fvSolution中设置代数求解器选项
TPBiCG
{
preconditionerDILU;
tolerance1e-06;
relTol0;
};
注意T方程形成的矩阵是非对称的，不要用PCG和DIC
6）运行
myicoFoam.mycavity
七.OpenFOAM程序开发――求解器的详细分析1
进入icoFoam目录
可以看到
createFields.HicoFoam.CicoFoam.depMake/
Make/为wmake编译所需的文件
IcoFoam.C为主程序文件，它包含createFields.H
编辑icoFoam.C
可以看到icoFoam.C首先引入的头文件为fvCFD.H。
所以你可以看到，在编译选项options中
EXE_INC=\
-I$（LIB_SRC）/finiteVolume/lnInclude//fvCFD.H的存放目录
EXE_LIBS=\
-lfiniteVolume//需要链接的库
找到fvCFD.H,编辑，可以看出这些是主程序必须的类库
#ifndeffvCFD_H
#definefvCFD_H
#include"parRun.H"
#include"Time.H"时间类
#include"fvMesh.H"网格类
#include"fvc.H"fvc类
#include"fvMatrices.H"fvMatrix类
#include"fvm.H"fvm类
#include"linear.H"
#include"calculatedFvPatchFields.H"
#include"fixedValueFvPatchFields.H"
#include"adjustPhi.H"
#include"findRefCell.H"
#include"mathematicalConstants.H"
#include"OSspecific.H"
#include"argList.H"
#ifndefnamespaceFoam
#definenamespaceFoam
usingnamespaceFoam;
#endif
#endif
再看看icoFoam的程序体,了解一下求解程序的结构
#include"fvCFD.H"――――――――――――――――（头文件）
通常位于main函数前，是程序所需的类的定义
//*************************************//
intmain（intargc,char*argv[]）
{
―――――――――――――――（包含文件）
#include"setRootCase.H"
#include"createTime.H"
#include"createMesh.H"
―――――――――――包含文件通常是程序片断，如创建时间、创建网格等
―――――――――――――――（求解器代码）――――――
#include"createFields.H"
需要根据应用，单独写的代码，如"createFields.H"和Main，以及Ueqn，pEqn等
―――――――――――――――――――――――――――――――――――――
#include"initContinuityErrs.H"
//*************************************//
。
。
。
。
。
}
八.OpenFOAM程序开发――求解器的详细分析2
a.场变量的定义
引用前面的温度场
Info<<"ReadingfieldT\n"<volScalarFieldT
（
IOobject
（
"T",
runTime.timeName（）,
mesh,
IOobject:
:
MUST_READ,
IOobject:
:
AUTO_WRITE
）,
mesh
）;
例如
volScalarFieldCO2
（
IOobject
（
"CO2",
runTime.timeName（）,
mesh,
IOobject:
:
READ_IF_PRESENT,
IOobject:
:
AUTO_WRITE
）,
mesh,
//Optionaldeclaration,thiscanbedonebyaccessingafilein"case/0/"，量纲可在文件中读
//dimensionedScalar（"zero",dimensionSet（1,-1,-3,0,0,0,0）,value）
）;
b.场变量的使用
场变量有定义在内部cell上的值，也有边界上的值
例如给组分限值
Exampleofamassfractionlimiterusedinthisproject:
//InitializethevariableY_iforuseinaloop
scalarField&CO2Internal=CO2.internalField（）;引用内部点
//Loopforallmeshpoints遍历内部点
forAll（CO2,celli）
{
//Limitsthemassfractiontoapositivenumber
if（CO2Internal[celli]<0.0）
{
CO2Internal[celli]=0.0;
}
//Limitsthemassfractiontomax1.0
if（CO2Internal[celli]>1.0）
{
CO2Internal[celli]=1.0;
}
}
c.定义输运方程
OpenFOAM定义方程时要选择一种类型的fvMatrix，有fvScalarMatrix和fvVectorMatrix
离散格式在case/system/fvSchemes.中设定
//DefineaScalarMatrixasaobject
fvScalarMatrixCO2Eqn定义系数矩阵
（
fvm:
:
div（phi,CO2）对流项离散fvm:
:
div
-fvm:
:
laplacian（turbulence->nuEff（）,CO2）扩散项离散fvm:
:
div
==S_CO2源项
）;
//Applyunderrelaxationtotheequation
//Underrelaxationfactorsdefinedinfile:
fvSolution
CO2Eqn.relax（）;松弛
CO2Eqn.solve（）;求解
升级会员 