OpenFOAM 场(field)的操作和运算

Field Operation And Manipulation

field 相关类的结构

几个常见的类:

volScalarField
volVectorField
surfaceScalarField
surfaceVectorField

其实它们都是别名,定义如下:

typedef GeometricField<scalar, fvPatchField, volMesh> volScalarField;//src\finiteVolume\fields\volFields\volFieldsFwd.H
typedef GeometricField<vector, fvPatchField, volMesh> volVectorField;
typedef GeometricField<scalar, fvsPatchField, surfaceMesh> surfaceScalarField;//src\finiteVolume\fields\surfaceFields\surfaceFieldsFwd.H
typedef GeometricField<vector, fvsPatchField, surfaceMesh> surfaceVectorField;

发现,它们实际上都是 GeometricField,不过是提供的模板不同。
下面给出 GeometricField 类的 UML 类图。首先我们看到它需要三个模板:<Type,PatchField,GeoMesh>
结合 volScalarFieldsurfaceScalarField 的定义,我们不难发现:第一个 Type 可以是 scalar vector 等,第二个 PatchField 可以是 fvPatchField fvsPatchField(这里的 s 表示 surface) 等,第三个 GeoMesh 可以是 volMesh surfaceMesh 等。

GeometricField 类定义于 src\OpenFOAM\fields\GeometricFields\GeometricField

  • GeometricField – 场 + 网格,包含内部场及其边界场,边界场是场的场(FieldField)
  • DimensionedField – 带单位的场,只有内部场,没有边界场
  • Field – 数组(即 List)+ 代数操作

field 相关类的用法

构造

