4.5 NekMesh modules

NekMesh is designed to provide a pipeline approach to mesh generation. To do this, we break up tasks into three different types. Each task is called a module and a chain of modules specifies the pipeline.

• Input modules read meshes in a variety of formats;

• Processing modules modify meshes to aide in generation processes;

• Output modules write meshes in a variety of formats.

The figure below depicts how these might be coupled together to form a pipeline:

On the command line, we would define this as:

NekMesh -m process1 -m process2 input.msh output.xml

Process modules can also have parameters passed to them, that can take arguments, or not.

NekMesh -m process1:p1=123:booleanparam input.msh output.xml

To list all available modules use the -l command line argument:

Available classes: 
  Input: dat: 
    Reads Tecplot polyhedron ascii format converted from Star CCM (.dat). 
...

and then to see the options for a particular module, use the -p command line argument:

Options for module detect: 
        vol: Tag identifying surface to process.

4.5.1 Input modules

Input and output modules use file extension names to determine the correct module to use. Not every module is capable of reading high-order information, where it exists. The table below indicates support currently implemented.

Note that you can override the module used on the command line. For example, Semtex session files rarely have extensions. So for a session called pipe-3d we can convert this using the syntax

NekMesh pipe-3d:sem pipe-3d.xml

The NekMesh input module also has an option to re-process all composites. By default only composites are reprocessed. However when extracting a smaller mesh from a larger mesh definition by redefining the composite of volumetric elements you can force the edges and/or faces to be reprocessed removing the definition of any edges and/or faces not required in the smaller mesh. This option is called processall and has the syntax

NekMesh input.xml:xml:processall output.xml

Typically, mesh generators allow physical surfaces and volumes to contain many element types; for example a cube could be constructed from a mixture of hexes and prisms. In Nektar++, a composite can only contain a single element type. Whilst the converter will attempt to preserve the numbering of composites from the original mesh type, sometimes a renumbering will occur when a domain contains many element types. For example, for a domain with the tag 150 containing quadrilaterals and triangles, the Gmsh reader will print a notification along the lines of:

Multiple elements in composite detected; remapped: 
- Tag 150 => 150 (Triangle), 151 (Quadrilateral)

The resulting file therefore has two composites of IDs 150 and 151 respectively, containing the triangular and quadrilateral elements of the original mesh. We note there is one exception to this convention in three-dimensional meshes where a face composite can contain both triangular and quadrilateral elements.

Typically a NekMesh call requires both an input and output module to be called, however, by specifying the output file name or file extension as stdout no output file will be created. This option is typically used for mesh statistics processing or inspecting composite values etc. An example call would be:

NekMesh input.xml stdout

or

NekMesh input.xml name.stdout

4.5.2 Output modules

The following output formats are supported:

Note that for Gmsh, it is highly likely that you will need to experiment with the source code in order to successfully generate meshes since robustness is not guaranteed.

The default for xml and vtk is into binary data which has been converted into base64. If you wish to see an ascii output you need to specify the output module option uncompress. For the uncompressed xml the user can execute:

NekMesh Mesh.msh output.xml:xml:uncompress

If the user wishes to obtain the vtk output in the legacy format, the output module option legacy should be specified by executing:

NekMesh Mesh.xml output.vtk:vtk:legacy

Finally, both the Gmsh and Nektar++ output modules support an order parameter, which allows you to generate a mesh of a uniform polynomial order. This is used in the same manner as the above, so that the command

NekMesh Mesh.msh output.msh:msh:order=7

will generate an order 7 Gmsh mesh. In the rest of these subsections, we discuss the various processing modules available within NekMesh.

It is possible to use FieldConvert to extract a smaller region of a mesh from a larger mesh using the “-r xmin,xmax,ymin,ymax,zmin” range option, e.g.

FieldConvert -r xmin,xmax,ymin,ymax,zmin,zmax bigMesh.xml smallMesh.xml

However this will not provide a composite of he faces on the boundary if they were not part of the original boundary composites of the original bigMesh.xml. These can however be recovered by an output module option chkbndcomp in NekMesh. To do this one should specify

NekMesh smallMesh.xml newSmallMesh.xml:xml:chkbndcomp

