1 准备文件

本演示案例只需要一个源文件。利用下面的命令创建案例文件结构。

cd  $FOAM_RUN
mkdir demo4 && cd demo4
mkdir Make
touch demo4.C Make/files Make/options
cp -r $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity .

文件结构如下图所示。

这里的cavity文件夹用于程序测试。

修改files文件与options文件。

• files文件的内容

demo4.C
EXE = demo4

• options文件的内容

EXE_INC = 
        -I$(FOAM_SRC)/finiteVolume/lnInclude 
        -I$(FOAM_SRC)/meshTools/lnInclude
 
EXE_LIBS = 
        -lfiniteVolume 
        -lmeshTools

2 源文件

本案例演示从文件中获取场数据，然后再对读取的数据进行操作。

程序代码如下。

#include "fvCFD.H"
 
// 声明一个函数，后面有该函数的具体实现过程
// 输入时间t，空间坐标x，参考点x0以及缩放因子scale
scalar calculatedPressure(scalar t, vector x, vector x0, scalar scale);
 
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
 
    Info << "读取transportProperties文件数据" << endl;
    // 从transportProperties文件中读取数据
    // 先定义一个IOdictionary对象，其构造函数参数为一个IOobject对象
    IOdictionary transportProperties(
        IOobject(
            "transportProperties",           // 字典文件名
            runTime.constant(),              // 字典文件所在路径，这里为constant文件夹下
            mesh,                            // 一个objectRegistry类对象，这里没什么用
            IOobject::MUST_READ_IF_MODIFIED, // 如果文件被修改，则必须重新读取
            IOobject::NO_WRITE               //表示文件为只读
            ));
 
    // 定义一个dimensionedScalar变量nu
    // nu有量纲，其量纲dimViscosity等同于(0,2,-1,0,0,0,0),单位为m2/s
    dimensionedScalar nu(
        "nu",         //指定名称
        dimViscosity, //指定scalar的量纲
        // 旧版本或org版本的写法：transportProperties.lookup("nu")
        transportProperties // com新版本写法，自动根据名称在字典中搜索
    );
    Info << "读取到的nu为：" << nu << endl;
 
    //-------------从p文件中读取压力场p--------
    // 不需要查找关键字，因为整个文件都是关于压力场的数据
    Info << "读取压力场数据" << endl;
    // 定义一个标量场p，无需指定量纲，因为其量纲已经在相应的文件中指定了
    volScalarField p(
        IOobject(
            "p",                //指定名称
            runTime.timeName(), // 获取当前时间
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE),
        mesh);
 
    // -------从U文件中读取速度场---------
    Info << "读取速度场数据" << endl;
    // 定义一个向量场U
    volVectorField U(
        IOobject(
            "U",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE),
        mesh);
 
    // 定义一个场向量
    const vector originVector(0.05, 0.05, 0.005);
 
    // 计算向量originVector与各网格中心的距离中的最大值
    // 利用dimensionedVector将vector转换为具有长度量纲的向量
    // mag函数计算向量的模
    // value函数将数值转换为无量纲标量值
    const scalar rFarCell = max(
                                mag(dimensionedVector(
                                        "x0",
                                        dimLength,
                                        originVector) -
                                    mesh.C()))
                                .value();
 
    // 在案例迭代过程中进行场变量计算
    while (runTime.loop())
    {
        Info << "Time = " << runTime.timeName() << nl << endl;
        // 在所有的网格上循环
        for (label cellI = 0; cellI < mesh.C().size(); cellI++)
        {
            // 根据自定义的函数计算压力值
            p[cellI] = calculatedPressure(runTime.time().value(), mesh.C()[cellI],
                                          originVector, rFarCell);
        }
 
        // 计算压力梯度，并转换为速度，由于量纲不一致，所以需要乘以一个时间量纲进行转换
        // 注意在不可压缩求解器中，压力的量纲为[0 2 -2 0 0 0 0]，压力梯度量纲[0 1 -2 0 0 0 0]
        // 因此压力梯度乘上时间得到速度量纲[0 1 -1 0 0 0 0]
        U = fvc::grad(p) * dimensionedScalar("tmp", dimTime, 1.0);
 
        // 将数据写入到文件中，标记为AUTO_WRITE的场数据才可被写入
        runTime.write();
    }
 
    Info << "计算完成";
 
    return 0;
}
 
scalar calculatedPressure(scalar t, vector x, vector x0, scalar scale)
{
    scalar r(mag(x - x0) / scale);
    // 分母加1e-12是为了防止除以零
    scalar rR(1.0 / (r + 1e-12));
    scalar f(1.0);
    // 返回一个以x0为中心的正弦分布的压力
    return Foam::sin(2.0 * Foam::constant::mathematical::pi * f * t) * rR*100;
}
 

在demo4路径利用wmake编译程序。

3 测试程序

程序需要配合OpenFOAM案例运行。本案例已经预先准备了cavity案例，利用下面的命令运行程序：

cd  $FOAM_RUN/demo4/cavity
blockMesh
../demo4

运行结果如下图所示。

可以进入ParaView查看数据。