simpleFoam是一个稳态流动求解器，可以求解层流与湍流流动问题。但该求解器中并不包含温度场的计算。本文描述在simpleFoam求解器基础上中添加温度场求解功能，将其改造为新求解器TsimpleFoam的基本过程。

内容参考：https://www.cfd-online.com/Forums/openfoam-programming-development/84480-adding-temperature-simplefoam.html

1 求解器改造

• 利用下面的命令拷贝icoFoam源代码

cd $FOAM_SOLVERS/incompressible
mkdir -p $WM_PROJECT_USER_DIR/applications/solvers
cp -r icoFoam $WM_PROJECT_USER_DIR/applications/solvers/TsimpleFoam
cd $WM_PROJECT_USER_DIR/applications/solvers/TsimpleFoam
# 删除多余的求解器文件
rm -r SRFSimpleFoam porousSimpleFoam/
mv simpleFoam.C TsimpleFoam.C

处理完毕后，icoTempFoam文件夹中的文件组织结构如下图所示。

1. 修改Make/files文件

文件内容修改为：

TsimpleFoam.C
EXE = $(FOAM_APPBIN)/TsimpleFoam

1. 修改createFields.H文件

在createFields.H文件中添加必须的场数据读取。这里需要添加热扩散系数DT与温度场T。

添加以下代码增加DT：

// 添加热扩散系数DT
dimensionedScalar DT
{
    "DT",
    dimViscosity,
    transportProperties.lookup("DT")
};

添加以下代码增加温度场T的读写：

//===增加温度场T===================
Info<< "Reading field T\n" << endl;
volScalarField T
(
    IOobject
    (
        "T",
        runTime.timeName(),
        mesh,
        IOobject::MUST_READ,
        IOobject::AUTO_WRITE
    ),
    mesh
);
//==============================

1. 修改TsimpleFoam.C文件

增加以下温度场求解功能。

//============增加温度求解功能====
fvScalarMatrix TEqn
(
    fvm::ddt(T)  fvm::div(phi,T)- fvm::laplacian(DT,T)
);
TEqn.relax();
TEqn.solve();
//===============================

1. 编译求解器

wmake

如下图所示编译过程中没有出现错误提示的话，即表示编译成功。

2 测试求解器

利用simpleFoam官方算例pitzDaily进行测试。

run
cp -r $FOAM_TUTORIALS/incompressible/simpleFoam/pitzDaily .
mv pitzDaily testCase
cd testCase/

删除多余的文件，最终剩下的文件结构如下所示。

2.1 修改constant文件夹

1. 修改transportProperties文件

文件内容修改为：

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "constant";
    object      transportProperties;
}
// * * * * * * * * * * * * * * * * * * //
 
transportModel  Newtonian;
 
nu              [0 2 -1 0 0 0 0] 1e-05;
DT              [0 2 -1 0 0 0 0] 3e-6;

这里的DT为热扩散系数：

其中，为热传导系数；为定压比热容。

1. 修改momentumTransport文件

采用层流计算。

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "constant";
    object      momentumTransport;
}
// * * * * * * * * * * * * * * * * * //
// 采用层流计算
simulationType laminar;
 
RAS
{
    // Tested with kEpsilon, realizableKE,
    // kOmega, kOmegaSST, v2f,
    // ShihQuadraticKE, LienCubicKE.
    model           kEpsilon;
    turbulence      off;
    printCoeffs     on;
}

2.2 设置0文件夹

在0文件夹下增加一个T文件。

cp 0/p 0/T

修改T文件内容，如下所示。

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      T;
}
// * * * * * * * * * * * * * * * * * * * //
// 注意修改温度的单位
dimensions      [0 0 0 1 0 0 0];
internalField   uniform 300;
 
boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform 300;
    }
 
    outlet
    {
        type            zeroGradient;
    }
 
    upperWall
    {
        type            fixedValue;
        value           uniform 350;
    }
 
    lowerWall
    {
        type            fixedValue;
        value           uniform 320;
    }
 
    frontAndBack
    {
        type            empty;
    }
}

2.3 修改system文件夹

利用命令foamGet residuals往system文件夹中添加residuals文件。

1. residuals文件

在residuals文件中添加温度T的输出。

fields (p U T);

1. controlDict文件

文件内容如下所示。

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * //
 
application     TsimpleFoam;
startFrom       startTime;
startTime       0;
stopAt          endTime;
endTime         1000;
deltaT          1;
writeControl    timeStep;
writeInterval   100;
purgeWrite      0;
writeFormat     ascii;
writePrecision  6;
writeCompression off;
timeFormat      general;
timePrecision   6;
runTimeModifiable true;
functions
{
    #includeFunc residuals
}

1. fvSchemes文件

fvSchemes文件中指定离散方法。文件内容如下所示。

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * //
 
ddtSchemes
{
    default         steadyState;
}
 
gradSchemes
{
    default         Gauss linear;
}
 
divSchemes
{
    default         none;
    div(phi,U)      bounded Gauss linearUpwind grad(U);
    div(phi,k)      bounded Gauss limitedLinear 1;
    div(phi,epsilon) bounded Gauss limitedLinear 1;
    div(phi,omega)  bounded Gauss limitedLinear 1;
    div(phi,v2)     bounded Gauss limitedLinear 1;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
    div(nonlinearStress) Gauss linear;
 
    div((nuEff*dev(grad(U).T()))) Gauss linear;   
    div(phi,T) Gauss upwind;
}
 
laplacianSchemes
{
    default         Gauss linear corrected;
    laplacian(DT,T) Gauss linear corrected;
}
 
interpolationSchemes
{
    default         linear;
}
 
snGradSchemes
{
    default         corrected;
}
 
wallDist
{
    method meshWave;
}

1. fvSolutions文件

该文件修改为：

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * //
 
solvers
{
    p
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.01;
        smoother        GaussSeidel;
    }
 
    T
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-07;
        relTol          0.01;
    }
 
    U
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0.01;
    }
}
 
SIMPLE
{
    nNonOrthogonalCorrectors 1;
    consistent      yes;
 
    residualControl
    {
        p               1e-3;
        U               1e-3;
        T               1e-6;    
    }
}
 
relaxationFactors
{
    equations
    {
        T               0.9;
        U               0.7; // 0.9 is more stable but 0.95 more convergent
        p               0.3; // 0.9 is more stable but 0.95 more convergent
    }
}

2.4 迭代计算

使用前面编译的求解器TsimpleFoam进行计算。

blockMesh
TsimpleFoam

计算迭代869次后达到收敛。

计算时可以查看残差。

foamMonitor -l postProcessing/residuals/0/residuals.dat

残差如下图所示。

2.5 计算结果

• 温度场分布

• 速度场分布

声明：原创文章，欢迎留言与我讨论，如需转载留言