this will then add a composite for any face that is on the boundary and not part of an existing boundary condition and put it into a composite with id=9999. This capability is currently only set up for 3D meshes

4.5.2.1 HDF5 format

NekMesh and all solvers within Nektar++ - along with subsequent FieldConvert modules - also support the HDF5 format. This allows for faster loading of geometries and meshes within each solver - and is a significant improvement over the XML format. HDF5 is recommended input format for any larger cases.

Converting from XML to HDF5 is a simple task that only requires the one NekMesh command:

NekMesh XMLMesh.xml HDF5Mesh.nekg

This will create two files HDF5Mesh.xml and HDF5Mesh.nekg which are both needed in the same directory to run the simulation. An additional flag in the session file is required, ensuring it is placed before the expansion list being:

1    <GEOMETRY DIM="3" SPACE="3" HDF5FILE="HDF5Mesh.nekg" />

HDF5 also has the additional advantage of ensuring the mesh and session file are split - which allows for easy ammending of the session file - whilst allowing for use of FieldCovnert modules that require only 1 XML input file - rather than having to concatenate the session and mesh XML files. Solvers and any FieldConvert modules can be run by referencing only the session file after the GEOMETRY tag is included.

4.5.3 Extract surfaces from a mesh

Often one wants to visualise surfaces of a 3D mesh, or extract the values of variables on the surface and visualise them. To support this, NekMesh can extract two-dimensional surfaces which can be visualised using FieldConvert in order to extract the value of a 3D field on a given surface.

As an example, we can extract composite surfaces 2 and 3-5 from a mesh using the extract module:

NekMesh -m extract:surf=2,3-5 Mesh.xml output.xml

If you also wish to have the boundaries of the extracted surface detected add the detectbnd option

NekMesh -m extract:surf=2,3-5:detectbnd Mesh.xml output.xml

which will produce new composites for the extracted boundary.

4.5.4 Negative Jacobian detection and histogram visualisation

To detect elements with negative Jacobian determinant, use the jac module:

NekMesh -m jac Mesh.xml output.xml

To get a detailed list of elements which have negative Jacobians, one may use the list option:

NekMesh -m jac:list Mesh.xml output.xml

To extract the elements whose Jacobian is under a specific value for the purposes of visualisation within the domain, use the extract boolean parameter. The value should be a number below 1.0 and is 0.0 by default. This means if extract is set by user with no value, the output mesh file will only have invalid mesh elements:

NekMesh -m jac:extract=value Mesh.xml MeshWithNegativeElements.xml

To show the histogram of Jacobian on the screen, use the histo parameter.

NekMesh -m jac:histo=value1,value2,value3:quality:detail:histofile=filename Mesh.xml output.xml

The histo takes 3 input values: The maximum value shown on the histogram, the number of positive bins on the histogram and the number of negative bins. Each value is separated with a comma. By default the values are: 1.0,10,1. This means by default the histogram will have 10 positive bins between (0.0, 1.0) with an interval 0.1, and will have 1 negative bins which shows all the invalid elememnts. An example can be found below in figure 4.3 (a).

This command histo also generates a text file called "jacs.txt" which contains the information of each element in the histogram. The name of the file can be changed by using the command histofile. The information in this text file are: Element ID, Jacobian value, element type, if this element is a boundary element, the boundary edge/face ID, vertex coordinates and neighbor element ID will also be included. Note that only elements whose Jacobian is smaller than "value1" will be counted, and they will be outputted in the histogram and text file. If command detail is added, the composite ID and name of each element will also be outputted into the text file. An example is shown in figure 4.3 (c).

The command quality is an option which show the percentage of each bin of the histogram. The results of quality always match the histogram. An example is shown in figure 4.3 (b).

An example with the following command is shown in figure 4.3. The results of quality and text file always correspond to the histogram.

NekMesh -m jac:histo=0.9,3,2:quality:detail Mesh.xml output.xml

To turn off curvature associated with negative jacobians one can try to use the linearise module:

NekMesh -m linerise:invalid  Mesh.xml output.xml

This option will remove all high order curvature on all element types with singular jacobians.

4.5.5 Spherigon patches

