本文边幅哄骗fvMesh类获取网格信息的一些基本用法。

OpenFOAM应用进程中需要芜俚地与网格进行交互，比如获取网格面信息，得到网格中心坐标等。哄骗fvMesh类中的相应成员函数很容易已毕这类需求。

fvMesh类是一个极为复杂的类，其继承自polyMesh、lduMesh、fvSchemes、surfaceInterpolation、fvSolution及data类，因此不仅不错造访网格，还不错造访与网格干系的各式算法形势。其中对网格信息的造访功能（如造访网格面、得到网格限度信息等）主要继承自polyMesh类。

图片

与前边的案例相似，本案例只需一个源文献及Make文献夹，哄骗底下的呐喊创建文献结构。

cd $FOAM_RUNmkdir demo3cd demo3touch demo3.Cmkdir Makecd Maketouch filestouch optionscd ..

文献结构如下图所示。

图片

demo3.CEXE = demo3

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

本案例触及到网格操作，因此需要准备网格文献。简单起见，使用一个OpenFOAM案例库中的文献进行测试，如使用cavity。

cp -r $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity .cd cavityblockMesh

践诺blockMesh生成网格后不错稽察网格信息，如下图所示。

图片

看到网格数目为400，网格面数为1640。

哄骗fvMesh类的C()及Cf()函数不错获取网格体中心坐标列表及里面面中心列表。

裁剪源文献demo3.C，输入底下的代码：

#include "fvCFD.H" int main(int argc, char *argv[]){#include "setRootCase.H"#include "createTime.H"// 在头文献createMesh.H中创建了一个名为mesh的fvMesh对象#include "createMesh.H"     // 输出刻下网格信息    // 哄骗mesh.C()得到网格单位体心坐标列表，mesh.Cf()得到网格里面面中心坐标列表    Info << "刻下时辰文献：" << runTime.timeName() << nl         << "网格数目：" << mesh.C().size() << nl         << "网格面数目：" << mesh.Cf().size() << nl << nl;     // 不错用轮回输出一些网格坐标信息    for (int cellI = 0; cellI < mesh.C().size(); cellI++)    {        // 由于网格数目太多，这里作念一些不祥,远离20输出一个        if (cellI % 20 == 0)        {            Info << "Cell[" << cellI << "]的中心坐标为：" << mesh.C()[cellI] << endl;        }    }    Info << endl;    return 0;}

编译上头的轮番并运转。

图片

详确需要到测试旅途（本案例为cavity文献夹）下践诺此轮番。运转后果如下图所示。

图片

不错看到，轮番正确地获取了网格数目与网格面的数目，并输出了网格中心坐标。

若念念要得到网格面的owner及neighbour网格索引，不错使用mesh.owner()及mesh.neighbour()函数。

在前边的代码中添加底下施行：

    // 不错使用mesh.owner()及mesh.neighbour()得到owner及neightbour的网格编号    for (label faceI = 0; faceI < mesh.owner().size(); faceI++)    {        if (faceI % 100 == 0)        {            Info << "里面面[" << faceI << "]的中心坐标为：" << mesh.C()[faceI] << nl                 << "owner[" << faceI << "]为：" << mesh.owner()[faceI] << nl                 << "neighbour[" << faceI << "]为：" << mesh.neighbour()[faceI] << nl;        }    }    Info << endl;

编译践诺，如下图所示。

图片

哄骗mesh.boundaryMesh()函数不错获取限度面上的信息。

添加以下代码。

// 哄骗mesh.boundaryMesh()函数获取限度面信息// forAll是OpenFOAM界说的轮回宏forAll(mesh.boundaryMesh(), patchI){    Info << "Patch[" << patchI << "]：" << mesh.boundary()[patchI].name() << "包括"        << mesh.boundary()[patchI].Cf().size() << "个网格面" << nl        << "肇端面为：" << mesh.boundary()[patchI].start() << endl;}Info << endl;

编译并运转，后果如下图所示。

图片

除了frontAndBack的网格面数外，其他限度的信息齐是莫得问题的。至于为什么会有问题，因为该限度类型为empty导致，管制轮番在文末。

若念念要获取限度面上某个网格的面积，不错使用函数sf()，如底下的代码：

