OpenFOAM小代码

李东岳

这些是我感觉对初学者有用的一些代码，收藏在此，持续更新。有问题回帖提问，我会对我贡献的代码负责回答。

计算Patch的面积

    // 计算某个patch的面积 http://www.cfd-china.com/topic/3498
    {
        label patchID = mesh.boundaryMesh().findPatchID("inlet");

        // Create a polyPatch for looping
        const polyPatch& myPatch = mesh.boundaryMesh()[patchID];
        
        // Initialize patchArea
        scalar patchArea = 0.0;
        
        // Loop trhough all faces on the polyPatch, adding their magnitude surface
        // area vectors
        forAll(myPatch, faceI)
        {
            patchArea += mesh.magSf().boundaryField()[patchID][faceI];
        } 
    }

获得Patch的名字

Info << "patch name is " << U.boundaryField()[patch].patch().name() << nl;

网格相关量

const scalarField& V = mesh.V();                            // 网格体积
const surfaceVectorField& Sf = mesh.Sf();                   // 网格面矢量
const surfaceScalarField& magSf = mesh.magSf();             // 网格面积的模
const surfaceVectorField& N = Sf/magSf;                     // 网格面法向

更改壁面一层网格值

        label patchID = mesh.boundaryMesh().findPatchID("wall");
        const scalarField TfaceCell = 
            T.boundaryField()[patchID].patchInternalField();

        k.boundaryFieldRef()[patchID] == sqrt(TfaceCell);

#include "wallFvPatch.H"
forAll(nu.boundaryField(),patchi)
{
    if (isA<wallFvPatch>(mesh.boundary()[patchi]))
    {
        nu.boundaryFieldRef()[patchi] == 0.5*nu.boundaryField()[patchi].patchInternalField();
    }
}

给定patch调用网格

如果p是一个fvPatch

const fvMesh& mesh = p.boundaryMesh().mesh();

类内调用时间步长

假定mesh_是类内可以调用的成员变量。

const scalar deltaT = mesh_.time().deltaT().value();

如果mesh_不可获取，U_可获取：

const scalar deltaT = U_.mesh().time().deltaT().value();

声明一个tmp变量

tmp<volScalarField> deltaLambdaT   
(
    volScalarField::New
    (
        "deltaLambdaT",
        mesh_,
        dimensionedScalar(dimless, 1.0)
    )   
);

简写边界场

volScalarField::Boundary& meanUb = meanU_.boundaryFieldRef();

forAll(meanUb, patchi)
{
    fvPatchScalarField test = meanUb[patchi];
}

判断是否是固定值边界、processor边界

    forAll(mesh_.boundary(), P)                                              
    {                                                                        
        if
        (
            isA<fixedValueFvPatchScalarField>
            (
                T.boundaryField()[P]
            ) 
         ||
            isA<processorPolyPatch>(mesh_.boundaryMesh()[P])
        )
        {
        }
    }

codedFixedValue

inlet
    {
        type            codedFixedValue;
        value           uniform (0 0 0); //default value
        name           linearTBC; //name of new BC type
        code
        #{
            //const fvPatch& boundaryPatch = patch();
            //const fvBoundaryMesh& boundaryMesh = boundaryPatch.boundaryMesh();
            //const fvMesh& mesh = boundaryMesh.mesh();
            //const scalar &t = this->db().time().deltaTValue();
            
            scalar inletS = 0.0;
            const fvPatch& inletPatch = this->patch();
            
            forAll(inletPatch, faceI)
            {
                inletS += inletPatch.magSf()[faceI];
            } 

            scalar outletS = 0.0019304;
            scalar alphain = 0.1;

            vectorField& vf = *this; 

            forAll(vf, i)
            { 
                vf[i].y() = 0.01*outletS/inletS/alphain;
                vf[i].x() = 0.0; 
                vf[i].z() = 0.0; 
            }
        #};
    }

努塞尔数

下列代码可以添加到controlDict中执行