Where high-order information is not available (e.g. when using meshes from imaging software), various techniques can be used to apply a smoothing to the high-order element. In NekMesh we use spherigons, a kind of patch used in the computer graphics community used for efficiently smoothing polygon surfaces.

Spherigons work through the use of surface normals, where in this sense ‘surface’ refers to the underlying geometry. If we have either the exact or approximate surface normal at each given vertex, spherigon patches approximate the edges connecting two vertices by arcs of a circle. In NekMesh we can either approximate the surface normals from the linear elements which connect to each vertex (this is done by default), or supply a file which gives the surface normals.

To apply spherigon patches on two connected surfaces 11 and 12 use the following command:

NekMesh -m spherigon:surf=11,12 \ 
    MeshWithStraighEdges.xml MeshWithSpherigons.xml

If the two surfaces "11" and "12" are not connected, or connect at a sharp edge which is C⁰ continuous but not C¹ smooth, use two separate instances of the spherigon module.

NekMesh -m spherigon:surf=11 -m spherigon:surf=12 \ 
    MeshWithStraighEdges.xml MeshWithSpherigons.xml

This is to avoid the approximated surface normals being incorrect at the edge.

If you have a high-resolution mesh of the surfaces 11 and 12 in ply format it can be used to improve the normal definition of the spherigons. Run:

NekMesh -m spherigon:surf=11,12:usenormalfile=Surf_11-12_Mesh.ply \ 
    MeshWithStraighEdges.xml MeshWithSpherigons.xml

This can be useful, for example, when meshing the Leading edge of an airfoil. Starting from a linear mesh (left figure) the spherigon patches curve the surface elements producing leading edge closer to the underlying geometry:

4.5.6 Periodic boundary condition alignment

When using periodic boundary conditions, the order of the elements within the boundary composite determines which element edges are periodic with the corresponding boundary composite.

To facilitate this alignment, NekMesh has a periodic alignment module which attempts to identify pairs of mutually periodic edges. Given two surfaces surf1 and surf2, which for example correspond to the physical surface IDs specified in Gmsh, and an axis which defines the periodicity direction, the following command attempts to reorder the composites:

NekMesh -m peralign:surf1=11:surf2=12:dir=y \ 
    -m peralign:surf1=13:surf2=14:dir=z Mesh.xml Mesh_aligned.xml

Here the surfaces with IDs 11 and 12 will be aligned normal to the y-axis and the surfaces 13 and 14 will be aligned normal to the z-axis.

Note that this command cannot perform magic – it assumes that any given edge or face lying on the surface is periodic with another face on the opposing surface, that there are the same number of elements on both surfaces, and the corresponding edge or face is the same size and shape but translated along the appropriate axis.

When using periodic boundary conditions that are rotationally aligned the following rotational options should be applied:

NekMesh -m peralign:surf1=11:surf2=12:dir=x:rot=PI/6 \ 
              Mesh.xml Mesh_aligned.xml

where rot specifies the rotation angle in radians from surf1 to surf2 about the axis specified by dir (i.e. the “x” axis in this example).

The rotation/translation is assumed to be exact within a relative tolerance. An optional factor, which is used to scale the tolerance, tolfact can also be specified. The default tolerance factor is 4, and it needs to be tolfact ≥ 1. For example:

NekMesh -m peralign:surf1=11:surf2=12:dir=x:rot=PI/6:tolfact=100 \ 
              Mesh.xml Mesh_aligned.xml

In 3D, where prismatic or tetrahedral elements are connected to one or both of the surfaces, additional logic is needed to guarantee connectivity in the XML file. In this case we append the orient parameter:

NekMesh -m peralign:surf1=11:surf2=12:dir=y:orient input.dat output.xml

4.5.7 Boundary layer splitting

Often it is the case that one can generate a coarse boundary layer grid of a mesh. NekMesh has a method for splitting prismatic and hexahedral elements into finer elements based on the work presented in [34] and [35]. You must have a prismatic mesh that is O-type – that is, you can modify the boundary layer without modifying the rest of the mesh.

Given n layers, and a ratio r which defines the relative heights of elements in different layers, the method works by defining a geometric progression of points

in the standard segment [-1,1]. These are then projected into the coarse elements to construct a sequence of increasingly refined elements, as depicted in figure 4.5.

