GMesh 初识

taotao_2016 2023-07-19 发布于辽宁

展开全文

GMesh是三维有限元网格生成器（自然支持1D，2D网格），内含自有CAD内核。GMesh有四大模块：几何、网格、求解器、后处理。可以使用C、C++，Python、Fortran等语言调用其API。本文结合Tutorials对GMesh做简要梳理。

1.几何

GMesh几何采用BRep表示法，支持自有（build-in）和OCC两种几何内核。GMesh不做几何转化，亦不将其他内核的数据转化为自己的几何数据。各种几何操作都是使用原有的几何内核处理。每种几何都有特定标号（tag），点、线、面、体的tag互不影响。标号从1开始，可以连续或不连续。GMesh读入IGS、step使用OCC 几何内核。GMesh对读入的离散数据（STL）使用“reparametrization” procedure”处理，由离散面片恢复原有几何，重新参数化，再进行网格划分。1.1 build-in 几何内核Build-in几何内核创建基本和BRep一致，遵循Point->Curve->Loop->Surface。所有几何元素都有其唯一标识，用integer类型标识。各个类型的标号直接没有关系，如：可以有Point 1，curve 1和surface 1。在创建完几何之后，必须将几何同步到GMesh，才可以划分网格，使用gmsh::model::geo::synchronize()。因为几何内核和GMesh使用的数据结构是不同步的。只有同步之后，才会创建相应的Gmesh Data Structure。

创建模型：GMesh的几何模型使用 gmsh::model::add('Name')添加，如果此命令没有调用，则默认使用“unNamed”。

创建点：GMesh的几何操作函数在命名gmsh::model::geo下，创建点使用gmsh::model::geo::addPoint(),接受四个参数，前三个参数必须提供，后两个是可选参数。前三个是点的坐标。第四个参数设置此点附近的网格大小，第五个参数是点的编号。返回创建点的编号。int p4 = gmsh::model::geo::addPoint(0, .3, 0, 0.01);

创建直线：使用gsh::model::geo::addLine();接受3个参数，前两个必选参数是点的编号，后一个可选参数是线段编号。返回线段编号。创建圆弧：gmsh::model::geo::addCircleArc(),接受7个参数。前三个参数分别为，起始点、圆心、终止点的编号。第四个参数是圆弧编号。后三个参数标识圆弧所在面，可选参数。

创建Loop：使用gmsh::model::geo::addCurveLoop(),接受2个参数，第一个必选参数，由线段编号组成的向量，负号表示和线段方向相反。第二个是可选参数，loop编号。返回Loop编号。

surface：平面由gmsh::model::geo::addPlaneSurface(),接受2个参数，第一个由Loop编号组成的向量。第二个是面的编号。返回面的编号。

Shell:gmsh::model::geo::addSurfaceLoop(),接受2个参数。第一个整数向量组参数，表面组成shell 的面。第二个参数是shell的编号。

solid：gmsh::model::geo::addVolume();接受2个参数。第一个参数是组成solid的shell编号组。第二个参数的solid的编号。

获取点的坐标：先将集合同步到model，再使用gmsh::model::getValue()获得。如下，获得0维5号点的坐标。gmsh::model::geo::synchronize();std::vectorxyz;gmsh::model::getValue(0, 5, {}, xyz);

平移：gmsh::model::geo::translate(), 接受四个参数。第一个参数整数对组，表明将那个元素平移。后三个参数表示平移距离，dx,dy,dz. 如 gmsh::model::geo::translate({ {0, 5} }, -0.02, 0, 0); 表示将0维5号点，沿x负方向平移0.02。

旋转：gmsh::model::geo::rotate(),接受8个参数。第一个整数对组，表明将那个元素选择。接下来三个参数指明旋转中心。再三个参数指明旋转轴。最后一个参数是旋转角度。如 gmsh::model::geo::rotate({ {0, 5} }, 0, 0.3, 0, 0, 0, 1, -M_PI / 4);将0维5号点，绕着过点{0,0.3,0}的Z轴旋转-PI / 4。

复制：gmsh::model::geo::copy();接受2个参数。第一个参数是整数对组，表明要被复制的元素。第二个参数是复制后对新元素的整数对组。如：将0维3号点，复制。会给复制后的元素自动编号，假如模型中点的最大编号是1，此时ov={0,2}。然后将生成的点沿y轴平移0.05。

std::vectorov;

gmsh::model::geo::copy({ {0, 3} }, ov);

gmsh::model::geo::translate(ov, 0, 0.05, 0);也可以同时copy多个元素，如，复制2维，1号面，11号面。此时的ov包含复制后2个面的对组，因为std::vector是动态数组。再将复制后的2个面沿x方向平移0.12.gmsh::model::geo::copy({ {2, 1}, {2, 11} }, ov);gmsh::model::geo::translate(ov, 0.12, 0, 0);