functions
{
    coded
    {
        functionObjectLibs ( "libutilityFunctionObjects.so" );
        enabled         true;
        type            coded;
        redirectType    printMinU;
        writeControl   timeStep;
        writeInterval  1;

        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{
            const volScalarField& T
            (
                mesh().lookupObject<volScalarField>("T")
            );

            const fvPatchList& patches = mesh().boundary();

            forAll(patches, patchi)
            {
                const fvPatch& currPatch = patches[patchi];

                if (mesh().time().value() == 30001)
                {
                    if (currPatch.name() == "downwall")
                    {
                        fvPatchScalarField nus = T.boundaryField()[patchi];

                        scalar L = 0.014231;
                        scalar Tref = 309.1;
                        scalar Timp = 347.6;
                        
                        nus = T.boundaryField()[patchi].snGrad()*L/(Timp - Tref);

                        forAll(T.boundaryField()[patchi], facei)
                        {
                            Pout << mesh().C().boundaryField()[patchi][facei].x()/0.014231
                                 << " " << nus[facei] << nl;
                        }
                    }
                }
            } 
       #};
    }
}

获取某个函数的计算时间

#include <chrono>

auto start = std::chrono::steady_clock::now();                                                   
// Functions here
auto end = std::chrono::steady_clock::now();                                                     
auto diff = end - start;
Info<< "Calculate nodes and weights, using " 
    << std::chrono::duration <double, std::milli> (diff).count()                                 
    << " ms" << endl;

在某个文件夹下输出scalarField

by @Samuel-Tu

    IOField<scalar> utau
    (
        IOobject
        (
            "utau",
             runTime.constant(),
            "../postProcessing",
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE 
        ),
        scalarField(totalFSize,0.0)
    );

BlockMesh-simpleGrading

处理非均一网格

blocks
(
    hex (0 1 2 3 4 5 6 7) (100 300 100)
    simpleGrading
    (
        1                  // x-direction expansion ratio
        (
            (0.2 0.3 4)    // 20% y-dir, 30% cells, expansion = 4
            (0.6 0.4 1)    // 60% y-dir, 40% cells, expansion = 1
            (0.2 0.3 0.25) // 20% y-dir, 30% cells, expansion = 0.25 (1/4)
        )
        3                  // z-direction expansion ratio
    )
);

输出yplus

functions
{
       yplus
      {
         type                   yPlus;	 
         functionObjectLibs      ("libutilityFunctionObjects.so")；
         outputControl           outputTime；
         enabled                 true;
       }
}

by @Sloan

在流场中添加颗粒失踪

controlDict添加下述内容，同时需要在constant文件夹添加kinematicProperties字典文件

particles
    {
        libs        ("liblagrangianFunctionObjects.so");
        type        particles;
    }

在流场中添加欧拉标量传输

controlDict添加下述内容，同时需要在0文件夹下添加s标量场，以及fvScheme、fvSolution都要做适配

functions
{
    scalarTransport1
        {
            type            scalarTransport;
            libs            ("libsolverFunctionObjects.so");
    
            field           s;
            D               1e-4;
            fvOptions
            {
                s
                {
                    type            semiImplicitSource;
                    active          true;
                    selectionMode   points;
                    points          
                    (
                        //(0.57 0.001 0.325)//w/h=5
                        //(0.62 0.001 0.325)//w/h=5
                        (0.67 0.001 0.325)//w/h=5
                        //(0.72 0.001 0.325)//w/h=5
                        //(0.77 0.001 0.325)//w/h=5
                    );
                    volumeMode      absolute;
                        
                    sources
                    {
                        s
                        {
                            explicit 0.0001;
                            implicit 0;
                        }
                    }
                }
            }
        }
}

对网格进行赋值

forAll(T, C)//遍历网格，T为场
{
    T[C] = min(nu[C], 100);//对T的第C个网格，赋值第C个网格的nu与100的最小值
}

声明多个场PtrList

PtrList<surfaceScalarField> abPosCoeff(4);                                                           
                                                                                                     
forAll(abPosCoeff, i)
{                                                                                                    
    abPosCoeff.set
    (                                                                                                
        i,                                                                                           
        new surfaceScalarField                                                                       
        (                                                                                            
            IOobject                                                                                 
            (                                                                                        
                "abPosCoeff" + name(i),                                                              
                mesh.time().timeName(),                                                              
                mesh,                                                                                
                IOobject::NO_READ,                                                                   
                IOobject::NO_WRITE                                                                   
            ),                                                                                       
            mesh,                                                                                    
            dimensionedScalar(inv(dimLength), 0.0)                                                   
        )                                                                                            
    );                                                                                               
}

OpenFOAM的二维数组

    PtrList<PtrList<scalar>> te(4);
                            
    forAll(te, i)
    {            
        te.set(i, new PtrList<scalar>(4));
                              
        forAll(te[i], j)   
        {             
            te[i].set(j, new scalar (0.0));
        }                   
    }

输出最大值最小值

functions
{
    minMaxp
    {
        type        fieldMinMax;
        functionObjectLibs ("libfieldFunctionObjects.so");
        fields
        (
             U
        );
        location        no;
        writeControl    timeStep;
        writeInterval   1;
    }
}

dimensionedScalar rho(dimDensity, 1.0);
    dimensionedScalar U(dimVelocity, 1.0);
    dimensionedScalar Ain(dimArea, 0.000481);
    dimensionedScalar sin(dimless, 1);
    dimensionedScalar min = rho*U*Ain*sin;

    tmp<volScalarField> tDEF
    (
        volScalarField::New
        (
            type(),
            mesh_,
            dimensionedScalar(dimless, 0)
        )
    );

    tmp<volScalarField> mwj
    (
        volScalarField::New
        (
            type(),
            mesh_,
            dimensionedScalar(dimMass/dimTime, 0)
        )
    );

    const volScalarField& s = lookupObject<volScalarField>("s");
    const volScalarField::Boundary& sBf = s.boundaryField();

    volScalarField::Boundary& mwjBf = mwj.ref().boundaryFieldRef();

    scalar mw = 0.0;

    forAll(mesh_.boundary(), patchi)
    {
        if (mesh_.boundary()[patchi].name() == "wall")
        {
            scalarField Aj = mesh_.boundary()[patchi].magSf();

            mwjBf[patchi] = -rho.value()*Aj*sBf[patchi].snGrad()*1e-5;

            forAll(mesh_.boundaryMesh()[patchi], facei)
            {
	        mw += mwjBf[patchi][facei];
            }
        }
    }

    volScalarField::Boundary& DEFBf = tDEF.ref().boundaryFieldRef();

    forAll(mesh_.boundary(), patchi)
    {
        if (mesh_.boundary()[patchi].name() == "wall")
        {
            scalar area = gSum(mesh_.magSf().boundaryField()[patchi]);
	    Info<< "Wall area is " << area << nl;

            scalarField Aj = mesh_.boundary()[patchi].magSf();

            //DEFBf[patchi] = mwjBf[patchi]/Aj/(mw/A);
            DEFBf[patchi] = mwjBf[patchi]/Aj;
        }
    }

    return tDEF;

单独一个文件输出直径信息

meanDiameter
    {
        type            coded;
        libs            ("libutilityFunctionObjects.so");
        name            error;
    
        codeExecute
        #{
            const volScalarField& d =
                mesh().lookupObject<volScalarField>("d.alpha.oil");

            scalar d32 = d.weightedAverage(mesh().V()).value();
    
            if (Pstream::master())
            {
                std::ofstream file;
                file.open ("d32", std::ofstream::out | std::ofstream::app);
                file << mesh().time().timeName() << " " << d32 << "\n";
                file.close();
            }
        #};
    }

拉格朗日输出alpha

cloudFunctions 
{ 
    f1
    {
        type voidFraction; 
    } 
}

yfclark

补充一个CodeStream实现充分发展的湍流入口，请东岳老师不要建议
$\begin{array}{l}
\delta_{T}>\delta_{w} \rightarrow u=u_{\infty}\left(\frac{d_{w}}{\delta_{T}}\right)^{1 / 7} \
\delta_{T}<\delta_{w} \rightarrow u=u_{\infty}
\end{array}
\delta_{T}=0.37\left(v / u_{\infty}\right)^{0.2} L_{T}^{0.8}
$
dw为到壁面的最短距离

OpenFoam 7：

gasinlet
{
    type            fixedValue;
    value           #codeStream
    {
        codeInclude
        #{
            #include "fvCFD.H"
            #include "DynamicList.H"
        #};

        codeOptions
        #{
            -I$(LIB_SRC)/finiteVolume/lnInclude \
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

	    //libs needed to visualize BC in paraview
	    codeLibs
	    #{
    	    -lmeshTools \
    	    -lfiniteVolume
	    #};

        code
        #{
            const IOdictionary& d = static_cast<const IOdictionary&>
			(
                dict.parent().parent()
            );
            const fvMesh& mesh = refCast<const fvMesh>(d.db());
            const label id = mesh.boundary().findPatchID("gasinlet");
            const fvPatch& patch = mesh.boundary()[id];

            vectorField U(patch.size(), vector(0, 0, 0));

            const scalar U_inf=93.5;
            const scalar LT=0.2;
            const scalar visco=2.3183558929634845e-06;
            const scalar deltaT=0.37*pow(visco/U_inf,0.2)*pow(LT,0.8);
            scalar dw=0.0;
            forAll(U, i)
            {
                const scalar x = patch.Cf()[i][0];
                const scalar y = patch.Cf()[i][1];
                const scalar z = patch.Cf()[i][2];
                DynamicList<scalar>  distance(4);
                distance.append(mag(y));
                distance.append(mag(y-30e-6));
                distance.append(mag(z-23e-6));
                distance.append(mag(z+23e-6));
                dw=min(distance);
                distance.clearStorage();
                scalar Ux=0.0;
                if(deltaT>dw)
                {
                    Ux=U_inf*pow(dw/deltaT,1.0/7.0);
                }
                else
                {
                    Ux=U_inf;
                }
	            U[i] = vector(Ux, 0, 0);
            }

            writeEntry(os, "", U);
        #};
    };
}

Zhong-combustion

@yfclark 谢谢贡献代码。有一个疑问还请教，这里没有脉动量，只有平均速度的径向分布。在LES或者DNS中，脉动量比较重要，所以关于脉动量的话有什么看法么？还是说这个入口可以作为一个入口，然后在一个长管道内比较容易自然发展成充分发展的湍流？代码里的U_inf可以理解为平均速度么？

马乔

wordHashSet unMatchedkeys(dict.toc());
forAll(names, nameI)
{
    const word& name = names[nameI];
    const entry* ePtr = dict.findEntry(name, keyType::REGEX);
    if(ePtr)
       unMatchedKeys.erase(ePtr->keyWord());
}

马乔

class base
{
private:
    bool ownedByRegistry_;
public:
    template<class Type>
    static void store(Type* p)
    {
	//p->base::ownedByRegistry_ =true;
	p->ownedByRegistry_ = true;
	p->printInfo();
    }
};

class derived
:public base
{
public:
    derived() {}
    void printInfo()
    {
	Info << "derived class" << endl;
    }
};

int main(int argc, char *argv[])
{
...
    derived d;
    base::store(&d);
...
}

D.Benjamin

@东岳在 OpenFOAM小代码中说：

东岳老师，请教几个小问题：

在编译过程中报错：'mesh_' was not declared，所以mesh_是不可以直接使用的，但是我不太理解什么是类内可以调用的成员变量，感觉速度场U，温度场T等都是类内可以调用的成员变量？
如果mesh_不可获取，U_可获取。上面说mesh_是类内可以调用的成员变量，然而U_又是什么？
声明一个tmp变量，mesh_表示什么意思呢?在声明U、T、p等场的时候，mesh_的位置常常写的是mesh

期待回复，祝好！

李东岳

@D-Benjamin

我不太理解什么是类内可以调用的成员变量

我说的这个东西，就是类的成员变量。如果存在成员变量mesh, 它在OpenFOAM类内都被定义为mesh_，在类外那肯定是mesh，后面有没有_只是写法风格，OpenFOAM类里面的变量都带_，类外面的变量都没有_。

OpenFOAM这面类继承会涉及到好几层。你可以一层一层的找有没有mesh_的声明，你可以直接用mesh_看看是否报错，没报错就是声明了，报错了就是没生命。

霜染丹枫

@Zhong-combustion 你好，我目前也在用LES做加热管流的湍流研究，关于湍流入口的设置看了一些资料，如有兴趣，可以加个qq:1191364394，共同学习探讨。

Zhy2022

东岳老师您好，请问如果是两相流中相分数（比如alpha.particle）随时间变化，它在codeFixedValue中应该怎么定义呢？

李东岳

@Zhy2022 codedFixedValue上面不是有嘛，你举一反三一下咩

Zhy2022

@东岳东岳老师，没能找出报错原因，只能再来麻烦一下您：

INLET1
    {
        type            codedFixedValue;
        value           $internalField;
        name            stepInjection;
        redictType      variedValue;
        code
        #{
             scalar alphapar = this->patch().fluid().phase();

             scalar t = this->db().time().value();
       
             if(t >= 0 && t<=5)
             {
                  alphapar = 0;
             }
             
             else
	     {
	          alphapar = 0.15;
	     }

	     operator==(alphapar)

         #};
        
    }

Zhy2022

@东岳东岳老师，问题已解决！无法编译与文件路径有关

veen

fvOption小工具

根据GitHub上一个python脚本改写了一个生成fvOption的小工具 https://github.com/Veenxz/fvSchemes
用Foam没有多久，各位大佬可以多提提意见。

# version 2.0
# Based on https://github.com/Carlopasquinucci/fvSchemes
# Generation by Veenxz
# relaesed under license GPL GNU 3.0

steady = True
pseudo_transient = False
precision = 2    # First order 1 or Second order 2
unbounded = False
LUST = False
secondorder = True

maxOrtho = 80
maxSkew = 20

# Header and Footer
h = [
    "/*--------------------------------*- C++ -*----------------------------------*|",
    "| =========                 |                                                 |",
    "| \\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |",
    "|  \\\    /   O peration     | Website:  https://openfoam.org                  |",
    "|   \\\  /    A nd           | Version:  7                                     |",
    "|    \\\/     M anipulation  |                                                 |",
    "\*---------------------------------------------------------------------------*/",
    "FoamFile", "{", "    version     2.0;", "    format      ascii;",
    "    class       dictionary;", '    location    "system";', "    object      fvSchemes;",
    "}",
    "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //",
    ""
]
footer = '\n// ************************************************************************* //'

nonOrthogonalCorrectors = 1
if (maxOrtho) > 80:
    print(
        'Warning: mesh is not so nice. Use 3 nonOrthogonalCorrectors in fvSolution file'
    )
    nonOrthogonalCorrectors = 3
if (maxSkew) > 8:
    print('Warning: mesh is not so nice')

if (maxOrtho) > 80:

    gradSchemes = (
        '{\n    default          cellLimited Gauss linear 0.5;\n'
           '    grad(U)          faceLimited Gauss linear 1.0;\n}'
    )
    divSchemes = (
        '{\n    div(phi,U)       Gauss linearUpwind grad(U);\n'
           '    div(phi,omega)   Gauss upwind;\n'
           '    div(phi,k)       Gauss upwind;\n'
           '    div(phi,e)       Gauss upwind;\n'
           '    div((nuEff*dev(T(grad(U))))) Gauss linear;\n}'
    )
    laplacianSchemes = (
        '{\n    default          Gauss linear limited 0.333;\n}'
    )
    snGradSchemes = (
        '{\n    default          Gauss linear limited 0.333;\n}'
    )

    blending = 0.2
    nonOrthogonalCorrectors = 3

if (maxOrtho) > 70:

    gradSchemes = (
        '{\n    default          cellLimited Gauss linear 0.5;\n'
           '    grad(U)          cellLimited Gauss linear 1.0;\n}'
    )
    divSchemes = (
        '{\n    div(phi,U)       Gauss linearUpwind grad(U);\n'
           '    div(phi,omega)   Gauss linearUpwind grad(omega);\n'
           '    div(phi,k)       Gauss linearUpwind grad(k);\n'
           '    div(phi,e)       Gauss linearUpwind grad(e);\n'
           '    div((nuEff*dev(T(grad(U))))) Gauss linear;\n}'
    )
    laplacianSchemes = (
        '{\n    default          Gauss linear limited 0.5;\n}'
    )
    snGradSchemes = (
        '{\n    default          Gauss linear limited 0.5;\n}'
    )

    blending = 0.5
    nonOrthogonalCorrectors = 3

if (maxOrtho) > 60:

    gradSchemes = (
        '{\n    default          cellMDLimited Gauss linear 0.5;\n'
           '    grad(U)          cellMDLimited Gauss linear 0.5;\n}'
    )
    divSchemes = (
        '{\n    div(phi,U)       Gauss linearUpwind grad(U);\n'
           '    div(phi,omega)   Gauss linearUpwind grad(omega);\n'
           '    div(phi,k)       Gauss linearUpwind grad(k);\n'
           '    div(phi,e)       Gauss linearUpwind grad(e);\n'
           '    div((nuEff*dev(T(grad(U))))) Gauss linear;\n}'
    )
    laplacianSchemes = (
        '{\n    default          Gauss linear limited 0.777;\n} '
    )
    snGradSchemes = (
        '{\n    default          Gauss linear limited 0.777;\n} '
    )

    blending = 0.7
    nonOrthogonalCorrectors = 2

if (maxOrtho) > 0:

    gradSchemes = (
        '{\n    default          cellMDLimited Gauss linear 0;\n'
           '    grad(U)          cellMDLimited Gauss linear 0.333;\n}'
    )
    divSchemes = (
        '{\n    div(phi,U)       Gauss linearUpwind grad(U);\n'
           '    div(phi,omega)   Gauss linearUpwind grad(omega);\n'
           '    div(phi,k)       Gauss linearUpwind grad(k);\n'
           '    div(phi,e)       Gauss linearUpwind grad(e);\n'
           '    div((nuEff*dev(T(grad(U))))) Gauss linear;\n}'
    )
    laplacianSchemes = (
        '{\n    default          Gauss linear limited 0.95;\n}'
    )
    snGradSchemes = (
        '{\n    default          Gauss linear limited 0.95;\n}'
    )

    blending = 0.8
    nonOrthogonalCorrectors = 1

if (LUST):
    gradSchemes = (
        '{\n    default          Gauss linear;\n'
           '    grad(U)          cellMDLimited leastSquares 1;\n}'
    )
    divSchemes = (
        '{\n    div(phi,U)       Gauss LUST grad(U);\n'
           '    div(phi,omega)   Gauss LUST grad(omega);\n'
           '    div(phi,k)       Gauss LUST grad(k);\n'
           '    div(phi,e)       Gauss LUST grad(e);\n'
           '    div((nuEff*dev(T(grad(U))))) Gauss linear;\n}'
    )
    laplacianSchemes = (
        '{\n    default          Gauss linear corrected;\n}'
    )
    snGradSchemes = (
        '{\n    default          Gauss linear corrected;\n}'
    )
    
    blending = 0.9
    nonOrthogonalCorrectors = 1

if (steady):
    ddtSchemes = (
        '{\n    default          steadyState;\n}'
    )
    if (pseudo_transient):
                ddtSchemes = (
                    '{\n    default          localEuler;\n}'
                )
else:
    ddtSchemes = (
        '{\n    default           CrankNicolson ' + str(blending) + ' ;\n}'
    )
    if precision == 1:
        ddtSchemes = (
            '{\n    default           Euler;\n}'
        )
    if (unbounded):
        ddtSchemes = (
            '{\n    default           backward;\n}'
        )


wallDist = (
    "{\n    method           meshWave;\n}"
)

#open fvSchemes and write inside

f = open("fvSchemes", "w")

for i in h:
	f.write(i + "\n")

f.write("ddtSchemes" + "\n")
f.write(ddtSchemes + "\n")
f.write("\n")
f.write("gradSchemes" + "\n")
f.write(gradSchemes + "\n")
f.write("\n")
f.write("divSchemes" + "\n")
f.write(divSchemes + "\n")
f.write("\n")
f.write("laplacianSchemes" + "\n")
f.write(laplacianSchemes + "\n")
f.write("\n")
f.write("snGradSchemes" + "\n")
f.write(snGradSchemes + "\n")
f.write("\n")
f.write("wallDist" + "\n")
f.write(wallDist + "\n")
f.write(footer)

f.close()

print('File fvSchemes created')

李东岳

@veen 厉害厉害

veen

@veen 这个代码还有些格式没法用……我回头测试了再放到github大家看仓库里边的代码就好了。

hongjiewang

@李东岳在 OpenFOAM小代码中说：

IOField<scalar> utau
(
IOobject
(
"utau",
runTime.constant(),
"../postProcessing",
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
scalarField(totalFSize,0.0)
);

请问老师，如果我想要在postProcessing中输出一个标量，他只是一个数，并不是场量，我应该怎么定义他的类型。上面这种方法是不是只适应于场量的输出。我想输出的是下面c的数值。

const volScalarField& b = mesh().lookupObject<volScalarField>("alpha.liquid");
scalar c= b.weightedAverage(mesh().V()).value();

李东岳

@hongjiewang

    IOList<scalar> utau
    (
        IOobject
        (
            "utau",
            runTime.constant(),
            "../postProcessing",
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE 
        ),
        1
    );
    
    utau[0] = b.weightedAverage(mesh().V()).value();

hongjiewang

@李东岳

  totalLiquid
  {        
        libs            (utilityFunctionObjects);
        type            coded;
        name            totalLiquid;
        enabled         true;
        writeControl    timeStep;
        writeInterval   1;
        
        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{

            const volScalarField& b =
                mesh().lookupObject<volScalarField>("alpha.liquid");

            IOList<scalar> liquidFraction
            (
                IOobject
                (
                    "liquidFraction",
                    mesh().time().constant(),
                    "../postProcessing",
                    mesh(),
                    IOobject::NO_READ,
                    IOobject::AUTO_WRITE
                ),
                1
            );
            liquidFraction[0] = b.weightedAverage(mesh().V()).value();
        #};
    }

东岳老师，我这样添加到controlDict里，运行后postProcessing里并没有出现liquidFraction的值。
1618213472(1).png

hongjiewang

@hongjiewang 在 OpenFOAM小代码中说：

@李东岳

  totalLiquid
  {        
        libs            (utilityFunctionObjects);
        type            coded;
        name            totalLiquid;
        enabled         true;
        writeControl    timeStep;
        writeInterval   1;
        
        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{

            const volScalarField& b =
                mesh().lookupObject<volScalarField>("alpha.liquid");

            IOList<scalar> liquidFraction
            (
                IOobject
                (
                    "liquidFraction",
                    mesh().time().constant(),
                    "../postProcessing",
                    mesh(),
                    IOobject::NO_READ,
                    IOobject::AUTO_WRITE
                ),
                1
            );
            liquidFraction[0] = b.weightedAverage(mesh().V()).value();
        #};
    }

东岳老师，我这样添加到controlDict里，运行后postProcessing里并没有出现liquidFraction的值。
1618213472(1).png

需要修改两个部分，
1.在下方添加 .write（）
2.将codeExecute改为codeWrite
之后就可以得到想要的相含量结果。

hongjiewang

@李东岳

FoamFile
{
    version     2.0;
    format      ascii;
    class       scalarField;
    location    "constant/../postProcessing";
    object      totalLiquid;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

1(0.2)

麻烦老师再帮我看一下，为什么每次都只会输出一个时刻的值，不会随时间增加输出的

李东岳

meanDiameter
    {
        type            coded;
        libs            ("libutilityFunctionObjects.so");
        name            error;
    
        codeExecute
        #{
            const volScalarField& d =
                mesh().lookupObject<volScalarField>("d.alpha.oil");

            scalar d32 = d.weightedAverage(mesh().V()).value();
    
            if (Pstream::master())
            {
                std::ofstream file;
                file.open ("d32", std::ofstream::out | std::ofstream::app);
                file << mesh().time().timeName() << " " << d32 << "\n";
                file.close();
            }
        #};
    }

我最近也凑巧要用这个，我用的上面这种方式。

五好青年

@李东岳
请问东岳老师

(1) 东岳老师，您上面写的：codedFixedValue
(2) yfclark朋友，上面写的：CodeStream

在自定义边界条件方面，这两个有啥区别啊，我研究了几天没想明白

五好青年

@yfclark
朋友，您好：

codedFixedValue和CodeStream

你知道在自定义边界条件方面，这两个有啥区别啊，我研究了几天没想明白

李东岳

codedFixedValue更简单，直接用

code
        #{
            operator==(min(10, 0.1*this->db().time().value()));
        #};

就可以赋值。CodeStream更普适性，我见过一些代码用CodeStream对内部场进行初始化。比如建立U的内部场与空间的函数。

hongjiewang

@李东岳目前我需要将温度边界和压力边界定义为随相含量变化的一个场。我只可以做到变化之后的也是同一类边界条件，但是像下图中的气相区压力是第一类边界条件，液相区温度是第二类边界条件。请问这种变边界类型的边界应该怎么定义~ 1632279639(1).png

hongjiewang

@hongjiewang 在 OpenFOAM小代码中说：

@李东岳目前我需要将温度边界和压力边界定义为随相含量变化的一个场。我只可以做到变化之后的也是同一类边界条件，但是像下图中的气相区压力是第一类边界条件，液相区温度是第二类边界条件。请问这种变边界类型的边界应该怎么定义~

已解决压力的~但是温度的第二类和第三类的转换还没有解决~

疏影横斜水清浅

你好，请问const polyPatch& cPatch = mesh.boundaryMesh()[patchID]；与const fvPatch& cPatch = mesh.boundary()[patchID];有什么不同？

Dingcy

@李东岳在 OpenFOAM小代码中说：

//const fvBoundaryMesh& boundaryMesh = boundaryPatch.boundaryMesh();
//const fvMesh& mesh = boundaryMesh.mesh();

我想用codedFixedValue,写一个壁面的边界条件，壁面都是wall，而不是patch，所以这里想用这两个命令，但是编译过程中总是提示我是不是应该声明点啥呢？我也查了网上的，发现全是在inlet使用这个功能，然后入口是patch

李东岳

@dingcy

const polyPatch& cPatch = mesh.boundaryMesh()[patchID]

类似这个？

Dingcy

@李东岳好像不是这个又显示了mesh没有声明....

vbcwl

@李东岳能请教下这个具体怎么使用么？我没能输出utau

李东岳

@vbcwl 上面那些代码么。上面那些代码都需要自己写到代码里面去。还不是简单的使用。需要有编程基础，这些代码起到一定的提示作用。

2019201300

李老师您好，请问输出努塞尔数代码是直接放进去吗？还是里边的部分数据需要进行修改呢微信截图_20221112113853.png ？

李东岳

@2019201300 这个我太久之前写的。已经忘了。但是这个代码只是提供案例。具体数值大概率需要自己改一下。直接放在controlDict下面。

2019201300

@李东岳好的，谢谢李老师，我再认真研究研究。