网格体积

[zhanghan]

请问有哪位大神提取过网格体积吗？？？？

[赵一铭]

mesh.V()

可以输出网格体积

[飞火流星jyj]

@赵一铭 赵老师，那网格的坐标呢，x，y，z

[队长别开枪]

@飞火流星jyj

    //class: fvMesh
    Info << "\n-Class: fvMesh--------" << endl;

    //- Return the object registry - resolve conflict polyMesh/lduMesh.
        //  Type: virtual const objectRegistry &
    Info << mesh.thisDb() << endl;

    //- Return reference to name.
        //  Type: const word &
    Info << mesh.name() << endl;

    //- Return reference to boundary mesh.
        //  Type: const fvBoundaryMesh &
    mesh.boundary();

    //- Internal face owner.
        // Type: const labelUList &
    labelList owners(mesh.owner());

    //- Internal face neighbour.
        //  Type: const labelUList &
    labelList neighbours(mesh.neighbour());

    //- Return cell volumes.
        //  Type: const DimensionedField< scalar, volMesh > &
    Info << mesh.V() << endl;

    //- Return old-time cell volumes.
        //  Type: const DimensionedField< scalar, volMesh > &
    Info << mesh.V0() << endl;

    //- Return old-old-time cell volumes.
        //  Type: const DimensionedField< scalar, volMesh > &
    Info << mesh.V00() << endl;

    //- Return sub-cycle cell volumes.
        //  Type: tmp< DimensionedField< scalar, volMesh > >
    Info << mesh.Vsc() << endl;

    //- Return sub-cycle old-time cell volumes.
        //  Type: tmp< DimensionedField< scalar, volMesh > >
    Info << mesh.Vsc0() << endl;

    //- Return cell face area vectors.
        //  Type: const surfaceVectorField &
    Info << mesh.Sf() << endl;

    //- Return cell face area magnitudes.
        //  Type: const surfaceScalarField &
    Info << mesh.magSf() << endl;

    //- Return cell face motion fluxes.
        //  Type: const surfaceScalarField &
    Info << mesh.phi() << endl;

    //- Return cell centres as volVectorField.
        //  Type: const volVectorField &
    Info << mesh.C() << endl;

    //- Return face centres as surfaceVectorField.
        //  Type: const surfaceVectorField &
    Info << mesh.Cf() << endl;

    //- Return face deltas as surfaceVectorField.
        //  Type: tmp< surfaceVectorField >
    Info << mesh.delta() << endl;  

    Info << "----------------------\n" << endl;

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    //class: fvMesh

    pointField points(mesh.points());
    faceList faces(mesh.faces());
    cellList cells(mesh.cells());    
    


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    //Class: point = vector
    Info << "\n-Class: point---------" << endl;
    point &pt(points[3]);
    Info << "pt = " << pt << endl;
    //- Return x component
    Info << "pt.x() = " << pt.x() << endl;
    //- Return y component
    Info << "pt.y() = " << pt.y() << endl;
    //- Return z component
    Info << "pt.z() = " << pt.z() << endl;
    Info << "----------------------\n" << endl;


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    //class: edge
    Info << "\n-Class: edges----------" << endl;
    edge eg(faces[1].faceEdge(1)), eg_(faces[1].faceEdge(3));
    label pt_1(eg.start()), &pt_2(eg.end());
    Info << "eg = " << eg << endl;
    //- Return start vertex label
    Info << "eg.start() = " << pt_1 << endl;
    //- Return end vertex label
    Info << "eg.end() = " << pt_2 << endl;
    //- Given one vertex, return the other
    Info << "eg.otherVertex(eg.end()) = " << eg.otherVertex(pt_2) << endl;
    //- Return common vertex
        //  - -1: no common vertex
    Info << "eg.commonVertex(eg_) = " << eg.commonVertex(eg_) << endl;
    //- Return reverse edge
    Info << "eg.reverseEdge() = " << eg.reverseEdge() << endl;
    //- Return centre (centroid)
    Info << "eg.centre(points) = " << eg.centre(points) << endl;
    //- Return the vector (end - start)
    Info << "eg.vec(points) = " << eg.vec(points) << endl;
    //- Return scalar magnitude
    Info << "eg.mag(points) = " << eg.mag(points) << endl;
    //- Return edge line
    Info << "eg.line(points) = " << eg.line(points) << endl;
    //- compare edges
        //  Returns:
        //  -  0: different
        //  - +1: identical
        //  - -1: same edge, but different orientation
    Info << "Foam::edge::compare(eg, eg.reverseEdge()) = " 
         << Foam::edge::compare(eg, eg.reverseEdge()) << endl;
    Info << "----------------------\n" << endl;


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    //class: face
    Info << "\n-Class: face----------" << endl;
    face &fe(faces[4]);
    Info << "fe = " << fe << endl;
    //- Return true if the face is empty
    Info << "fe.empty() = " << fe.empty() << endl;
    //- Return No. of points corresponding to this face
    Info << "fe.size() = " << fe.size() << endl;
    //- Return first point
    Info << "fe.first() = " << fe.first() << endl;
    //- Return last point
    Info << "fe.last() = " << fe.last() << endl;
    //- Return n-th point
    Info << "fe.operator[](0) = " << fe.operator[](0) << endl;
    //- Return the points corresponding to this face
    Info << "fe.points(points) = " << fe.points(points) << endl;
    //- Centre point of face
    Info << "fe.centre(points) = " << fe.centre(points) << endl;
    //- Calculate average value at centroid of face
    Info << "fe.average(points, points) = " 
         << fe.average(points, points) << endl;
    //- Magnitude of face area
    Info << "fe.mag(points) = " << fe.mag(points) << endl;
    //- Vector normal; magnitude is equal to area of face
    Info << "fe.normal(points) = " << fe.normal(points) << endl;
    //- Return face with reverse direction
        //  The starting points of the original and reverse face are identical.
    Info << "fe.reverseFace() = " << fe.reverseFace() << endl;
    //- Which vertex on face (face index given a global index)
        //  returns -1 if not found
    Info << "fe.which(1966) = " << fe.which(1966) << endl;
    //- Next vertex on face
    Info << "fe.nextLabel(1) = " << fe.nextLabel(1) << endl;
    //- Previous vertex on face
    Info << "fe.prevLabel(1) = " << fe.prevLabel(1) << endl;
    //- Return number of edges
    Info << "fe.nEdges() = " << fe.nEdges() << endl;
    //- Return edges in face point ordering,
        //  i.e. edges()[0] is edge between [0] and [1]
    Info << "fe.edges() = " << fe.edges() << endl;
    //- Return n-th face edge
    Info << "fe.faceEdge(1) = " << fe.faceEdge(1) << endl;
    //- Return the edge direction on the face
        //  Returns:
        //  -  0: edge not found on the face
        //  - +1: forward (counter-clockwise) on the face
        //  - -1: reverse (clockwise) on the face
    Info << "fe.edgeDirection(fe.faceEdge(1)) = " 
         << fe.edgeDirection(fe.faceEdge(1)) << endl;
    //- compare faces
        //   0: different
        //  +1: identical
        //  -1: same face, but different orientation
    Info << "Foam::face::compare(fe, fe) = " 
         << Foam::face::compare(fe, fe) << endl;
    Info << "----------------------\n" << endl;


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    //class: cell
    Info << "\n-Class: cell----------" << endl;
    cell &cl(cells[100]);
    Info << "cl = " << cl << endl;
    //- Return true if the cell is empty
    Info << "cl.empty() = " << cl.empty() << endl;
    //- Return No. of faces corresponding to this cell
    Info << "cl.size() = " << cl.size() << endl;
    //- Return first face
    Info << "cl.first() = " << cl.first() << endl;
    //- Return last face
    Info << "cl.last() = " << cl.last() << endl;
    //- Return n-th face
    Info << "cl.operator[](0) = " << cl.operator[](0) << endl;
    //- Return number of faces
    Info << "cl.nFaces() = " << cl.nFaces() << endl;
    //- Return labels of cell vertices
    Info << "cl.labels(faces) = " << cl.labels(faces) << endl;
    //- Return the cell vertices
    Info << "cl.points(faces, points) = " << cl.points(faces, points) << endl;
    //- Return cell edges
    Info << "cl.edges(faces) = " << cl.edges(faces) << endl;
    //- Returns cell centre
    Info << "cl.centre(points, faces) = " << cl.centre(points, faces) << endl;
    //- Returns cell volume
    Info << "cl.mag(points, faces) = " << cl.mag(points, faces) << endl;
    Info << "----------------------\n" << endl;