// 哄骗函数faceCells()得到网格面相邻的网格// 哄骗函数sf()不错得到网格面的面积向量，该向量的模即为面积//界说一个索引，仅为测试使用，要是要得到几何面积，则需要将该限度上的统统网格面积乞降label patchFaceI(0);forAll(mesh.boundaryMesh(), patchI){    Info << "Patch[" << patchI << "]的网格面" << patchFaceI << "与网格["        << mesh.boundary()[patchI].patch().faceCells()[patchFaceI] <<"]相邻,"        << "其面积向量：" << mesh.boundary()[patchI].Sf()[patchFaceI] << ",其面积为："        << mag(mesh.boundary()[patchI].Sf()[patchFaceI])        << endl;     }Info << endl;

轮番运转后果如下图所示。

图片

关于里面面，不错凯旋在mesh实例上调用Sf()轮番。若要获取网格面的面积，不错使用函数magSf()，其等同于mag(Sf())，复返一个标量面积值。

关于里面面，法向量从owner指向neightbour，何况owner的索引小于neighbour。关于限度面，法向量老是指向缱绻域外部。

不错更详确地稽察组成每个网格面的点。最初，咱们通过获取对网格中各个对象的援用来界说一些变量。这些变量被界说为const，因为咱们不贪图以任何样子变嫌网格。详确：这些列表指的是网格的物理界说，因此包括限度面。不错使用mesh.edge()[patchi].Cf().size()和mesh.edge()[patchi].start()轮番搜检网格面是位于里面已经位于限度上。

如底下的代码：

// 稽察网格细节// 将变量界说为const，则此变量只可读取，而弗成修改const faceList &fcs = mesh.faces();            // 得到统统的网格面列表const List<point> &pts = mesh.points();        //得到统统的网格节点列表const List<point> &cents = mesh.faceCentres(); //得到网格面心点列表 // 遍历索引能被80整除的网格面forAll(fcs, faceI) if (faceI % 80 == 0){    if (faceI < mesh.Cf().size())    {        Info << "里面面";    }    else    {        forAll(mesh.boundary(), patchI) if ((mesh.boundary()[patchI].start() <= faceI) &&                                            (faceI < mesh.boundary()[patchI].start() + mesh.boundary()[patchI].Cf().size()))        {            Info << "网格面在限度[" << patchI << "]上，网格面[";            break;        }    }     Info << faceI << "]的中心为：" << cents[faceI] << "，领有" << fcs[faceI].size() << "个节点：";     forAll(fcs[faceI], vertexI)    {        Info << " " << pts[fcs[faceI][vertexI]];    }    Info << endl;}Info << endl; return 0;}

运转 后果如下图所示。

图片

在本测试案例中，FrontAndBack限度被界说为empty类型。这是一种额外的限度，在使用上头的函数时会出现一些出东说念主预念念的问题，如调用Cf函数时会复返空列表，导致size()函数的值为零。为了幸免出现此类问题，在进行Cf()函数调用时，不错先搜检限度类型。

如底下的代码：

// 判断限度类型是否为empty类型forAll(mesh.boundaryMesh(), patchID){    const polyPatch &pp = mesh.boundaryMesh()[patchID];    if (isA<emptyPolyPatch>(pp))    {        Info << "patch[" << patchID <<"]:" << mesh.boundary()[patchID].name()            << "的限度类型为empty" << endl;    }    else    {        Info << "patch[" << patchID << "]:"<< mesh.boundary()[patchID].name()            <<"的限度类型不是empty" << endl;    }}

运转后果如下图所示。

图片

也不错通过patch的称号哄骗函数findPatchID获取其ID，如底下的代码：

// 哄骗限度称号得到其IDword patchName("movingWall");label patchID = mesh.boundaryMesh().findPatchID(patchName);Info << "patch name:" << patchName << nl << "patchID:" << patchID << endl;

运转后果如下图所示。

图片

本案例的十足轮番代码如下所示。