To split a prism boundary layer on surface 11 into 3 layers with a growth rate of 2 and 7 integration points per element use the following command:

NekMesh -m bl:surf=11:layers=3:r=2:nq=7 MeshWithOnePrismLayer.xml \ 
      MeshWith3PrismsLayers.xml

4.5.8 High-order cylinder generation

Generating accurate high-order curved geometries in Gmsh is quite challenging. This module processes an existing linear cylindrical mesh, with axis aligned with the z-coordinate axis, to generate accurate high-order curvature information along the edges.

NekMesh -m cyl:surf=2:r=1.0:N=5 LinearCylinder.xml HighOrderCylinder.xml

The module parameters are:

• surf: Surface on which to apply curvature. This should be the outer surface of the cylinder.

• r: Radius of the cylinder.

• N: Number of high-order points along each element edge.

4.5.9 Linearisation

The ability to remove all the high-order information in a mesh can be useful at times when mesh generation is tricky or impossible in the presence of curvature

To do this in NekMesh use the command:

NekMesh -m linearise:all high-order-mesh.xml linear-mesh.xml

The output will contain only the linear mesh information, all curved information is removed. Alternatively

NekMesh -m linearise:invalid high-order-mesh.xml linear-mesh.xml

attempts to remove curvature from elements only where necessary. This is a simple algorithm that removes curvature from invalid elements and repeats until all elements are valid. Either all or invalid must be specified.

• all: remove curvature from all elements.

• invalid: remove curvature from invalid elements.

• prismonly: consider only prisms when removing curvature. This is useful in the presence of a prismatic boundary layer.

4.5.10 Extracting interface between tetrahedra and prismatic elements

When the mesh is three-dimensional and comprised of a prismatic boundary layer with tetrahedra in the interior of the domain, this module extracts the prismatic elements only, and constructs a boundary region for the interface between the tetrahedra and prisms. This is useful in, for example, the study of aortic flows, where the prismatic boundary layer can be extracted and refined to study unsteady advection-diffusion problems on a more refined grid inside the boundary layer.

To use this module you therefore use the command:

NekMesh -m extracttetprisminterface input.xml output.xml

There are no configuration options for this module, as it is highly specific to a certain class of meshes.

4.5.11 Boundary identification

Some mesh formats lack the ability to identify boundaries of the domain they discretise. NekMesh has a rudimentary boundary identification routine for conformal meshes, which will create a composite of edges (2D) or faces (3D) which are connected to precisely one element. This can be done using the detect module:

NekMesh -m detect volume.xml volumeWithBoundaryComposite.xml

4.5.12 Scalar function curvature

This module imposes curvature on a surface given a scalar function z = f(x,y). For example, if on surface 1 we wish to apply a surface defined by a Gaussian z = exp[-(x² + y²)] using 7 quadrature points in each direction, we may issue the command

NekMesh -m scalar:surf=1:nq=7:scalar=exp\(x*x+y*y\) mesh.xml deformed.xml

4.5.13 Link Checking

It is quite possible that a mesh contains some sort of hanging entity or element connectivity error. The check link module is a fast check that, a) elements are correctly connected and b) the boundary entities (composites) match the interior domain:

NekMesh -m linkcheck mesh.xml mesh2.xml

This module should be added to the module chain if the user suspects there may be a mesh issue. The module will print a warning if there is a connectivity error.

4.5.14 2D mesh extrusion

This module allows a 2D mesh, quads, triangles or both, to be extruded in the z direction to make a simple 3D mesh made of prisms and hexahedra. It is also capable of extruding the high-order curvature within the 2D mesh. The module requires two parameters:

NekMesh -m extrude:layers=n:length=l 2D.xml 3D.xml

length which determines how long the z extrusion will be and layers, the number of elements in the z direction.

4.5.15 2D mesh revolution

This module allows a 2D mesh, quads, triangles or both, to be revolved around the y axis to make a simple 3D mesh made of prisms and hexahedra. It is also capable of revolving the high-order curvature within the 2D mesh. The module requires two parameters:

NekMesh -m revolve:layers=n:angle=a 2D.xml 3D.xml