拉伸：gmsh::model::geo::extrude(),接受5个参数。第一个参数是被拉伸元素的对组表示。第2到4个参数，表示拉伸方向。最后参数表示拉伸得到元素的对组表示。

将外部几何文件加入到Gmesh Model：try { gmsh::merge('../t13_data.stl');}catch (...) { gmsh::logger::write('Could not load STL mesh: bye!'); return;}

1.2 OCC内核OCC的操作与原始的OCC基本一致，在命名空间gmsh::model::occ下，此处不详细说明。

使用OCC导入BREP, STEP or IGES格式文件：gmsh::model::occ::importShapes(); 下式导入step文件，并将最终的几何导入到v中。std::vectorv;try {gmsh::model::occ::importShapes('../t20_data.step', v);//此处可以设置只将solid导入，或者将solid下的点、线、面、shell都导入。}catch (...) {gmsh::logger::write('Could not load STEP file: bye!');gmsh::finalize();return 0;}

几何同步到Model：gmsh::model::occ::synchronize();

以下使用OCC定义六面体，gmsh::model::add('t16');gmsh::model::occ::addBox(0, 0, 0, 1, 1, 1, 1);//前三个参数是起始角点坐标，后三个是在xyz方向上的增量，最后一个是编号。

创建球面：gmsh::model::occ::addSphere(x,y,z,r,编号)//参数依次为，圆心坐标、半径、编号。

其他几何命令：gmsh::model::occ::addCircle(); gmsh::model::occ::addCurveLoop(); gmsh::model::occ::addThruSections();

2. Mesh

2.1unstructured算法选择

MeshAdapt：长边被分割，短边被合并。考虑swap交换，网格更优。

Delaunay：插入点在外接圆圆心处，考虑各向异性（anisotropic Delaunay criterion）。

Frontal-Delaunay：网格质量优

Frontal-Delaunay for Quads

BAMG：四边形

MeshAdapt划分复杂曲面是最稳定的，Frontal-Delaunay划分的单元网格质量好。Delaunay划分大的平面，速度快；处理网格尺寸变化大的速度优于Frontal-Delaunay。当Frontal-Delaunay和Delaunay算法失败时，自动切换为MeshAdapt。平面默认使用“Delaunay”，曲面默认使用“MeshAdapt”。

2.2 网格尺寸的设置方式

设置某些点的尺寸：当命令Mesh.MeshSizeFromPoints and Mesh.MeshSizeExtendFromBoundary使用时，可以设置区域中几个点处的期望网格大小，其他区域的网格大小，通过这几个位置插值得到。

通过曲率控制：Mesh.MeshSizeFromCurvature使用时，设置每2π弧度上的网格个数。

通过场控制：

box field，在Box区域内是设置的网格参数。ANSYS有类似功能。

Distance field，距某一特征特定位置区域，使用设置的网格大小。

MathEval，使用数学方程。

PostView，与已经存在的网格的计算结果作为Background size控制尺寸。

编辑一个Scalar Triangle作为Background Mesh，存入.pos文件：

View 'background mesh' {ST(0,0,0,10,0,0,5,5,0){0.5,0.1,0.01};ST(10,0,0,10,10,0,5,5,0){0.1,1.0,0.01};ST(10,10,0,0,10,0,5,5,0){1,0.2,0.01};ST(0,10,0,0,0,0,5,5,0){0.2,0.5,0.01};};

主要命令说明：

gmsh::model::geo::mesh::setSize()，设置某些点处的网格大小。第一个参数是点的对组表示，第二个参数是网格大小。如：将0维编号103点处的网格大小设置为0.1；gmsh::model::geo::mesh::setSize({ {0, 103}},0.1);

gmsh::option::setNumber('Mesh.MeshSizeFactor', 0.1);将网格大小设置为当前定义的0.1倍。

gmsh::model::geo::mesh::setTransfiniteCurve(2, 20); 设置2号线上，20个节点，包含端点。还可以设置节点分布的方式，如，在3号曲线上布置30个节点，线段长度按1.2渐变。相邻2段直线长度之比为1.2，故点越来越密。

gmsh::model::geo::mesh::setTransfiniteCurve(3, 30, 'Progression', 1.2);

gmsh::model::geo::mesh::setRecombine(2, 1);将2D平面1划分为四边形网格。将三角网格合并为四边形，仅2维可用。阀值默认45度。

合并所有2D网格：gmsh::option::setNumber('Mesh.RecombineAll', 1);

gmsh::option::setNumber('Mesh.Smoothing', 100);光滑处理。

使用数学函数作为网格大小控制：要使用setAsBackgroundMesh才能生效

gmsh::model::mesh::field::add('MathEval', 1);

gmsh::model::mesh::field::setString(1, 'F', '0.01*(1.0+30.*(y-x*x)*(y-x*x) + (1-x)*(1-x))');

gmsh::model::mesh::field::setAsBackgroundMesh(1); //参数是field的编号

修改point周围的网格尺寸：gmsh::model::geo::mesh::setSize({ {0, 1}, {0, 2}, {0, 3}, {0, 4} }, lc);将点1,2,3,4处的网格尺寸修改为lc.

假如硬点、硬线到面、体上。如：将line 5 嵌入到 face 1上，在这条直线上会生成带点。

gmsh::model::geo::addLine(1, 3, 5);

gmsh::model::geo::synchronize();//嵌入之前，先同步

gmsh::model::mesh::embed(1, { 5 }, 2, 1);

设置最大最小网格控制：

gmsh::option::setNumber('Mesh.MeshSizeMin', 0.001);

gmsh::option::setNumber('Mesh.MeshSizeMax', 0.3);

设置按曲率控制网格大小：

gmsh::option::setNumber('Mesh.MeshSizeFromCurvature', 20);//(here with a target of 20 elements per 2*Pi radians。

3. 文件格式

此处仅介绍GMesh中最主要的两种格式，MSH文件存储几何、网格信息。pos文件存储计算结果。

3.1 MSH文件格式

msh文件按模块存储，模块名称以“$”标识，如，表示是节点信息。

$Nodes内容$EndNodes

一个msh文件至少要包含三个模块：$MeshFormat、$Nodes 、$Elements。注释使用$Comments/$EndComments。编号从1开始，0作为系统保留。除此之外，还有$PhysicalNames，$Entities，$PartitionedEntities，$Periodic，$NodeData等。本处仅介绍$MeshFormat，$Entities，$Nodes和$Elements四个模块。以Tutorial 1的四边形为例，介绍.msh的文件格式。

$MeshFormat// 指定文件的格式，目前是4.14.1 0 8 //版本；0 for ASCII mode, 1 for binary mode$EndMeshFormat$PhysicalNames1//name的个数2 6 'My surface'//Physical的维数；编号；名称$EndPhysicalNames

$Entities //存储几何信息4 4 1 0 //点的个数；线的个数；面的个数；体的个数1 0 0 0 0 //点的编号；x坐标；y坐标；z坐标；PhysicalGroup的数量（没有Group，0标记）2 0.1 0 0 0 3 0.1 0.3 0 0 4 0 0.3 0 0 1 0 0 0 0.1 0 0 1 5 2 1 -2 //线的编号；minX；minY；minZ；maxX；maxY；maxZ；PhysicalGroup的数量（不在任何Group，0标记);Group编号；点的数量；点的编号2 0.1 0 0 0.1 0.3 0 1 5 2 3 -2 3 0 0.3 0 0.1 0.3 0 0 2 3 -4 // 3号线不在任何Group，故标识为0 04 0 0 0 0 0.3 0 1 5 2 4 -1 // 1,2,4号线在Group 5中，1 0 0 0 0.1 0.3 0 1 6 4 4 1 -2 3 //面的编号；minX；minY；minZ；maxX；maxY；maxZ；PhysicalGroup的数量；Group编号；边界线的数量；线的编号$EndEntities

$Nodes8 10 1 10//块个数；节点个数；最小编号；最大编号0 1 0 1//节点依附的Entity维数（节点1-4都在Point的位置，维数为0）；EntityTag（对应于$Entities的信息）；parametric；此块中的节点数量10 0 00 2 0 120.1 0 00 3 0 130.1 0.3 00 4 0 140 0.3 01 1 0 0 //line 1，3上没有节点，但是line1 在Group 5中，所以line1导出，line3 不导出。1 2 0 2 //节点5-8在直线上，所以维数是1，此Block中包含两个节点，先列出节点号，再列出坐标。560.1 0.2000000000002564 00.1 0.1000000000002544 01 4 0 2780 0.1757875301054114 00 0.07760584538134013 02 1 0 2 //节点9,10在面上，维数是2，面对应$Entities中的1号面。9100.05000000000000001 0.03880292269067007 00.05 0.2439468825264169 0

$EndNodes

$Elements4 17 1 17//块个数（以下共存在四个Block单元组）；单元总数；最小单元编号；最大单元编号1 1 1 1 //单元维数；EntityTag；单元类型（P347-P348）；此块中的单元数量1 1 2 1 2 1 32 3 5 3 5 6 4 6 2 1 4 1 35 4 7 6 7 8 7 8 1 2 1 2 108 8 6 7 9 6 5 7 10 4 7 10 11 1 2 9 12 6 8 9 13 7 5 10 14 2 6 9 15 3 4 10 16 5 3 10 17 8 1 9 $EndElements

$Entities的示意图：