#include "fvCFD.H" int main(int argc, char *argv[]){#include "setRootCase.H"#include "createTime.H"// 在头文献createMesh.H中创建了一个名为mesh的fvMesh对象#include "createMesh.H"     // 输出刻下网格信息    // 哄骗mesh.C()得到网格单位体心坐标列表，mesh.Cf()得到网格里面面中心坐标列表    Info << "刻下时辰文献：" << runTime.timeName() << nl         << "网格数目：" << mesh.C().size() << nl         << "网格面数目：" << mesh.Cf().size() << nl << nl;     // 不错用轮回输出一些网格坐标信息    for (int cellI = 0; cellI < mesh.C().size(); cellI++)    {        // 由于网格数目太多，这里作念一些不祥,远离20输出一个        if (cellI % 20 == 0)        {            Info << "Cell[" << cellI << "]的中心坐标为：" << mesh.C()[cellI] << endl;        }    }    Info << endl;     // 不错使用mesh.owner()及mesh.neighbour()得到owner及neightbour的网格中心坐标列表    for (label faceI = 0; faceI < mesh.owner().size(); faceI++)    {        if (faceI % 100 == 0)        {            Info << "里面面[" << faceI << "]的中心坐标为：" << mesh.C()[faceI] << nl                 << "owner[" << faceI << "]为：" << mesh.owner()[faceI] << nl                 << "neighbour[" << faceI << "]为：" << mesh.neighbour()[faceI] << nl;        }    }    Info << endl;     // 哄骗mesh.boundaryMesh()函数获取限度面信息    // forAll是OpenFOAM界说的轮回宏    forAll(mesh.boundaryMesh(), patchI)    {        Info << "Patch[" << patchI << "]：" << mesh.boundary()[patchI].name() << "包括"             << mesh.boundary()[patchI].Cf().size() << "个网格面" << nl             << "肇端面为：" << mesh.boundary()[patchI].start() << endl;    }    Info << endl;     // 哄骗函数faceCells()得到网格面相邻的网格    // 哄骗函数sf()不错得到网格面的面积向量，该向量的模即为面积    //界说一个索引，仅为测试使用，要是要得到几何面积，则需要将该限度上的统统网格面积乞降    label patchFaceI(0);    forAll(mesh.boundaryMesh(), patchI)    {        Info << "Patch[" << patchI << "]的网格面" << patchFaceI << "与网格["             << mesh.boundary()[patchI].patch().faceCells()[patchFaceI] << "]相邻,"             << "其面积向量：" << mesh.boundary()[patchI].Sf()[patchFaceI] << ",其面积为："             << mag(mesh.boundary()[patchI].Sf()[patchFaceI])             << endl;    }    Info << endl;     // 稽察网格细节    // 将变量界说为const，则此变量只可读取，而弗成修改    const faceList &fcs = mesh.faces();            // 得到统统的网格面列表    const List<point> &pts = mesh.points();        //得到统统的网格节点列表    const List<point> &cents = mesh.faceCentres(); //得到网格面心点列表     // 遍历索引能被80整除的网格面    forAll(fcs, faceI) if (faceI % 80 == 0)    {        if (faceI < mesh.Cf().size())        {            Info << "里面面";        }        else        {            forAll(mesh.boundary(), patchI) if ((mesh.boundary()[patchI].start() <= faceI) &&                                                (faceI < mesh.boundary()[patchI].start() + mesh.boundary()[patchI].Cf().size()))            {                Info << "网格面在限度[" << patchI << "]上，网格面[";                break;            }        }         Info << faceI << "]的中心为：" << cents[faceI] << "，领有" << fcs[faceI].size() << "个节点：";         forAll(fcs[faceI], vertexI)        {            Info << " " << pts[fcs[faceI][vertexI]];        }        Info << endl;    }    Info << endl;     // 判断限度类型是否为empty类型    forAll(mesh.boundaryMesh(), patchID)    {        const polyPatch &pp = mesh.boundaryMesh()[patchID];        if (isA<emptyPolyPatch>(pp))        {            Info << "patch[" << patchID << "]:" << mesh.boundary()[patchID].name()                 << "的限度类型为empty" << endl;        }        else        {            Info << "patch[" << patchID << "]:" << mesh.boundary()[patchID].name()                 << "的限度类型不是empty" << endl;        }    }    Info << endl;     // 哄骗限度称号得到其ID    word patchName("movingWall");    label patchID = mesh.boundaryMesh().findPatchID(patchName);    Info << "patch name:" << patchName << nl << "patchID:" << patchID << endl;     Info << endl;     return 0;} 

（本文完结）