angle which determines the angle of revolution (θ) around the y axis and layers, the number of elements in the θ direction.

The angle can also be set to "full" for a full revolution of 2π, although this is the default setting:

NekMesh -m revolve:layers=n:angle=full 2D.xml 3D.xml

The angle must be between 0 and 2π. The 2D mesh must not cross or touch the y axis i.e. x > 0 for all points. By default the mesh is revolved through the full 2π.

4.5.16 Variational Optimisation

This module can correct invalid and improve the quality of elements in high-order meshes by applying curvilinear deformation to the interiors of domains. It achieves this by solving a solid mechanics system which, using variational calculus has been cast is a non-linear energy optimsation problem. It is basis of the work in [48].

It works by considering the boundary (curved) mesh entities to be fixed and moving the interior nodes to a lower energy configuration. This new configuration in most scenarios is a higher quality mesh. The energy is evaluated depending on which functional is chosen. We find hyperleastic to be the most reliable but it can also model the mesh and a linearelastic solid as well as functionals based on the Winslow equation and the distortion method proposed by Roca et al. [16].

There are a large number of options which can be viewed using the help function but the basic usage is:

NekMesh -m varopti:type inital.xml optimised.xml

where type can be hyperelastic, linearelastic, winslow or roca.

4.5.17 r-adaptation: interpolation file

This module can deform an existing mesh by using the variational optimiser presented above. A file must be provided that contains a list of points and a scaling value for each of them. This scaling factor is then used to target an element size based on the initial size of the element. Scaling values are interpolated throughout the domain based on the interpolation method of the main library. The file should look like

0 0 0 2.0 
0 1 0 2.0 
1 0 0 0.5 
1 1 0 0.5

where the first three columns are x, y, z and the last column is the scaling factor.

The call is identical to the variational optimisation module above:

NekMesh -m varopti:type:scalingfile=file.txt:subiter=x inital.xml adapted.xml

where subiter is an additional parameter to the variational optimiser that defines the frequency at which individual elements update their target scaling based on their latest location in the domain. subiter should be a scalar and is the number of steps between updates. It is often recommended to run r-adaptation on a linear mesh for stability and performance reasons. Note also that the mesh must have CAD information in order for nodes to slide on curves and surfaces.

4.5.18 r-adaptation: CAD curves

Another option for r-adaptation can be defined where instead of providing an interpolation file, the user can provide a list of CAD curves towards which the scaling takes place. With this option, the r-adaption uses a fixed scaling factor for all elements which are on or within a specified distance from the curve.

The call is similar to the r-adaption with interpolation file, with the difference of replacing the file with a list of CAD curve IDs, the scaling factor and an optional distance from the curve:

NekMesh -m varopti:type:radaptcurves=CurveIDs:radaptscale=scale:radaptrad=rad:subiter=x inital.xml adapted.xml

where radaptcurves, radaptscale and radaptscale are additional containing the list of CAD curve IDs (integers separated by commas, or with a dash for ranges, i.e. 1-3,6 would select curves with ID 1, 2, 3, and 6), the scaling factor as a float and the distance from the CAD curves as a float, respectively. As with the r-adaption with interpolation file, subiter is to define the frequency at which individual elements update their target scaling based on their current proximity to the selected CAD curves.

4.5.19 Mesh projection

This module can take any linear mesh, providing that it is a close representation of the CAD and project the boundary of the mesh onto the CAD. This will curve the surface of the mesh. The method has a number of failsafes ensuring that even bad CAD or poor linear meshes should be able to be curved to some degree. If the method encounters an issue, such as the linear mesh being a large distance from the CAD, it will simply leave that element straight sided. A well made CAD model and accurate linear mesh should be curved with little issue.

The module needs to be informed of the CAD file to project the mesh to and the order at which to curve the surface:

NekMesh -m projectcad:file=cadfile.step:order=x inital.xml optimised.xml

4.5.20 Mesh combine

This module can combine two mesh files into one output file. This will read all information from two mesh files and generate the combine meish file with all these information. The process will duplicate composite IDs from mesh 1 detected in mesh 2 and will remap into the output file.

NekMesh -m combine:file=mesh1.xml mesh2.xml combine.xml