[队长别开枪]

@飞火流星jyj 我贴了一些我自己的测试代码:laughing:

[飞火流星jyj]

@队长别开枪 十分感谢，我回去好好看看:laughing:

[李东岳]

@队长别开枪 在 网格体积 中说：

上半年年度最佳代码 :cheeky:

[zhanghan]

@赵一铭

能指点我一下吗，实在不知道该怎么做

[李东岳]

@zhanghan

volScalarField volume
{
    IOobject
        (
            "volume",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        mesh.V()
}

把上面的代码复制进去你的代码测试一下，我没测试。可以添加在createFields.H中。

[飞火流星jyj]

@李东岳 李老师我测试了下

 Info<< "Reading volume\n" << endl;
volScalarField volume
(
    IOobject
        (
            "volume",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        mesh.V()
);

提示错误是：

 no matching function for call to ‘Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::GeometricField(Foam::IOobject, const Foam::DimensionedField<double, Foam::volMesh>&)’
 );

没有这种类型

[李东岳]

volScalarField volume
(
    IOobject
	(
		"volume",
		runTime.timeName(),
		mesh,
		IOobject::NO_READ,
		IOobject::AUTO_WRITE
	),
	mesh,
	scalar(0.0)
);
volume.internalField() = mesh.V(); 

试试这个

[李东岳]

如果是OpenFOAM-4.x，应该是：

volume.ref() = mesh.V(); 

[飞火流星jyj]

@李东岳 李老师可以了，3.0.x没问题，生成的数据是怎么排列的呢？

[飞火流星jyj]

@队长别开枪 请问您测试过吗，我测试了在createField.H加入

Info << "\n-Class: cell----------" << endl;
    cell &cl(cells[100]);
Info << "cl.centre(points, faces) = " << cl.centre(points, faces) << endl;

报错如下：

createFields.H:33:14: error: ‘cells’ was not declared in this scope
     cell &cl(cells[100]);
              ^
createFields.H:34:52: error: ‘points’ was not declared in this scope
 Info << "cl.centre(points, faces) = " << cl.centre(points, faces) << endl;
                                                    ^
createFields.H:34:60: error: ‘faces’ was not declared in this scope
 Info << "cl.centre(points, faces) = " << cl.centre(points, faces) << endl;
                                                            ^

怎么修改就好了

[队长别开枪]

@飞火流星jyj 我今晚或者明天把整个测试代码发上来吧:happy:

[飞火流星jyj]

@队长别开枪 好的，十分感谢

[程迪]

@飞火流星jyj
别用cl.centre()

这个东西和mesh.C()似乎不一样，OF有两种中心，几何重心和FVM意义上的中心。二者在三维网格有warp的情况下有一些微小的差别。

[飞火流星jyj]

@程迪 好的，一个是polymesh的一个是扩展的FvMesh，那我该用哪个表示出某网格的x，y，z，这该怎么写啊

[程迪]

@飞火流星jyj
GeometricField类型有一个component函数，返回的是分量，参数是direction，实际上类型是char/uint_8，你试试'x','y','z'或者1,2,3就好。

[飞火流星jyj]

@程迪 请问怎么查看是0,1,2还是1,2,3这样的编号:confused:

[赵一铭]

@飞火流星jyj

例如你要输出U的x方向，就是U.x();
不同的版本好像不一样，例如有时候可能是U.component(0);

:big_mouth:

[队长别开枪]

@飞火流星jyj 不好意思，这两天太忙没顾得上，刚刚把代码发上来。链接地址：meshInfoFoam，它会将体单元的体积、中心和面单元的面积、中心、面积矢量、单位矢量全部写入初始文件夹0，具体结果可以参考链接meshInfoFoam output。代码里还测试了一些涉及单个网格元素操作的成员函数。

[zhanghan]

@队长别开枪
您好！我试过您的meshInfoFoam了，非常棒的求解器。
我有一个疑问，请问我怎样才能将导出的faceArea与我的cellCenter一一对应上呢？？？

[程迪]

@zhanghan
不可能对上，两者数目都不一样，一个是面场，一个是体心场。

[队长别开枪]

@zhanghan 你是指获取属于某一体单元的面单元的信息么？

[飞火流星jyj]

@队长别开枪 请问如果在forAll(,celli)循环里，mesh.c()是指循环到的网格的中心嘛，还是所有的，如果不是怎样写循环到的celli的中心坐标:crying: 求大神指教

[队长别开枪]

@飞火流星jyj mesh.C()返回的是所有体单元的中心，要获取特定体单元中心的话可以这样写

forAll(mesh.cells(), cellI)
{
    Info << mesh.C()[cellI] << endl;
    // or
    Info << mesh.C().operator[](cellI) << endl;
}

[飞火流星jyj]

@队长别开枪 十分感谢，对于
多个网格的U[celli]的保存，我想用定义一个速度单位的数组存起来，dimensionedvector 好像只是保存一个值吧，这样的数组我不知道该写好，请大神指教

[鲸落]

@队长别开枪 老师您好，请问要是想输出动网格的体积等信息应该作何修改呢？

[队长别开枪]

@鲸落 哪种动网格？overset，AMR？

[鲸落]

@队长别开枪 抱歉回复晚了，动网格文件的代码是这样的

dynamicFvMesh           dynamicMotionSolverFvMesh;           
solver                  velocityLaplacian;                     
velocityLaplacianCoeffs 
{
  diffusivity          uniform;                              
}
v0                      v0 [ 0 1 -1 0 0 0 0 ]   (0 0 0);

[桎梏]

@队长别开枪

感谢大佬提供的工具。但是用了以后为何总是自动停止呢，0文件夹里什么也没有，

请问是因为我的网格数太多（800w），虚拟机的内存（4g）不够么？

[李东岳]

@桎梏

killed不一定总是内存的问题。但是你这个一定是内存的问题。内存太小。老铁买个笔记本吧。现在16G笔记本都有点小了。要上32G

[桎梏]

@李东岳 图方便，我只带了一个商务本，，，周一去实验室试试
笔记本买不啊，这头的实在是贵