volScalarField A
(
IOobject
(
"A",
mesh.time().timeName(),
mesh, // 字典注册对象
IOobject::MUST_READ, // MUST_READ_IF_MODIFIED NO_READ READ_IF_PRESENT
IOobject::AUTO_WRITE, // NO_WRITE
true // 默认是注册
),
mesh
)
volScalarField B
(
IOobject
(
"B",
mesh.time().timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar(dimensionSet(1, -3, -1, 0, 0, 0, 0), 0.)
// 量纲还可以:dimless,或者 dimMass/dimVolume
)

上述两种构造函数对应 GeometricField 源代码中的:

GeometricField(const IOobject&,const Mesh&);
GeometricField(const IOobject&,const Mesh&,const dimensioned<Type>&,const word& patchFieldType=PatchField<Type>::calculatedType());

第一种是从文件中读取,内部场和边界场的初始值都由文件给定,边界条件也由文件给定。
第二种不需要从文件中读取,内部场和边界场的初始值都是这里给定的 0,边界条件默认是 calculated。也可以指定边界条件类型,如:

volScalarField C
(
IOobject
(
"C",
mesh.time().timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar(dimensionSet(1, -3, -1, 0, 0, 0, 0), 0.),
zeroGradientFvPatchScalarField::typeName
)

下面介绍第三种构造方式,以一个量为蓝本,构造第二个

volScalarField D
(
IOobject
(
"D",
mesh.time().timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
C
)

这个用法很有意思,当你提供了 D 时,它就会从文件读取;没提供时,它就会从 C 复制。
本人亲测。但是这个用法背后对应的代码有待挖掘。
它使用的是这个构造函数:


template<class Type, template<class> class PatchField, class GeoMesh>
Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
(
const IOobject& io,
const GeometricField<Type, PatchField, GeoMesh>& gf
)
:
Internal(io, gf),
timeIndex_(gf.timeIndex()),
field0Ptr_(nullptr),
fieldPrevIterPtr_(nullptr),
boundaryField_(*this, gf.boundaryField_)
{
if (debug)
{
InfoInFunction
<< "Constructing as copy resetting IO params"
<<endl << this->info() << endl;
}

if (!readIfPresent() && gf.field0Ptr_)
{
field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
(
io.name() + "_0",
*gf.field0Ptr_
);
}
}

其中 Internal 的是 DimensionedField 的别名,调用的这个构造函数:

template<class Type, class GeoMesh>
DimensionedField<Type, GeoMesh>::DimensionedField
(
const IOobject& io,
const DimensionedField<Type, GeoMesh>& df
)
:
regIOobject(io),
Field<Type>(df),
mesh_(df.mesh_),
dimensions_(df.dimensions_)
{}

需要注意一个场景,当你尝试这样使用时,是没问题的。

volScalarField B
(
IOobject
(
"B",
mesh.time().timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
A
);

但是如果你想的是继承 A 的边界类型,和量纲,但是数值想让它等于 0 的话。
你尝试这样做:

volScalarField B
(
IOobject
(
"B",
mesh.time().timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
0.*A
);

就会出现问题,有些边界类型会变成 calculated。
原因如下:
0.*A 调用的是定义于 GeometricFieldFunctionsM.C 中的这个函数

TEMPLATE                                                                       \
tmp<GeometricField<ReturnType, PatchField, GeoMesh>> operator Op \
( \
const GeometricField<Type1, PatchField, GeoMesh>& gf1 \
) \
{ \
tmp<GeometricField<ReturnType, PatchField, GeoMesh>> tRes \
( \
GeometricField<ReturnType, PatchField, GeoMesh>::New \
( \
#Op + gf1.name(), \
gf1.mesh(), \
Dfunc(gf1.dimensions()) \
) \
); \
\
Foam::OpFunc(tRes.ref(), gf1); \
\
return tRes; \
} \

它返回一个 tmp 的场,通过这个 New 函数创建的(GeometricField.H)。注意,它的边界类型使用的默认 calculated!

//- Return a temporary field constructed from name, mesh, dimensionSet
// and patch type.
static tmp<GeometricField<Type, PatchField, GeoMesh>> New
(
const word& name,
const Mesh&,
const dimensionSet&,
const word& patchFieldType=PatchField<Type>::calculatedType() // 这里有默认参数
);

template<class Type, template<class> class PatchField, class GeoMesh>
Foam::tmp<Foam::GeometricField<Type, PatchField, GeoMesh>>
Foam::GeometricField<Type, PatchField, GeoMesh>::New
(
const word& name,
const Mesh& mesh,
const dimensionSet& ds,
const word& patchFieldType // PatchField<Type>::calculatedType()
)
{
return tmp<GeometricField<Type, PatchField, GeoMesh>>
(
new GeometricField<Type, PatchField, GeoMesh>
(
IOobject
(
name,
mesh.thisDb().time().timeName(),
mesh.thisDb(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh,
ds,
patchFieldType // PatchField<Type>::calculatedType()
)
);
}

// 调用了这个构造函数:
template<class Type, template<class> class PatchField, class GeoMesh>
Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
(
const IOobject& io,
const Mesh& mesh,
const dimensionSet& ds,
const word& patchFieldType // PatchField<Type>::calculatedType()
)
:
Internal(io, mesh, ds, false),
timeIndex_(this->time().timeIndex()),
field0Ptr_(nullptr),
fieldPrevIterPtr_(nullptr),
boundaryField_(mesh.boundary(), *this, patchFieldType)
{
if (debug)
{
InfoInFunction <<"Creating temporary" << endl << this->info() << endl;
}

readIfPresent();
}

// 构造 boundaryField_ 时,调用的是:
template<class Type, template<class> class PatchField, class GeoMesh>
Foam::GeometricField<Type, PatchField, GeoMesh>::Boundary::
Boundary
(
const BoundaryMesh& bmesh,
const DimensionedField<Type, GeoMesh>& field,
const word& patchFieldType // PatchField<Type>::calculatedType()
)
:
FieldField<PatchField, Type>(bmesh.size()), // Boundary 的基类
bmesh_(bmesh) // 成员变量 const BoundaryMesh& bmesh_;
{
if (debug)
{
InfoInFunction << endl;
}

forAll(bmesh_, patchi)
{
this->set
(
patchi,
PatchField<Type>::New
(
patchFieldType, // PatchField<Type>::calculatedType()
bmesh_[patchi],
field
)
);
}
}


template<class Type>
Foam::tmp<Foam::fvPatchField<Type>> Foam::fvPatchField<Type>::New
(
const word& patchFieldType, // PatchField<Type>::calculatedType()
const fvPatch& p,
const DimensionedField<Type, volMesh>& iF
)
{
return New(patchFieldType, word::null, p, iF);
}



template<class Type>
Foam::tmp<Foam::fvPatchField<Type>> Foam::fvPatchField<Type>::New
(
const word& patchFieldType, // PatchField<Type>::calculatedType()
const word& actualPatchType, // word::null
const fvPatch& p, // 网格边界类型,比如 wall patch cyclic symmetry
const DimensionedField<Type, volMesh>& iF
)
{
if (debug)
{
InfoInFunction
<< "patchFieldType =" << patchFieldType
<<":" << p.type()
<< endl;
}

typename patchConstructorTable::iterator cstrIter =
patchConstructorTablePtr_->find(patchFieldType);

if (cstrIter == patchConstructorTablePtr_->end())
{
FatalErrorInFunction
<< "Unknown patchField type"
<< patchFieldType << nl << nl
<< "Valid patchField types are :" << endl
<<patchConstructorTablePtr_->sortedToc()
<<exit(FatalError);
}

typename patchConstructorTable::iterator patchTypeCstrIter =
patchConstructorTablePtr_->find(p.type()); // 查找网格的边界类型,wall patch 不在这个列表里边

if
(
actualPatchType == word::null
|| actualPatchType != p.type()
)
{
if (patchTypeCstrIter != patchConstructorTablePtr_->end())
{
return patchTypeCstrIter()(p, iF); // 如果 cyclic,symmetry,则调用 cyclicFvPatchField, symmetryFvPatchField 的构造函数
}
else
{
return cstrIter()(p, iF); // 如果 wall,patch,则调用 calculatedFvPatchField 的构造
}
}
else
{
tmp<fvPatchField<Type>> tfvp = cstrIter()(p, iF);

// Check if constraint type override and store patchType if so
if ((patchTypeCstrIter != patchConstructorTablePtr_->end()))
{
tfvp.ref().patchType() = actualPatchType;
}
return tfvp;
}
}

这里给出一个解决方案,换一个构造函数即可:

volScalarField B
(
IOobject
(
"B",
mesh.time().timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
0*A,
A.boundaryField()
);

使用的是这个构造函数:

template<class Type, template<class> class PatchField, class GeoMesh>
Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
(
const IOobject& io,
const Internal& diField,
const PtrList<PatchField<Type>>& ptfl
)
:
Internal(io, diField),
timeIndex_(this->time().timeIndex()),
field0Ptr_(nullptr),
fieldPrevIterPtr_(nullptr),
boundaryField_(this->mesh().boundary(), *this, ptfl)
{
if (debug)
{
InfoInFunction
<<"Constructing from components" << endl << this->info() << endl;
}

readIfPresent();
}

除了常规的 IOobject 构造,还有复制构造:

tmp<volScalarField> tmagGradP = mag(fvc::grad(p));
volScalarField normalisedGradP
(
"normalisedGradP",// 如果不指定名字,默认是什么?
tmagGradP()/max(tmagGradP())
);

获取

  • const

    • const Internal& internalField() const, 有量纲的内部场
    • const typename Internal::FieldType& primitiveField() const, 无量纲的内部场
    • const Boundary& boundaryField() const, 无量纲的边界场
  • non-const

    • Internal& ref(), 有量纲的内部场
    • typename Internal::FieldType& primitiveFieldRef(), 无量纲的内部场
    • Boundary& boundaryFieldRef(), 无量纲的边界场

forAll 来遍历的时候,遍历的其实就是那个 field<Type>
所以总体是有量纲的,用 forAll 对每个 cell 或 face 遍历的时候,就失去量纲了。

使用举例:

Info<<"T============"<<T<<endl; // 有量纲的 内部 + 边界
Info<<"T.internalField()============"<<T.internalField()<<endl; // 有量纲的 内部 + 边界
Info<<"T.primitiveField()============"<<T.primitiveField()<<endl; // 无量纲的 内部
Info<<"T.boundaryField()============"<<T.boundaryField()<<endl; // 无量纲的 边界
Info<<"T.ref()============"<<T.ref()<<endl; // 有量纲的 内部 + 边界
Info<<"T.primitiveFieldRef()============"<<T.primitiveFieldRef()<<endl; // 无量纲的 内部
Info<<"T.boundaryFieldRef()============"<<T.boundaryFieldRef()<<endl; // 无量纲的 边界

T 文件:

dimensions [0 0 0 1 0 0 0];

internalField uniform 900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

输出:

T============dimensions      [0 0 0 1 0 0 0];

internalField uniform 900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

T.internalField()============dimensions [0 0 0 1 0 0 0];

internalField uniform 900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

T.primitiveField()============8{900}
T.boundaryField()============
1
(
type fixedValue;
value uniform 1800;

)

T.ref()============dimensions [0 0 0 1 0 0 0];

internalField uniform 900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

T.primitiveFieldRef()============8{900}
T.boundaryFieldRef()============
1
(
type fixedValue;
value uniform 1800;

)

这里有一个奇怪的现象,即 ref 指的是内部场,但这里输出的确实内部 + 边界。不要被这里迷惑了,其它地方仍指的是内部场:

dimensionedScalar A(dimTemperature, 100.);
T.ref() = A;
Info<<"T============"<<T<<endl;

输出:

T============dimensions      [0 0 0 1 0 0 0];

internalField uniform 100;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

场的操作和运算 (Field Operation And Manipulation)

IOobject 中定义的函数

  • writeOpt
    使用举例:
mu_.writeOpt() = IOobject::AUTO_WRITE;

DimensionedField 中定义的函数

  • 平均
p.average() = gAverage(p)
  • 加权平均
p.weightedAverage(weightField) = gSum(weightField*p)/gSum(weightField)

使用举例:

Info<<"average p is"<<p.weightedAverage(mesh.V()).value()<<endl;
// 全场压力的体积加权平均值
  • dimensions
    使用举例:
Z_.dimensions().reset(dimless);

GeometricField 中定义的函数

  • max, min
    取最大值,包括内部场和边界场。并行计算时,Info 输出的是 master 的还是全局的?未测试
    使用举例:
Info<<"T==="<<max(T)<<endl;
  • writeMinMax
    输出最小最大值,只包括内部场。
template<class Type, template<class> class PatchField, class GeoMesh>
void Foam::GeometricField<Type, PatchField, GeoMesh>::writeMinMax
(
Ostream& os
) const
{
os <<"min/max(" << this->name() << ") ="
<<Foam::min(*this).value() << ","
<<Foam::max(*this).value()
<< endl;
}

使用举例:

Z_.writeMinMax(Info);

下面是详细的测试:
T 文件:

dimensions [0 0 0 1 0 0 0];

internalField uniform 900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

使用举例:

Info<<"max(T)==="<<max(T)<<endl; // 内部 + 边界
Info<<"min(T)==="<<min(T)<<endl; // 内部 + 边界
Info<<"max(T.internalField)==="<<max(T.internalField())<<endl; // 内部
Info<<"min(T.internalField)==="<<min(T.internalField())<<endl; // 内部
Info<<"max(T.ref)==="<<max(T.ref())<<endl; // 内部
Info<<"min(T.ref)==="<<min(T.ref())<<endl; // 内部
Info<<"max(T.boundaryField)==="<<max(T.boundaryField())<<endl; // 边界
Info<<"min(T.boundaryField)==="<<min(T.boundaryField())<<endl; // 边界
T.writeMinMax(Info); // 内部
T.max(1000); // 如有小于 1000 的,令其等于 1000,对内部和边界都起作用
Info<<"T==="<<T<<endl;
T.max(1900);
Info<<"T==="<<T<<endl;

输出:

max(T)=== max(T) [0 0 0 1 0 0 0] 1800
min(T)=== min(T) [0 0 0 1 0 0 0] 900
max(T.internalField)=== max(T) [0 0 0 1 0 0 0] 900
min(T.internalField)=== min(T) [0 0 0 1 0 0 0] 900
max(T.ref)=== max(T) [0 0 0 1 0 0 0] 900
min(T.ref)=== min(T) [0 0 0 1 0 0 0] 900
max(T.boundaryField)=== 1800
min(T.boundaryField)=== 1800
min/max(T) = 900, 900
T=== dimensions [0 0 0 1 0 0 0];

internalField uniform 1000;

boundaryField
{
walls
{
type fixedValue;
value uniform 1800;
}
}

T=== dimensions [0 0 0 1 0 0 0];

internalField uniform 1900;

boundaryField
{
walls
{
type fixedValue;
value uniform 1900;
}
}

上边的 T.max(1000) 调用的是这个函数

template<class Type, template<class> class PatchField, class GeoMesh>
void Foam::GeometricField<Type, PatchField, GeoMesh>::max
(
const dimensioned<Type>& dt
)
{
Foam::max(primitiveFieldRef(), primitiveField(), dt.value());
Foam::max(boundaryFieldRef(), boundaryField(), dt.value());
}

我们之所以可以直接给它一个 scalar 实参(即上边的 1000),是因为有这样一个构造函数的存在:

template<class Type>
Foam::dimensioned<Type>::dimensioned(const Type& t)
:
name_(::Foam::name(t)),
dimensions_(dimless),
value_(t)
{}

上边的 Foam::max 定义于 FieldFunctions.C 中的 BINARY_TYPE_FUNCTION(Type, Type, Type, max)
这个宏在 FieldFunctionsM.H 中

#define BINARY_TYPE_FUNCTION_SF(ReturnType, Type1, Type2, Func)                \
\
TEMPLATE \
void Func \
( \
Field<ReturnType>& res, \
const Type1& s1, \
const UList<Type2>& f2 \
) \
{ \
TFOR_ALL_F_OP_FUNC_S_F \
( \
ReturnType, res, =, ::Foam::Func, Type1, s1, Type2, f2 \
) \
} \
\
TEMPLATE \
tmp<Field<ReturnType>> Func \
( \
const Type1& s1, \
const UList<Type2>& f2 \
) \
{ \
tmp<Field<ReturnType>> tRes(new Field<ReturnType>(f2.size())); \
Func(tRes.ref(), s1, f2); \
return tRes; \
} \
\
TEMPLATE \
tmp<Field<ReturnType>> Func \
( \
const Type1& s1, \
const tmp<Field<Type2>>& tf2 \
) \
{ \
tmp<Field<ReturnType>> tRes = reuseTmp<ReturnType, Type2>::New(tf2); \
Func(tRes.ref(), s1, tf2()); \
tf2.clear(); \
return tRes; \
}


#define BINARY_TYPE_FUNCTION_FS(ReturnType, Type1, Type2, Func) \
\
TEMPLATE \
void Func \
( \
Field<ReturnType>& res, \
const UList<Type1>& f1, \
const Type2& s2 \
) \
{ \
TFOR_ALL_F_OP_FUNC_F_S \
( \
ReturnType, res, =, ::Foam::Func, Type1, f1, Type2, s2 \
) \
} \
\
TEMPLATE \
tmp<Field<ReturnType>> Func \
( \
const UList<Type1>& f1, \
const Type2& s2 \
) \
{ \
tmp<Field<ReturnType>> tRes(new Field<ReturnType>(f1.size())); \
Func(tRes.ref(), f1, s2); \
return tRes; \
} \
\
TEMPLATE \
tmp<Field<ReturnType>> Func \
( \
const tmp<Field<Type1>>& tf1, \
const Type2& s2 \
) \
{ \
tmp<Field<ReturnType>> tRes = reuseTmp<ReturnType, Type1>::New(tf1); \
Func(tRes.ref(), tf1(), s2); \
tf1.clear(); \
return tRes; \
}


#define BINARY_TYPE_FUNCTION(ReturnType, Type1, Type2, Func) \
BINARY_TYPE_FUNCTION_SF(ReturnType, Type1, Type2, Func) \
BINARY_TYPE_FUNCTION_FS(ReturnType, Type1, Type2, Func)

TFOR_ALL_F_OP_FUNC_S_F 和 TFOR_ALL_F_OP_FUNC_F_S 定义于 FieldM.H

// member function : this f1 OP fUNC s, f2

#define TFOR_ALL_F_OP_FUNC_S_F(typeF1, f1, OP, FUNC, typeS, s, typeF2, f2) \
\
/* check the two fields have same Field<Type> mesh */ \
checkFields(f1, f2, "f1" #OP ""#FUNC"(s, f2)"); \
\
/* set access to f1, f2 and f3 at end of each field */ \
List_ACCESS(typeF1, f1, f1P); \
List_CONST_ACCESS(typeF2, f2, f2P); \
\
/* loop through fields performing f1 OP1 f2 OP2 f3 */ \
List_FOR_ALL(f1, i) \
List_ELEM(f1, f1P, i) OP FUNC((s), List_ELEM(f2, f2P, i)); \
List_END_FOR_ALL \


// member function : this f1 OP fUNC f2, s

#define TFOR_ALL_F_OP_FUNC_F_S(typeF1, f1, OP, FUNC, typeF2, f2, typeS, s) \
\
/* check the two fields have same Field<Type> mesh */ \
checkFields(f1, f2, "f1" #OP ""#FUNC"(f2, s)"); \
\
/* set access to f1, f2 and f3 at end of each field */ \
List_ACCESS(typeF1, f1, f1P); \
List_CONST_ACCESS(typeF2, f2, f2P); \
\
/* loop through fields performing f1 OP1 f2 OP2 f3 */ \
List_FOR_ALL(f1, i) \
List_ELEM(f1, f1P, i) OP FUNC(List_ELEM(f2, f2P, i), (s)); \
List_END_FOR_ALL

接下来调用了 doubleScalar.H 中的 MAXMINPOW(Scalar, Scalar, Scalar) 因为 max 定义于这个宏里边。

#define MAXMINPOW(retType, type1, type2)          \
\
MAXMIN(retType, type1, type2) \
\
inline double pow(const type1 s, const type2 e) \
{ \
return ::pow(Scalar(s), Scalar(e)); \
}

MAXMINPOW(Scalar, Scalar, Scalar)
MAXMINPOW(Scalar, Scalar, int)
MAXMINPOW(Scalar, int, Scalar)
MAXMINPOW(Scalar, Scalar, long)
MAXMINPOW(Scalar, long, Scalar)
MAXMINPOW(Scalar, Scalar, float)
MAXMINPOW(Scalar, float, Scalar)

#undef MAXMINPOW

MAXMIN 定义于 int.H

#define MAXMIN(retType, type1, type2)              \
\
inline retType max(const type1 s1, const type2 s2) \
{ \
return (s1> s2)? s1: s2; \
} \
\
inline retType min(const type1 s1, const type2 s2) \
{ \
return (s1 < s2)? s1: s2; \
}

FieldFunctions.C 中还定义了几个诸如此类的函数

#define TMP_UNARY_FUNCTION(ReturnType, Func)                                   \
\
template<class Type> \
ReturnType Func(const tmp<Field<Type>>& tf1) \
{ \
ReturnType res = Func(tf1()); \
tf1.clear(); \
return res; \
}

template<class Type>
Type max(const UList<Type>& f)
{
if (f.size())
{
Type Max(f[0]);
TFOR_ALL_S_OP_FUNC_F_S(Type, Max, =, max, Type, f, Type, Max)
return Max;
}
else
{
return pTraits<Type>::min;
}
}

TMP_UNARY_FUNCTION(Type, max)

那么 max(T) 调用的应该是这里。

  • gMax
    max 的并行版,但是返回的是无量纲的值?

FieldFunctions

一些底层的函数,一般用不到

DimensionedFieldFunctions

pow
sqr
magSqr
mag
cmptAv
gMax
gMin
gSum
gSumMag
gAverage
+
-
* 外积
/
^ 叉乘
& 点乘,内积
&& 双内积

GeometricFieldFunctions

pow
sqr
magSqr
mag
cmptAv
gMax 返回的是无量纲的?
gMin
gSum
gSumMag
gAverage
+
-
* 外积
/
^ 叉乘
& 点乘,内积
&& 双内积

fvc 中的操作

volumeIntegrate

volumeIntegrate 的功能是对每个网格的某个物理量,乘以其网格体积,然后形成一个新的场。

volumeIntegrate(df) = df.mesh().V()*df.field();

domainIntegrate

domainIntegrate 的功能是对某个物理量,乘以其网格体积,然后对所有网格求和。

domainIntegrate(df) = gSum(fvc::volumeIntegrate(df));

这里的 gSum 是对全场求和。

其它 fvc 中的操作

surfaceIntegrate
surfaceSum
Su
Sp
SuSp
snGrad
reconstruct
laplacian
grad
flux
div
DDt
curl
average
cellReduce
文章作者: Yan Zhang
文章链接: https://openfoam.top/fields/
版权声明: 本博客所有文章除特别声明外,均采用 CC BY-NC-SA 4.0 许可协议。转载请注明来自 OpenFOAM 成长之路
您的肯定会给我更大动力~