[Gmsh] c++ API - Import mesh vertices and elements
Ryan Durscher
dursch at ufl.edu
Sat Aug 22 06:00:07 CEST 2015
Simone,
Please see below for sample code. I have tried commenting out those line,
but it resulted in the same issue. The error I am seeing is :
"Initializing memorypools.
tetrahedron per block: 8188.
Size of a point: 128 bytes.
Size of a tetrahedron: 96 (96) bytes.
Size of a shellface: 216 (216) bytes.
Error: The point set is trivial."
These are tetgen outputs with the "The point set is trivial" being a
result of the 0 nodes being input into it. Incidentally if I call
refineMesh it properly refines the surface mesh I am inputting into the
gmshModel if I comment out the region stuff. From the GUI in order to get a
volume mesh I have to use the elementary entities -> add -> volume then
select the surface mesh (generated from step *** in the code). This creates
a *.geo file which I cant figure out how to replicate though the API. The
form of the *.geo file is:
Discrete Face(1) = {};
Surface Loop(3) = {1};
Volume(3) = {3};
Physical Surface(1) = {1};
I have been staring at this for a while now and any help would be greatly
appreciated.
Thanks,
Ryan
// Gmsh headers
#include "Gmsh.h"
#include "GModel.h"
#include "GEntity.h"
#include "discreteFace.h"
#include "discreteRegion.h"
#include "MTriangle.h"
#include "MFace.h"
// *verts is a 3*numVerts array of coordinates verts(3*i + 0) = x,
verts(3*i + 1) = y, and verts(3*i + 3) = z
// *conn is a 3*numTris array - define the triangles connectivity
int gmsh_Volume_Mesh(int numVerts,
int numTris,
double *verts,
int *conn)
{
// Initialize variables
int i,j;
printf("USING GMSH\n");
GmshInitialize();
GmshSetOption("Mesh","Algorithm",1.);
GmshSetOption("General","Verbosity", 100.);
//Create Gmsh model and set the factory
GModel *modelGmsh = new GModel;
modelGmsh->setFactory("Gmsh");
// Create a discrete face
GFace *interFace = new discreteFace(modelGmsh,1);
// Pack triangle verts into vertex vector
std::vector<MVertex*> vertex;
vertex.resize((int) numVerts);
for (i = 0 ; i < numVerts ; i ++) {
vertex[i] = new MVertex(verts[3*i + 0],
verts[3*i + 1],
verts[3*i + 2],
interFace,
i+1); // Want index to start at 1
// Add vertex to face
interFace->addMeshVertex(vertex[i]);
}
// Pack connectivity of triangles into face
for (i = 0; i < numTris; i++) {
interFace->addTriangle(new MTriangle(vertex[conn[3*i + 0]-1],
vertex[conn[3*i + 1]-1],
vertex[conn[3*i + 2]-1]));/*,
i+1,2));*/
}
printf("Number of elements in face =
%d\n",interFace->getNumMeshElements());
// Add face to model
modelGmsh->add(interFace);
// ****** modelGmsh->writeGEO("test.geo");
// ****** modelGmsh->writeMSH("test.msh");
// Create face loop - only 1 surface
std::vector<GFace*> interfaceFace; //in my case the surface is given in
one piece only
interfaceFace.push_back(interFace);
std::vector <std::vector<GFace*> > interfaceLoop;
interfaceLoop.push_back(interfaceFace);
// Create a physical region
GRegion *regionInterface = modelGmsh->addVolume(interfaceLoop);
//Add a volume physical entity
regionInterface->addPhysicalEntity(1);
//Add the region to the model
modelGmsh->add(regionInterface);
modelGmsh->writeGEO("test.geo");
modelGmsh->writeMSH("test.msh");
//Mesh the model
//modelGmsh->refineMesh(1); // Refining the surface mesh works fine if
I comment out the region stuff
modelGmsh->mesh(3);
/*printf("Num Regions = %d\n",gmshmodel->getNumRegions());
printf("Num Faces = %d\n",gmshmodel->getNumFaces());
printf("Num Edges = %d\n",gmshmodel->getNumEdges());
printf("Num Verts = %d\n",gmshmodel->getNumVertices());*/
delete modelGmsh;
GmshFinalize();
printf("Done meshing\n");
}
On Fri, Aug 21, 2015 at 11:15 PM, Simone GREMMO [531257] <
Simone.GREMMO at umons.ac.be> wrote:
> Hi Ryan,
> I used the addPhysicalEntity(1) so that when saving the mesh, it will
> write to file the volume mesh and not the surfaces used to bound the volume.
> As far as I could understand the addPhysicalEntity should not interfere
> with the mesh generation, but can used to save in the msh just a part of
> the mesh: so maybe you generate the mesh, but do not write to file. Is it
> possible?
> Have you tried to comment the lines
> regionInterface->addPhysicalEntity(1);
> modelGmsh->add(regionInterface);
> and see if in the resulting msh file you get something different?
>
> Can you copy the relevant lines of your code? It may be easier to
> understand your problem.
>
> Simone
>
>
> ------------------------------
> *De :* rjdurscher at gmail.com [rjdurscher at gmail.com]
> *Envoyé :* samedi 22 août 2015 1.37
> *À :* Gmsh
> *Cc :* gmsh at geuz.org; Simone GREMMO [531257]
> *Objet :* Re: [Gmsh] c++ API - Import mesh vertices and elements
>
> Hi Simone,
>
> I am trying to implement the same thing you describe here in my code with
> Gmsh. I have followed your code below exactly however I am unable to
> generate the volume mesh. Outputting a *.msh file (model->writeMSH()) at
> various stages the problem seems to come about during "
> regionInterface->addPhysicalEntity(1);" in which the previously defined
> surface mesh now disappears (according the *.msh file). Any help would be
> greatly appreciated. I am currently using the latest version of Gmsh
> (2.10.1). Furthermore any additional code you could share regarding your
> second post (meshing between two surface) would be greatly appreciated.
>
> Thanks,
> Ryan
>
>
> On Monday, February 9, 2015 at 6:41:49 AM UTC-5, Simone GREMMO [531257]
> wrote:
>>
>> Dear all,
>>
>> I think I have partially solved my problem. Here is the code that I have
>> written and let me to create the volume mesh (Thanks to Mikhail for his
>> suggestion in this answer:
>> http://geuz.org/pipermail/gmsh/2015/009536.html):
>>
>>
>>
>> GmshInitialize();
>>
>> GmshSetOption("Mesh","Algorithm",1.);
>>
>> //Create Gmsh model and set the factory
>>
>> GModel *modelGmsh = new GModel;
>>
>> modelGmsh->setFactory("Gmsh");
>>
>>
>>
>> //GFace is defined by a mesh and not through GVertex and GEdge
>>
>> GFace *interfaceGmsh = new discreteFace(modelGmsh,1);
>>
>> //Function that loops over mesh vertices and elements and import them
>> in a gmsh structure [see at the bottom of the main function]
>>
>> importSurfaceMeshintoGmsh(interface, interfaceGmsh);
>>
>> modelGmsh->add(interfaceGmsh);
>>
>>
>>
>> //Using the previously defined surface, define a surface loop, used
>> to enclose a volume
>>
>> std::vector<GFace*> interfaceFace; //in my case the surface is given
>> in one piece only
>>
>> interfaceFace.push_back(interfaceGmsh);
>>
>> std::vector <std::vector<GFace*> > interfaceLoop;
>>
>> interfaceLoop.push_back(interfaceFace);
>>
>>
>>
>> //Add the volume enclose by the surface loop to the model and create
>> the related GRegion
>>
>> GRegion *regionInterface = modelGmsh->addVolume(interfaceLoop);
>>
>> //Add a volume physical entity
>>
>> regionInterface->addPhysicalEntity(1);
>>
>> //Add the region to the model
>>
>> modelGmsh->add(regionInterface);
>>
>> //Mesh the model
>>
>> modelGmsh->mesh(3);
>>
>>
>> modelGmsh->writeGEO("/home/simoneg/Hexpress-Hybrid_projects/IncludeGmsh-HH/TESTGMSHAPI_keyword/interfaceVertGmsh.geo");
>>
>>
>> modelGmsh->writeMSH("/home/simoneg/Hexpress-Hybrid_projects/IncludeGmsh-HH/TESTGMSHAPI_keyword/interfaceVertGmsh.msh");
>>
>> delete modelGmsh;
>>
>> GmshFinalize();
>>
>>
>>
>> //Function to import the mesh in a Gmsh structure
>>
>> //inMesh is the input mesh, created with another tool
>>
>> //gmshMesh is the structure containing the inMesh, but that can be used
>> by Gmsh to generate the volume mesh
>>
>> bool importSurfaceMeshintoGmsh(const My_Mesh &inMesh, GFace *gmshMesh)
>>
>> {
>>
>> My_int numNodes = inMesh.numVerticies(); //get the number of vertices
>> from the inMesh
>>
>> My_int numCells = inMesh.numCells(); ////get the number of elements
>> from the inMesh
>>
>>
>>
>> //Get all the vertices from the mesh
>>
>> std::vector<MVertex*> mv;
>>
>> mv.resize((int)numNodes);
>>
>> for (MC_Index nodeCounter = 0; nodeCounter < numNodes; nodeCounter++)
>>
>> {
>>
>> My_Coord3D coordHH = inMesh.vertex(nodeCounter); //a 3 dimensions
>> vector with the vertex coordinates
>>
>> mv[nodeCounter] = new MVertex((double)coordHH[0],
>> (double)coordHH[1], (double)coordHH[2],gmshMesh, nodeCounter+1);
>>
>> gmshMesh->addMeshVertex(mv[nodeCounter]); //add the vertex to the
>> mesh
>>
>> }
>>
>>
>>
>>
>>
>> //Get all the elements from the mesh
>>
>> for (My_int cellCounter = 0; cellCounter < numCells; cellCounter++)
>>
>> {
>>
>> My_int numVertCell = inMesh.numElems(cellCounter); //for element
>> #cellCounter, get the number of vertices
>>
>> switch (numVertCell)
>>
>> {
>>
>> case 3: //triangle
>>
>> gmshMesh->addTriangle(new
>> MTriangle(mv[inMesh.elem(cellCounter,0)], //extract the vertexID for
>> element #cellCounter
>>
>>
>> mv[inMesh.elem(cellCounter,1)],
>>
>> mv[inMesh.elem(cellCounter,2)]));
>>
>> break;
>>
>> case 4: //quadrangle
>>
>> gmshMesh->addQuadrangle(new MQuadrangle(mv[inMesh.
>> elem(cellCounter,0)],
>>
>>
>>
>> mv[inMesh. elem(cellCounter,1)],
>>
>> mv[inMesh. elem(cellCounter,2)],
>>
>> mv[inMesh. elem(cellCounter,3)]));
>>
>> break;
>>
>> default:
>>
>> return false;
>>
>> break;
>>
>> }
>>
>> }
>>
>> return true;
>>
>> }
>>
>>
>>
>> This code generates a mesh in a volume without any hole inside. Does
>> anyone know how to add a hole inside the volume and mesh the region between
>> the two surfaces?
>>
>> The .geo commands looks like:
>>
>> *Merge "interfaceMesh.msh";*
>>
>> *Surface Loop(3) = {1};*
>>
>> *Merge "envelopeMovingMesh.msh";*
>>
>> *Surface Loop(4) = {2};*
>>
>> *Volume(1) = {3,4};*
>>
>> *Physical Volume(1) = {5};*
>>
>>
>>
>> Thanks for your help.
>>
>> Simone
>>
>>
>>
>> Dear all,
>>
>> I have included Gmsh as library in my c++ program and I need to use it to
>> mesh a volume region bounded by a surface mesh: in my program I have
>> defined a surface mesh and I need to import its vertices and elements into
>> a gmsh structure to be able to create a volume mesh starting from this
>> surface.
>>
>> Steps to do so, as far as I can understand, are:
>>
>> - Declare a GModel* that will contain the surface
>>
>> - Add a GFace* into the GModel
>>
>> - Add vertices and cells to the GFace* (get them from the mesh
>> created in my program). Actually define a mesh for this GFace object
>>
>> - Create the volume mesh
>>
>>
>>
>> Is this correct? Does anyone have any suggestion on how to write this
>> portion of code?
>>
>>
>>
>> PS I have already looked at the example in utils/api_demos/mainSimple.cpp
>> and have adapted it to my case, but I would rather avoid reading/writing
>> files to transfer meshes between the my program and gmsh to avoid
>> performance issues.
>>
>>
>>
>> Simone
>>
>
--
Durscher, Ryan J
dursch at ufl.edu
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.geuz.org/pipermail/gmsh/attachments/20150822/a536176b/attachment-0001.html>