6. Examples

This chapter contains example CEETRON Mesh programs demonstrating the use of the mesh generation modules. The source code for all examples resides in CeetronSAM/src/sam/vis/exam.

6.1. Example 41, Generate a Surface Mapped Mesh

This example illustrates mapped mesh generation using the MapMesh module. A simple beam with a length of 10. and a height of 1. is to be modeled with 8 node Serendipity quadrilateral elements in 2D. The beam geometry is input as 4 points forming a simple linear quadrilateral patch. A node association is defined along edge 4 of the patch so that the nodes generated along the edge may be identified, for example, to apply a clamped boundary restraint. The value chosen for the association, in this case 100, is arbitrary. A node association of 200 is assigned to point 2. This association could be possibly used to apply a concentrated load.

The msh_MapMeshSetPatchParami() function is used to specify that a 8 by 2 grid of quadrilateral elements is to be generated. The msh_MapMeshSetParami() function is used to specify parabolic elements. The nodes and finite elements are actually generated by calling msh_MapMeshGenerate(). The generated nodes and elements are entered into the Connect object which is passed in the argument list. Lastly, the generated nodes and elements are printed.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[4][3] = {{0., 0., 0.}, {10., 0., 0.}, {10., 1., 0.}, {0., 1., 0.}};

/*----------------------------------------------------------------------
                      Generate a Surface Mapped Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_MapMesh* mapmesh;
    vis_Connect* connect;

    Vint i;
    Vint pix[4];
    Vint numnp, numel;
    Vdouble x[3];
    Vint aid, aid1;
    Vint nix, ix[8];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create MapMesh object */
    mapmesh = msh_MapMeshBegin();
    msh_MapMeshDef(mapmesh, 4, 1);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define points */
    for (i = 0; i < 4; i++) {
        msh_MapMeshSetPoint(mapmesh, i + 1, coords[i]);
    }
    /* define patch connectivity */
    pix[0] = 1;
    pix[1] = 2;
    pix[2] = 3;
    pix[3] = 4;
    msh_MapMeshSetPatch(mapmesh, 1, VIS_SHAPEQUAD, 2, 0, 0, pix);

    /* set assoc of 100 at edge 4 */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 1, SYS_EDGE, 4, 100);

    /* set assoc of 200 at point 2 */
    msh_MapMeshSetPointAssoc(mapmesh, VIS_MISCID1, 2, 200);

    /* set mesh parameters */
    msh_MapMeshSetParami(mapmesh, VIS_MESH_MAXI, 3);
    msh_MapMeshSetPatchParami(mapmesh, 1, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMI, 8);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMJ, 2);

    /* generate */
    msh_MapMeshGenerate(mapmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        vis_ConnectNodeAssoc(connect, VIS_MISCID, 1, &i, &aid);
        vis_ConnectNodeAssoc(connect, VIS_MISCID1, 1, &i, &aid1);
        printf("id= %d  x= %f, y= %f, z= %f,  aid= %d, aid1= %d\n", i, x[0], x[1], x[2], aid, aid1);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5], ix[6], ix[7]);
    }
    /* end objects */
    vis_ConnectEnd(connect);
    msh_MapMeshEnd(mapmesh);
    return 0;
}

6.2. Example 41a, Generate a Multiple Patch Mapped Mesh

This example illustrates mapped mesh generation of quadrilateral shell elements using the MapMesh module with multiple patches. A shell structure with a circular section connected to a straight section is extruded in the z direction. The circular and straight sections are each modelled with a separate patch. The circular section at z=0. in the x-y plane is modelled with points 1 through 4 which describe a cubic approximation to the circle. Point 1 is located at the origin. The straight section connects point 4 on the circular section to point 5 which marks the other endpoint of the straight section. Points 6 through 10 describe the cross section at z=50. The first patch describes the circular shell, the second patch describes the straight shell. The two patches share points 4 and 9 where they are connected. Separate node associations are defined at the generated nodes on the z=0. and z=50. cross sections for boundary condition application or any other use. A 4 by 20 element mesh is generated for the first patch and a 10 by 20 element mesh for the second patch. The number of elements generated along the common J direction on each patch must be the same for all the generated nodes to be shared correctly between the two patches. The nodes and finite elements are generated by calling msh_MapMeshGenerate() and are entered into the Connect object which is passed in the argument list. Lastly, the generated nodes and elements are printed.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/*----------------------------------------------------------------------
                      Generate a Multiple Patch Mapped Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_MapMesh* mapmesh;
    vis_Connect* connect;

    Vdouble radius, angle, length, span, radian;
    Vdouble coords[10][3];
    Vint i;
    Vint pix[8];
    Vint numnp, numel;
    Vdouble x[3];
    Vint aid, aid1, pid;
    Vint nix, ix[4];

    radius = 3.;
    angle = 60.;
    length = 10.;
    span = 50.;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create MapMesh object */
    mapmesh = msh_MapMeshBegin();
    msh_MapMeshDef(mapmesh, 10, 2);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectDef(connect, 0, 0);

    /* compute points */
    /* cross section radius, compute 4 points */
    for (i = 0; i < 4; i++) {
        radian = i * angle / 3. * (3.14159 / 180.);
        coords[i][0] = radius - radius * cos(radian);
        coords[i][1] = radius * sin(radian);
        coords[i][2] = 0.;
    }
    /* cross section straight length, 1 additional point */
    coords[4][0] = coords[3][0] + length * sin(radian);
    coords[4][1] = coords[3][1] + length * cos(radian);
    coords[4][2] = 0.;
    /* span, 5 more points offset in z from first 5 points */
    for (i = 0; i < 5; i++) {
        coords[i + 5][0] = coords[i][0];
        coords[i + 5][1] = coords[i][1];
        coords[i + 5][2] = coords[i][2] + span;
    }
    /* define points */
    for (i = 0; i < 10; i++) {
        msh_MapMeshSetPoint(mapmesh, i + 1, coords[i]);
    }
    /* define patch point connectivity */
    /* patch 1 */
    pix[0] = 1;
    pix[1] = 2;
    pix[2] = 3;
    pix[3] = 4;
    pix[4] = 6;
    pix[5] = 7;
    pix[6] = 8;
    pix[7] = 9;
    msh_MapMeshSetPatch(mapmesh, 1, VIS_SHAPEQUAD, 4, 2, 0, pix);
    /* patch 2 */
    pix[0] = 4;
    pix[1] = 5;
    pix[2] = 9;
    pix[3] = 10;
    msh_MapMeshSetPatch(mapmesh, 2, VIS_SHAPEQUAD, 2, 2, 0, pix);

    /* set assoc of 1 at z=0. edge */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 1, SYS_EDGE, 1, 1);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 2, SYS_EDGE, 1, 1);

    /* set assoc of 2 at z=span edge */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 1, SYS_EDGE, 3, 2);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 2, SYS_EDGE, 3, 2);

    /* set property ids as element associations */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_PROPID, 1, SYS_ELEM, 0, 10);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_PROPID, 2, SYS_ELEM, 0, 20);

    /* set mesh parameters */
    msh_MapMeshSetParami(mapmesh, VIS_MESH_MAXI, 2);
    /* patch 1 */
    msh_MapMeshSetPatchParami(mapmesh, 1, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMI, 4);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMJ, 20);
    /* patch 2 */
    msh_MapMeshSetPatchParami(mapmesh, 2, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_MapMeshSetPatchParami(mapmesh, 2, MAPMESH_NUMELEMI, 10);
    msh_MapMeshSetPatchParami(mapmesh, 2, MAPMESH_NUMELEMJ, 20);

    /* generate */
    msh_MapMeshGenerate(mapmesh, connect);

    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam41a.bdf");

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        vis_ConnectNodeAssoc(connect, VIS_MISCID, 1, &i, &aid);
        vis_ConnectNodeAssoc(connect, VIS_MISCID1, 1, &i, &aid1);
        printf("id= %d  x= %f, y= %f, z= %f,  aid= %d, aid1= %d\n", i, x[0], x[1], x[2], aid, aid1);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        vis_ConnectElemAssoc(connect, VIS_PROPID, 1, &i, &pid);
        printf("id= %d  ix= %d %d %d %d,  pid= %d\n", i, ix[0], ix[1], ix[2], ix[3], pid);
    }

    /* end objects */
    vis_ConnectEnd(connect);
    msh_MapMeshEnd(mapmesh);
    return 0;
}

6.3. Example 42, Generate an Extruded Mesh

This example illustrates extruded mesh generation using the ExtMesh module. A cross section of 2 parabolic Serendipity triangles lying in the x-z plane is extruded about the z axis through 45 degrees using a cylindrical coordinate system. The resulting mesh contains 15 node parabolic Serendipity wedge elements.

The cross section geometry is input as 9 nodes in the x-z global coordinate plane connected by 2 parabolic triangular elements. The node coordinates are defined using msh_ExtMeshSetNode() and are always expressed in the global coordinate system. A CoordSys object is used to define a cylindrical coordinate system in which the points defining the extrusion path are defined. The two path endpoints are defined using msh_ExtMeshSetPoint(). The point coordinates are entered in the cylindrical coordinate system.

The msh_ExtMeshSetPathParami() function is used to specify that 4 planes of elements are to be extruded. Since the cross section contains 2 elements, the total number of elements generated will be 8. The function msh_ExtMeshGenerate() generates the nodes and elements. During the mesh generation process, the input node coordinates are transformed to the cylindrical coordinate system defined by CoordSys. All nodes are then extruded along the specified path in the cylindrical coordinate system. The generated nodes are then transformed back to the global coordinate system before they are entered into the Connect object which is passed in the argument list. Lastly, the generated nodes and elements are printed.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble xnode[9][3] = {{2., 0., 0.}, {3., 0., 0.}, {4., 0., 0.}, {2., 0., 1.}, {3., 0., 1.},
                              {4., 0., 1.}, {2., 0., 2.}, {3., 0., 2.}, {4., 0., 2.}};

static Vint ixtri[2][6] = {{1, 3, 7, 2, 5, 4}, {3, 9, 7, 6, 8, 5}};

static Vdouble origin[3] = {0., 0., 0.};
static Vdouble triad[3][3] = {{1., 0., 0.}, {0., 1., 0.}, {0., 0., 1.}};

/*----------------------------------------------------------------------
                      Generate an Extruded Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_ExtMesh* extmesh;
    vis_Connect* connect;

    Vint i;
    Vdouble xp[3];
    Vint pix[2];
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[15];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create ExtMesh object */
    extmesh = msh_ExtMeshBegin();
    msh_ExtMeshDef(extmesh, 9, 2);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectDef(connect, 0, 0);

    /* set coordinate system transformation */
    msh_ExtMeshSetOriginTriad(extmesh, SYS_CYLINDRICAL, origin, triad);

    /* define nodes */
    for (i = 0; i < 9; i++) {
        msh_ExtMeshSetNode(extmesh, i + 1, xnode[i]);
    }
    /* define parabolic triangle connectivity */
    for (i = 0; i < 2; i++) {
        msh_ExtMeshSetElem(extmesh, i + 1, VIS_SHAPETRI, 3, 0, ixtri[i]);
    }
    /* define path end points */
    xp[0] = 0.;
    xp[1] = 0.;
    xp[2] = 0.;
    msh_ExtMeshSetPoint(extmesh, 1, xp);
    xp[1] = 45.;
    msh_ExtMeshSetPoint(extmesh, 2, xp);
    pix[0] = 1;
    pix[1] = 2;
    msh_ExtMeshSetPath(extmesh, 2, pix);

    /* set mesh parameters */
    msh_ExtMeshSetPathParami(extmesh, 1, VIS_MESH_NUMELEMENTS, 4);

    /* generate */
    msh_ExtMeshGenerate(extmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %e, y= %e, z= %e\n", i, x[0], x[1], x[2]);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
        printf("           %d %d %d %d %d %d\n", ix[6], ix[7], ix[8], ix[9], ix[10], ix[11]);
        printf("           %d %d %d\n", ix[12], ix[13], ix[14]);
    }
    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam42.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_ExtMeshEnd(extmesh);
    return 0;
}

6.4. Example 42a, Generate an Extruded Mesh along Normals

This example illustrates extruded mesh generation along node normals using the ExtMesh module. This type of extrusion meshing is useful in creating boundary layer meshes. The original surface mesh to be extruded consists of 10 triangles which form 5 sides of a unit cube. The extrusion directions are defined for each node. The extrusion directions are normalized and scaled to the desired thickness of the boundary layer mesh, .25 in this case. The function msh_ExtMeshSetNodeNormal() is used to input the extrusion directions for each node.

The msh_ExtMeshSetParami() function is used to activate extrusion along node normals, specify the number of layers (4), flag the fact that the normals are in the opposite direction to the triangle connectivity sense and enable the assignment of the VIS_GEOFACE association to the element faces and nodes generated on the initial and final faces of the extruded mesh. Since the normals are opposite to the connectivity sense, the connectivity of the resulting wedge elements must be reversed to yield a positive volume.

The function msh_ExtMeshGenerate() generates the nodes and elements. Lastly, the generated nodes and elements are printed.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/*----------------------------------------------------------------------
                      Generate a Multiple Patch Mapped Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_MapMesh* mapmesh;
    vis_Connect* connect;

    Vdouble radius, angle, length, span, radian;
    Vdouble coords[10][3];
    Vint i;
    Vint pix[8];
    Vint numnp, numel;
    Vdouble x[3];
    Vint aid, aid1, pid;
    Vint nix, ix[4];

    radius = 3.;
    angle = 60.;
    length = 10.;
    span = 50.;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create MapMesh object */
    mapmesh = msh_MapMeshBegin();
    msh_MapMeshDef(mapmesh, 10, 2);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectDef(connect, 0, 0);

    /* compute points */
    /* cross section radius, compute 4 points */
    for (i = 0; i < 4; i++) {
        radian = i * angle / 3. * (3.14159 / 180.);
        coords[i][0] = radius - radius * cos(radian);
        coords[i][1] = radius * sin(radian);
        coords[i][2] = 0.;
    }
    /* cross section straight length, 1 additional point */
    coords[4][0] = coords[3][0] + length * sin(radian);
    coords[4][1] = coords[3][1] + length * cos(radian);
    coords[4][2] = 0.;
    /* span, 5 more points offset in z from first 5 points */
    for (i = 0; i < 5; i++) {
        coords[i + 5][0] = coords[i][0];
        coords[i + 5][1] = coords[i][1];
        coords[i + 5][2] = coords[i][2] + span;
    }
    /* define points */
    for (i = 0; i < 10; i++) {
        msh_MapMeshSetPoint(mapmesh, i + 1, coords[i]);
    }
    /* define patch point connectivity */
    /* patch 1 */
    pix[0] = 1;
    pix[1] = 2;
    pix[2] = 3;
    pix[3] = 4;
    pix[4] = 6;
    pix[5] = 7;
    pix[6] = 8;
    pix[7] = 9;
    msh_MapMeshSetPatch(mapmesh, 1, VIS_SHAPEQUAD, 4, 2, 0, pix);
    /* patch 2 */
    pix[0] = 4;
    pix[1] = 5;
    pix[2] = 9;
    pix[3] = 10;
    msh_MapMeshSetPatch(mapmesh, 2, VIS_SHAPEQUAD, 2, 2, 0, pix);

    /* set assoc of 1 at z=0. edge */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 1, SYS_EDGE, 1, 1);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 2, SYS_EDGE, 1, 1);

    /* set assoc of 2 at z=span edge */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 1, SYS_EDGE, 3, 2);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_MISCID, 2, SYS_EDGE, 3, 2);

    /* set property ids as element associations */
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_PROPID, 1, SYS_ELEM, 0, 10);
    msh_MapMeshSetPatchAssoc(mapmesh, VIS_PROPID, 2, SYS_ELEM, 0, 20);

    /* set mesh parameters */
    msh_MapMeshSetParami(mapmesh, VIS_MESH_MAXI, 2);
    /* patch 1 */
    msh_MapMeshSetPatchParami(mapmesh, 1, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMI, 4);
    msh_MapMeshSetPatchParami(mapmesh, 1, MAPMESH_NUMELEMJ, 20);
    /* patch 2 */
    msh_MapMeshSetPatchParami(mapmesh, 2, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_MapMeshSetPatchParami(mapmesh, 2, MAPMESH_NUMELEMI, 10);
    msh_MapMeshSetPatchParami(mapmesh, 2, MAPMESH_NUMELEMJ, 20);

    /* generate */
    msh_MapMeshGenerate(mapmesh, connect);

    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam41a.bdf");

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        vis_ConnectNodeAssoc(connect, VIS_MISCID, 1, &i, &aid);
        vis_ConnectNodeAssoc(connect, VIS_MISCID1, 1, &i, &aid1);
        printf("id= %d  x= %f, y= %f, z= %f,  aid= %d, aid1= %d\n", i, x[0], x[1], x[2], aid, aid1);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        vis_ConnectElemAssoc(connect, VIS_PROPID, 1, &i, &pid);
        printf("id= %d  ix= %d %d %d %d,  pid= %d\n", i, ix[0], ix[1], ix[2], ix[3], pid);
    }

    /* end objects */
    vis_ConnectEnd(connect);
    msh_MapMeshEnd(mapmesh);
    return 0;
}

6.5. Example 42b, Generate a TetMesh with an Extruded Boundary Layer

This example illustrates generating a boundary layer mesh from a surface geometry and filling the interior with tetrahedra. The process outlined in this example is typical of a complete mesh generation sequence for meshing the interior of a surface mesh in which a boundary layer mesh is required. The SurfMesh, ExtMesh and TetMesh objects are using in succession to generate the complete mesh.

A simple unit cube defined by 12 boundary triangles defines the surface geometry. These triangles are entered into the SurfMesh object for generating the initial surface mesh of triangles. The resulting triangles are entered into the ExtMesh object. The normals of the surface triangle mesh are averaged at each node and scaled by the desired boundary layer length of .25. Since the boundary layer mesh is extruded in the opposite direction of the sense of the boundary triangle connectivity, the parameter EXTMESH_REVERSE is set to create positive volume boundary layer elements. The extrusion distance is large enough that the node vectors begin to “intersect” resulting in inverted elements. The parameter EXTMESH_LOCALCHECK is enabled to check for inversion and recede the extrusion distance sufficiently to avoid inversion. The parameter EXTMESH_FACEMARK is enabled to mark all nodes and element face on the final extrusion layer with the association VIS_GEOFACE. All nodes and element faces on the initial extrusion layer are marked with a value of 1 and with a 2 on the final extrusion layer. This will facilitate passing the triangles on the final extrusion layer to the TetMesh module for the final tetrahedral meshing of the interior.

Finally, using the VIS_GEOFACE association, the nodes and triangles on the final extrusion layer are detected and passed to the TetMesh module. The VIS_MESH_INWARD parameter must be set since the sense of the triangle on the final extrusion layer will point to the inside of the geometry.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/* a simple unit cube */
static Vdouble xc[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                           {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};

static Vint tefl[12][3] = {{1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1}, {1, 0, 1}, {0, 1, 1},
                           {0, 1, 1}, {1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {0, 1, 1}, {1, 1, 0}};

/* surface triangle node normals */
static Vdouble norms[12][3][3] = {{{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}}, {{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}},
                                  {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},    {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},
                                  {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}}, {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}},
                                  {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},    {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},
                                  {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}}, {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}},
                                  {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}},    {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}}};

/*----------------------------------------------------------------------
         Cross product between two vectors
----------------------------------------------------------------------*/
static void
cross(Vdouble a[3], Vdouble b[3], Vdouble c[3])
{
    c[0] = a[1] * b[2] - a[2] * b[1];
    c[1] = a[2] * b[0] - a[0] * b[2];
    c[2] = a[0] * b[1] - a[1] * b[0];
}

/*----------------------------------------------------------------------
         Difference between two vectors
----------------------------------------------------------------------*/
static void
diff(Vdouble a[3], Vdouble b[3], Vdouble c[3])
{
    c[0] = a[0] - b[0];
    c[1] = a[1] - b[1];
    c[2] = a[2] - b[2];
}

/*----------------------------------------------------------------------
         Generate Extruded Boundary Layer with Interior TetMesh
----------------------------------------------------------------------*/
int
main()
{
    vis_IdTran* idtran;
    msh_SurfMesh* surfmesh;
    msh_TetMesh* tetmesh;
    msh_ExtMesh* extmesh;
    vis_Connect *connectsurf, *connectext, *connecttet, *connectall;
    vis_Group *ngroup, *fgroup;
    Vdouble(*normal)[3], x[3][3], u[3], v[3], w[3], mag;
    Vint i, j, nix, ix[10], numnpsurf, numelsurf, numnpext, numelext;
    Vint jx[3], numnp, numel, no, aid, num;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* perform surface meshing */
    /* create SurfMesh object */
    surfmesh = msh_SurfMeshBegin();
    msh_SurfMeshDef(surfmesh, 8, 12);

    /* define points */
    for (i = 0; i < 8; i++) {
        msh_SurfMeshSetPoint(surfmesh, i + 1, xc[i], 1);
    }
    /* define triangle connectivity */
    for (i = 0; i < 12; i++) {
        msh_SurfMeshSetTri(surfmesh, i + 1, tris[i], tefl[i]);
        msh_SurfMeshSetTriNorm(surfmesh, i + 1, norms[i]);
    }
    /* set mesh parameters */
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, 0.2);
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_MAXI, 2);

    /* create Connect object */
    connectsurf = vis_ConnectBegin();
    vis_ConnectPre(connectsurf, SYS_DOUBLE);
    vis_ConnectDef(connectsurf, 0, 0);
    /* generate */
    msh_SurfMeshGenerate(surfmesh, connectsurf);
    if (msh_SurfMeshError(surfmesh)) {
        printf("surf mesh generation error\n");
        exit(1);
    }
    else {
        printf("surf mesh generation complete\n");
    }

    /* write generated surf mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connectsurf, SYS_NASTRAN_BULKDATA, "exam42b_srf.bdf");

    /* perform extrusion meshing */
    /* compute averaged nodal normals */
    vis_ConnectNumber(connectsurf, SYS_ELEM, &numelsurf);
    vis_ConnectNumber(connectsurf, SYS_NODE, &numnpsurf);
    normal = (Vdouble(*)[3])malloc(3 * numnpsurf * sizeof(Vdouble));
    for (i = 0; i < numnpsurf; ++i) {
        normal[i][0] = normal[i][1] = normal[i][2] = 0.;
    }
    for (i = 0; i < numelsurf; ++i) {
        vis_ConnectElemNode(connectsurf, i + 1, &nix, ix);
        vis_ConnectCoordsdv(connectsurf, nix, ix, x);
        diff(x[1], x[0], u);
        diff(x[2], x[0], v);
        cross(u, v, w);
        mag = sqrt(w[0] * w[0] + w[1] * w[1] + w[2] * w[2]);
        w[0] /= mag;
        w[1] /= mag;
        w[2] /= mag;
        for (j = 0; j < 3; ++j) {
            normal[ix[j] - 1][0] += w[0];
            normal[ix[j] - 1][1] += w[1];
            normal[ix[j] - 1][2] += w[2];
        }
    }
    /* set normal length to layer size of 0.25 */
    for (i = 0; i < numnpsurf; ++i) {
        mag = sqrt(normal[i][0] * normal[i][0] + normal[i][1] * normal[i][1] + normal[i][2] * normal[i][2]);
        mag = 0.25 / mag;
        normal[i][0] *= -mag;
        normal[i][1] *= -mag;
        normal[i][2] *= -mag;
    }

    /* create ExtMesh object */
    extmesh = msh_ExtMeshBegin();
    msh_ExtMeshDef(extmesh, numnpsurf, numelsurf);

    /* enter generated surface mesh */
    for (i = 1; i <= numnpsurf; ++i) {
        vis_ConnectCoordsdv(connectsurf, 1, &i, &u);
        msh_ExtMeshSetNode(extmesh, i, u);
        msh_ExtMeshSetNodeNormal(extmesh, i, normal[i - 1]);
    }
    for (i = 1; i <= numelsurf; ++i) {
        vis_ConnectElemNode(connectsurf, i, &nix, ix);
        msh_ExtMeshSetElem(extmesh, i, SYS_SHAPETRI, 2, 0, ix);
    }
    /* set mesh parameters */
    msh_ExtMeshSetParami(extmesh, EXTMESH_NORMAL, SYS_ON);
    msh_ExtMeshSetParamd(extmesh, VIS_MESH_GROWTHRATE, 2.);
    msh_ExtMeshSetParami(extmesh, EXTMESH_LOCALRECEDE, SYS_ON);
    msh_ExtMeshSetParami(extmesh, EXTMESH_REVERSE, SYS_ON);
    msh_ExtMeshSetNormalPath(extmesh, 4, 1., SYS_OFF);
    /* mark faces on original layer with value of 1
       and faces on final layer with value of 2 */
    msh_ExtMeshSetParami(extmesh, EXTMESH_FACEMARK, VIS_GEOFACE);

    /* create Connect object and generate */
    connectext = vis_ConnectBegin();
    vis_ConnectPre(connectext, SYS_DOUBLE);
    vis_ConnectDef(connectext, 0, 0);
    msh_ExtMeshGenerate(extmesh, connectext);
    if (msh_ExtMeshError(extmesh)) {
        printf("extrusion mesh generation error\n");
        exit(1);
    }
    else {
        printf("extrusion mesh generation complete\n");
    }
    /* write generated ext mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connectext, SYS_NASTRAN_BULKDATA, "exam42b_ext.bdf");

    /* perform tet meshing of the remaining interior */
    /* create TetMesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, numnpsurf, numelsurf, 2);
    msh_TetMeshSetParami(tetmesh, VIS_MESH_INWARD, SYS_ON);

    /* fill TetMesh object with generated points on surface */
    vis_ConnectNumber(connectext, SYS_ELEM, &numelext);
    vis_ConnectNumber(connectext, SYS_NODE, &numnpext);
    numnp = 0;
    idtran = vis_IdTranBegin();
    /* find nodes on surface */
    /* first generate face group */
    fgroup = vis_GroupBegin();
    vis_GroupDef(fgroup, numelext, SYS_ELEM, SYS_FACE);
    vis_ConnectSetGroupParami(connectext, CONNECT_ASSOCTYPE, VIS_GEOFACE);
    vis_ConnectSetGroupParami(connectext, CONNECT_ASSOCID, 2);
    vis_ConnectKernel(connectext, 0);
    vis_ConnectFaceGroup(connectext, CONNECT_ASSOC, NULL, fgroup);
    /* second generate unique node group from face group */
    ngroup = vis_GroupBegin();
    vis_GroupDef(ngroup, numnpext, SYS_NODE, SYS_NONE);
    vis_ConnectNodeGroup(connectext, CONNECT_UNIQUE, fgroup, ngroup);
    /* set nodes */
    for (i = 1; i <= numnpext; i++) {
        if (!vis_GroupNode(ngroup, i))
            continue;
        numnp++;
        vis_IdTranSetId(idtran, numnp, i);
        vis_ConnectCoordsdv(connectext, 1, &i, &u);
        msh_TetMeshSetNode(tetmesh, numnp, u);
    }
    vis_GroupEnd(ngroup);
    vis_GroupEnd(fgroup);

    /* fill TetMesh object with generated faces on surface */
    numel = 0;
    for (i = 1; i <= numelext; i++) {
        vis_ConnectElemNum(connectext, SYS_FACE, i, &num);
        for (no = 1; no <= num; no++) {
            vis_ConnectElemEntAssoc(connectext, VIS_GEOFACE, SYS_FACE, i, no, &aid);
            if (aid != 2)
                continue;
            vis_ConnectElemCon(connectext, SYS_FACE, i, no, &nix, ix);
            vis_ConnectCoordsdv(connectext, nix, ix, x);
            vis_IdTranIndices(idtran, nix, ix, jx);
            numel++;
            msh_TetMeshSetTri(tetmesh, numel, jx);
        }
    }
    /* create Connect object to hold generated tet mesh */
    connecttet = vis_ConnectBegin();
    vis_ConnectPre(connecttet, SYS_DOUBLE);
    vis_ConnectDef(connecttet, 0, 0);

    /* generate mesh */
    msh_TetMeshGenerate(tetmesh, connecttet);
    if (msh_TetMeshError(tetmesh)) {
        printf("tet mesh generation error\n");
        exit(1);
    }
    else {
        printf("tet mesh generation complete\n");
    }
    /* write generated tet mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connecttet, SYS_NASTRAN_BULKDATA, "exam42b_tet.bdf");

    /* merge extruded and tet meshes */
    connectall = vis_ConnectBegin();
    vis_ConnectPre(connectall, SYS_DOUBLE);
    vis_ConnectAppend(connectall, connectext);
    vis_ConnectAppend(connectall, connecttet);
    vis_ConnectSetParamf(connectall, CONNECT_TOLERANCE, 0.0);
    vis_ConnectMerge(connectall, NULL);
    vis_ConnectWrite(connectall, SYS_NASTRAN_BULKDATA, "exam42b_all.bdf");

    /* end objects */
    free(normal);
    vis_IdTranEnd(idtran);
    vis_ConnectEnd(connectall);
    vis_ConnectEnd(connectext);
    vis_ConnectEnd(connectsurf);
    vis_ConnectEnd(connecttet);
    msh_ExtMeshEnd(extmesh);
    msh_SurfMeshEnd(surfmesh);
    msh_TetMeshEnd(tetmesh);
    return 0;
}

6.6. Example 42c, Generate an Extruded Mesh with Interface Layers

This example illustrates generating a normal-extruded mesh where the requested extrusion distances would lead to collisions between elements. These collisions can be avoided by enabling a global recession algorithm using the integer parameter EXTMESH_GLOBALRECEDE.

The example also illustrates the usage of 3 paths in the extrusion direction, each with its own fraction of the normal distance and number of layers. Interface elements are requested at the end of each extrusion path.

The parameter EXTMESH_FACEMARK is enabled to mark all element faces at the boundaries of each extrusion path.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble x[16][3] = {{0., 1., 1.}, {1., 1., 1.}, {2., 0., 1.}, {3., 0., 1.}, {4., 1., 1.}, {5., 1., 1.},
                           {0., 1., 0.}, {1., 1., 0.}, {2., 0., 0.}, {3., 0., 0.}, {4., 1., 0.}, {5., 1., 0.},
                           {0., 2., 0.}, {5., 2., 0.}, {0., 2., 1.}, {5., 2., 1.}};

static Vdouble xn[16][3] = {{0., 1., 0.},  {1., 1., 0.},  {1., 1., 0.},  {-1., 1., 0.}, {-1., 1., 0.}, {0., 1., 0.},
                            {0., 1., 0.},  {1., 1., 0.},  {1., 1., 0.},  {-1., 1., 0.}, {-1., 1., 0.}, {0., 1., 0.},
                            {0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}};

static Vint pix[6][4] = {{1, 2, 8, 7}, {2, 3, 9, 8}, {3, 4, 10, 9}, {4, 5, 11, 10}, {5, 6, 12, 11}, {13, 14, 16, 15}};

/*----------------------------------------------------------------------
                      Generate an Extruded Mesh with Interface Layers
----------------------------------------------------------------------*/
int
main()
{
    msh_ExtMesh* extmesh;
    vis_Connect* connect;
    Vint numnp, numel, i, j, nix, ix[8];
    Vdouble y[3];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create ExtMesh object */
    extmesh = msh_ExtMeshBegin();
    msh_ExtMeshDef(extmesh, 16, 6);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < 16; i++) {
        msh_ExtMeshSetNode(extmesh, i + 1, x[i]);
        msh_ExtMeshSetNodeNormal(extmesh, i + 1, xn[i]);
    }
    /* define linear quadrilateral connectivity */
    for (i = 0; i < 6; i++) {
        msh_ExtMeshSetElem(extmesh, i + 1, VIS_SHAPEQUAD, 2, 0, pix[i]);
    }
    /* set mesh parameters */
    msh_ExtMeshSetParami(extmesh, EXTMESH_NORMAL, SYS_ON);
    msh_ExtMeshSetParami(extmesh, EXTMESH_GLOBALRECEDE, SYS_ON);
    msh_ExtMeshSetParami(extmesh, EXTMESH_FACEMARK, VIS_GEOFACE);
    /* 3 sets of layers with an interface */
    msh_ExtMeshSetNormalPath(extmesh, 4, 0.25, SYS_ON);
    msh_ExtMeshSetNormalPath(extmesh, 5, 0.5, SYS_ON);
    msh_ExtMeshSetNormalPath(extmesh, 4, 0.25, SYS_ON);

    /* generate */
    msh_ExtMeshGenerate(extmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, &y);
        printf("id= %d  x= %e, y= %e, z= %e\n", i, y[0], y[1], y[2]);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix=", i);
        for (j = 0; j < nix; j++) {
            printf(" %d", ix[j]);
        }
        printf("\n");
    }
    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam42c.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_ExtMeshEnd(extmesh);
    return 0;
}

6.7. Example 43, Generate a 2D Planar Tri Mesh

This example illustrates generating a triangle mesh in a unit square with a rectangular hole, an embedded line loop and a short fixed interior edge. The input geometry consists of 12 nodes and 16 boundary lines. The number of input nodes and boundary lines are specified by calling msh_TriMeshDef(). Each node and boundary line are input using msh_TriMeshSetNode() and msh_TriMeshSetLine() respectively. Note that the embedded line loop (line segments 9 thru 14 connecting points 5, 6, 11 and 12) are defined by six back-to-back line segments. The hard edge (line segments 15 and 16) is specified using back-to-back line segments and disabling splitting during mesh generation using msh_TriMeshSetLineStat(). The default mesh sizing is specified by msh_TriMeshSetParamd(), and local mesh size at node 1 is specified using msh_TriMeshSetNodeSizing(). A possible addition to the local mesh sizing specification would be to compute the length of the hard edge and set that length as the local mesh size at the end nodes of the hard edge. This would tend to create equilateral triangles about the hard edge which would gradually increase in size (depending upon the growth rate) to the default mesh size.

Two varieties of line associations are set. The edge associations will be propagated into the final mesh as element edge entity associations on the generated triangle edges lying on the line. All nodes generated along the line will be given node associations. The line element associations will be propagated into the final mesh as element associations on all triangles bounded by the line and any other input lines.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 12
#define MAX_LINE 16

/* a square with a rectangular hole and fixed edge */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},  {.3, .3, 0.},   {.7, .3, 0.},
                                      {.3, .7, 0.}, {.7, .7, 0.}, {.8, .8, 0.}, {.91, .9, 0.}, {.15, .15, 0.}, {.85, .15, 0.}};

/* boundary lines */
static Vint lines[MAX_LINE][2] = {{1, 2},
                                  {2, 4},
                                  {4, 3},
                                  {3, 1},
                                  {5, 7},
                                  {7, 8},
                                  {8, 6},
                                  {6, 5},
                                  {11, 5},
                                  {5, 11},
                                  {12, 6},
                                  {6, 12},
                                  {11, 12},
                                  {12, 11},
                                  /* interior edge */
                                  {9, 10},
                                  {10, 9}};

/*----------------------------------------------------------------------
                      Generate a 2D Planar Tri Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;

    Vint i, k;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[6];
    Vint aid, num, aids[3];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* set line 16 as an undivided edge */
    msh_TriMeshSetLineStat(trimesh, 16, 0);

    /* tag the first line segment */
    msh_TriMeshSetLineAssoc(trimesh, VIS_MISCID, 1, SYS_EDGE, 1001);
    /* elements inside region with boundary edge 10 or 13 */
    msh_TriMeshSetLineAssoc(trimesh, VIS_MISCID1, 10, SYS_ELEM, 2000);
    msh_TriMeshSetLineAssoc(trimesh, VIS_MISCID1, 10, SYS_ELEM, 2002);
    msh_TriMeshSetLineAssoc(trimesh, VIS_MISCID1, 13, SYS_ELEM, 2001);
    /* tag the first node */
    msh_TriMeshSetNodeAssoc(trimesh, VIS_MISCID2, 1, 1002);

    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, 1);

    /* set mesh size */
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .1);

    /* refine about node 1 */
    msh_TriMeshSetNodeSizing(trimesh, 1, .02);

    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43.tri");
    /* generate */
    msh_TriMeshGenerate(trimesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    /* check for nodes generated on tagged line segment */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %13.10f, y= %13.10f, z= %13.10f\n", i, x[0], x[1], x[2]);
        vis_ConnectNodeAssoc(connect, VIS_MISCID2, 1, &i, &aid);
        if (aid) {
            printf("   VIS_MISCID2= %d\n", aid);
        }
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
        for (k = 1; k <= 3; k++) {
            vis_ConnectElemEntAssoc(connect, VIS_MISCID, SYS_EDGE, i, k, &aid);
            if (aid) {
                printf(" edge %d,   VIS_MISCID= %d\n", k, aid);
            }
        }
        vis_ConnectAllElemAssoc(connect, VIS_MISCID1, i, &num, aids);
        for (k = 0; k < num; k++) {
            if (aids[k]) {
                printf(" elem,   VIS_MISCID1= %d\n", aids[k]);
            }
        }
    }
    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    return 0;
}

6.8. Example 43a, Generate a 2D Planar Tri Mesh with Parabolic Lines

This example illustrates generating a triangle mesh in a unit square with a circular hole. The input boundary segments are parabolic. Note that in the output mesh all segment endpoint node numbers are retained. In general the midside node numbers will be changed.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define SQUARE_NODE 8
#define SQUARE_LINE 4

/* outside square */
static Vdouble coords[SQUARE_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {1., 1., 0.}, {0., 1., 0.},
                                         {.5, 0., 0.}, {1., .5, 0.}, {.5, 1., 0.}, {0., .5, 0.}};

/* square boundary lines */
static Vint lines[SQUARE_LINE][3] = {{1, 5, 2}, {2, 6, 3}, {3, 7, 4}, {4, 8, 1}};

/*----------------------------------------------------------------------
                      Generate a 2D Planar Tri Mesh with Parabolic Lines
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;

    Vint i;
    Vint numpts, numlns, lix[3], numcirclelines;
    Vint numnp, numel;
    Vdouble ang, rad, xc[3], x[3];
    Vint nix, ix[6];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    numcirclelines = 8;
    numpts = SQUARE_NODE + 2 * numcirclelines;
    numlns = SQUARE_LINE + numcirclelines;
    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, numpts, numlns, 3);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < SQUARE_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* generate circle nodes */
    rad = .25;
    ang = 0.;
    xc[0] = .5;
    xc[1] = .5;
    xc[2] = 0.;
    for (i = 0; i < 2 * numcirclelines; i++) {
        x[0] = rad * cos(ang) + xc[0];
        x[1] = rad * sin(ang) + xc[1];
        x[2] = 0. + xc[2];
        ang += 3.14159 / numcirclelines;
        msh_TriMeshSetNode(trimesh, SQUARE_NODE + i + 1, x);
    }
    /* define lines */
    for (i = 0; i < SQUARE_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* generate circle lines */
    for (i = 0; i < numcirclelines; i++) {
        if (i == numcirclelines - 1) {
            lix[0] = SQUARE_NODE + 1;
        }
        else {
            lix[0] = SQUARE_NODE + 2 * i + 3;
        }
        lix[1] = SQUARE_NODE + 2 * i + 2;
        lix[2] = SQUARE_NODE + 2 * i + 1;
        msh_TriMeshSetLine(trimesh, SQUARE_LINE + i + 1, lix);
    }
    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, 1);

    /* set mesh size */
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .05);

    /* write debugging file */
    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43a.tri");
    /* generate */
    msh_TriMeshGenerate(trimesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %f, y= %f, z= %f\n", i, x[0], x[1], x[2]);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
    }
    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43a.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    return 0;
}

6.9. Example 43b, Generate a 2D Planar Tri Mesh with Sizing Function

This example is an extension of Example 43 to illustrate using user-defined isotropic and anisotropic sizing functions to affect the sizing of the generated triangular mesh. Any sizing function set by the user will override any other sizing specifications set using functions such as msh_TriMeshSetParamd() or msh_TriMeshSetNodeSizing(). In addition, any anisotropic sizing function set by the user takes precedence over an isotropic sizing function which may also have been set. In this example the anisotropic sizing function, s_aniso, is set using the function msh_TriMeshSetFunction(). The anisotropic sizing function returns three orthogonal vectors, each vector is scaled to the desired size in that direction. The generated triangle aspect ratios will be 4 and inclined at an angle of 45 degrees.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 8
#define MAX_LINE 8

/* a square with a rectangular hole */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                                      {.2, .2, 0.}, {.8, .2, 0.}, {.2, .8, 0.}, {.8, .8, 0.}};

static Vint lines[MAX_LINE][2] = {{1, 2}, {2, 4}, {4, 3}, {3, 1}, {5, 7}, {7, 8}, {8, 6}, {6, 5}};

/* isotropic size function */
static void
s_iso(msh_TriMesh* trimesh, Vobject* usrobj, Vdouble x[3], Vdouble* s)
{
    printf("s_iso called\n");
    /* linear variation in y direction */
    *s = .05 + .05 * x[1];
    printf("x= %e %e %e, s= %e\n", x[0], x[1], x[2], *s);
}
/* anisotropic size function */
static void
s_aniso(msh_TriMesh* trimesh, Vobject* usrobj, Vdouble x[3], Vdouble s[3][3])
{
    Vdouble s1, s2, s3;
    printf("s_aniso called\n");
    /* 4 to 1 aspect ratio at 45 degree angle */
    s1 = .025;
    s[0][0] = s1 * .707;
    s[0][1] = s1 * .707;
    s[0][2] = 0.;

    s2 = .1;
    s[1][0] = s2 * -.707;
    s[1][1] = s2 * .707;
    s[1][2] = 0.;
    /* set the out of plane direction size to unity */
    s3 = 1.;
    s[2][0] = 0.;
    s[2][1] = 0.;
    s[2][2] = s3 * 1.;
    printf("x= %e %e %e, s1= %e, s2= %e\n", x[0], x[1], x[2], s1, s2);
}

/*----------------------------------------------------------------------
                      Generate a 2D Planar Tri Mesh with Sizing Function
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;

    Vint i, k;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[6];
    Vint aid;
    Vint anisoflag;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* tag the first line segment */
    msh_TriMeshSetLineAssoc(trimesh, VIS_MISCID, 1, SYS_EDGE, 1001);

    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, VIS_ON);

    /* use zero for isotropic sizing, one for anisotropic */
    anisoflag = 1;
    /* set sizing function */
    if (anisoflag == 0) {
        msh_TriMeshSetFunction(trimesh, TRIMESH_FUN_SIZING, (Vfunc*)s_iso, NULL);
    }
    else {
        msh_TriMeshSetFunction(trimesh, TRIMESH_FUN_ASIZING, (Vfunc*)s_aniso, NULL);
    }

    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43b.tri");
    /* generate */
    msh_TriMeshGenerate(trimesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
    }
    /* print element information */
    /* check for edges generated on tagged line segment */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
        for (k = 1; k <= 3; k++) {
            vis_ConnectElemEntAssoc(connect, VIS_MISCID, SYS_EDGE, i, k, &aid);
            if (aid) {
                printf(" edge %d,   VIS_MISCID= %d\n", k, aid);
            }
        }
    }
    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43b.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    return 0;
}
../_images/vistools-exam43b.gif

Figure 28-43b, Generated Anisotropic Triangle Mesh

6.10. Example 43c, Generate a 2D Planar Quad Mesh

This example illustrates using TriMesh to generate a planar triangulation which can then be used as the surface geometry for a subsequent mesh generated using SurfMesh.

The planar geometry consists of a unit square with a fixed interior edge. The desired goal of TriMesh is to generate a triangulation of the input geometry without generating any additional nodes on the boundary or interior. Do this by calling msh_TriMeshSetParami() and disable the parameters VIS_MESH_BOUNDREFINE and VIS_MESH_INTERREFINE. It is useful to know which edges of the generated triangles lie on input boundary lines. Use msh_TriMeshSetParami() and set line flag association to VIS_MISCID1. This will enable generation of the element association consisting of bitwise flags of triangle edges which lie on input boundary lines. The triangle mesh is generated using msh_TriMeshGenerate(). The generated triangle mesh is illustrated in the first figure below. Note that the number of nodes in the output Connect object, numnp is equal to the number of nodes (MAX_NODE) input to TriMesh. These nodes are set in SurfMesh using msh_SurfMeshSetPoint(). The triangles are set using msh_SurfMeshSetTri(). The edge flags are queried using vis_ConnectElemAssoc() with association type VIS_MISCID1. The edge flags are parsed bit by bit and the corresponding triangle edges are marked as preserved. Use msh_SurfMeshSetParami() and msh_SurfMeshSetParamd() to set meshing parameters. In this case linear quadrilateral elements are generated with a mesh size of .1. Note that the preserved edges appear in the final mesh. The generated quadrilateral mesh is illustrated in the second figure below.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 6
#define MAX_LINE 6

/* a square with a fixed interior edge */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.}, {.2, .3, 0.}, {.7, .8, 0.}};

/* boundary lines */
static Vint lines[MAX_LINE][2] = {{1, 2},
                                  {2, 4},
                                  {4, 3},
                                  {3, 1},
                                  /* interior edge */
                                  {5, 6},
                                  {6, 5}};

/*----------------------------------------------------------------------
                      Generate a 2D Planar Quad Surface Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect *connect, *connectsurf;
    msh_SurfMesh* surfmesh;

    Vint i, j;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[6], efl[3];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* generate boundary triangulation only */
    msh_TriMeshSetParami(trimesh, VIS_MESH_BOUNDREFINE, VIS_OFF);
    msh_TriMeshSetParami(trimesh, VIS_MESH_INTERREFINE, VIS_OFF);
    /* element associations with line flags */
    msh_TriMeshSetParami(trimesh, TRIMESH_LINEFLAG, VIS_MISCID1);

    /* generate */
    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43c.tri");
    msh_TriMeshGenerate(trimesh, connect);
    /* delete */
    msh_TriMeshEnd(trimesh);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("TriMesh, numnp= %d, numel= %d\n", numnp, numel);

    surfmesh = msh_SurfMeshBegin();
    msh_SurfMeshDef(surfmesh, numnp, numel);

    /* create connect object to receive surf mesh */
    connectsurf = vis_ConnectBegin();
    vis_ConnectPre(connectsurf, SYS_DOUBLE);
    vis_ConnectDef(connectsurf, 0, 0);

    /* load nodes and preserve them */
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        msh_SurfMeshSetPoint(surfmesh, i, x, 1);
    }
    /* load elements */
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        for (j = 0; j < 3; j++) {
            vis_ConnectElemEntAssoc(connect, VIS_MISCID1, SYS_EDGE, i, j + 1, &efl[j]);
        }
        msh_SurfMeshSetTri(surfmesh, i, ix, efl);
    }
    /* delete */
    vis_ConnectEnd(connect);

    /* set mesh parameters */
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_MAXI, 2);
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_SHAPE, VIS_SHAPEQUAD);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, .1);

    /* generate */
    msh_SurfMeshWrite(surfmesh, SYS_ASCII, "exam43c.srf");
    msh_SurfMeshGenerate(surfmesh, connectsurf);
    /* print generated nodes and elements */
    vis_ConnectNumber(connectsurf, SYS_NODE, &numnp);
    vis_ConnectNumber(connectsurf, SYS_ELEM, &numel);
    printf("SurfMesh, numnp= %d, numel= %d\n", numnp, numel);

    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connectsurf, SYS_NASTRAN_BULKDATA, "exam43c.bdf");

    /* end objects */
    vis_ConnectEnd(connectsurf);
    msh_SurfMeshEnd(surfmesh);
    return 0;
}
../_images/vistools-exam43c1.gif

Figure 28-43c1, Generated Triangle Mesh

../_images/vistools-exam43c2.gif

Figure 28-43c2, Generated Quadrilateral Mesh

6.11. Example 43d, Refine a 2D Planar Tri Mesh

This example illustrates using TriMesh to generate and subsequently refine a planar triangulation. The original mesh is generated with an edge length of .2. The refined mesh sets the edge length to .1 for the elements with any generated node in the upper right corner of the original mesh.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 8
#define MAX_LINE 8

/* a square with a rectangular hole */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                                      {.2, .2, 0.}, {.8, .2, 0.}, {.2, .8, 0.}, {.8, .8, 0.}};

static Vint lines[MAX_LINE][2] = {{1, 2}, {2, 4}, {4, 3}, {3, 1}, {5, 7}, {7, 8}, {8, 6}, {6, 5}};

/*----------------------------------------------------------------------
                      Refine a 2D Planar Tri Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;
    vis_State* state;

    Vint i;
    Vint numnp, numel;
    Vdouble x[3], v;
    Vint nix, ix[6];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, VIS_ON);
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .2);

    /* generate initial mesh */
    msh_TriMeshGenerate(trimesh, connect);
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("Generated original mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43d.bdf");

    /* create State object with new element edge lengths */
    state = vis_StateBegin();
    vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
    /* refine elements in upper right corner */
    for (i = 1; i <= numel; ++i) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        vis_ConnectCoordsdv(connect, 1, &ix[0], (Vdouble(*)[3])x);
        if (x[0] >= .8 && x[1] >= .8) {
            v = .1;
            vis_StateSetDatadv(state, i, &v);
        }
    }
    /* clear Connect object */
    vis_ConnectDef(connect, 0, 0);
    /* generate refined mesh */
    msh_TriMeshRefine(trimesh, state, connect);
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43d.ref.bdf");

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("Generated refined mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);

    /* print node information */
    /* check for nodes generated on tagged line segment */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
    }

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    vis_StateEnd(state);
    return 0;
}

6.12. Example 43e, Refine a 2D Planar Tri Mesh as an Orphan Mesh

This example is a modification of Example 43d using TriMesh to generate a mesh but then for purposes of illustration, the original mesh is erased using TriMeshDef and the original mesh is reentered as an orphan mesh for subsequent refinement. The original mesh is generated with an edge length of .2. The refined mesh sets the edge length to .1 for the elements with any generated node in the upper right corner of the original mesh.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 8
#define MAX_LINE 8

/* a square with a rectangular hole */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                                      {.2, .2, 0.}, {.8, .2, 0.}, {.2, .8, 0.}, {.8, .8, 0.}};

static Vint lines[MAX_LINE][2] = {{1, 2}, {2, 4}, {4, 3}, {3, 1}, {5, 7}, {7, 8}, {8, 6}, {6, 5}};

/*----------------------------------------------------------------------
                      Refine a 2D Planar Tri Mesh as on Orphan Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;
    vis_State* state;

    Vint i;
    Vint numnp, numel;
    Vdouble x[3], v;
    Vint nix, ix[6];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, VIS_ON);
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .2);

    /* generate initial mesh */
    msh_TriMeshGenerate(trimesh, connect);
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("Generated original mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43e.bdf");

    /* redefine trimesh to setup as taking orphan mesh */
    msh_TriMeshDef(trimesh, numnp, MAX_LINE, 2);
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, &x);
        msh_TriMeshSetNode(trimesh, i, x);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* set trimesh as background mesh */
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        msh_TriMeshSetElem(trimesh, i, SYS_SHAPETRI, 3, ix);
    }
    /* create State object with new element edge lengths */
    state = vis_StateBegin();
    vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
    /* refine elements in upper right corner */
    for (i = 1; i <= numel; ++i) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        vis_ConnectCoordsdv(connect, 1, &ix[0], (Vdouble(*)[3])x);
        if (x[0] >= .8 && x[1] >= .8) {
            v = .1;
            vis_StateSetDatadv(state, i, &v);
        }
    }
    /* clear Connect object */
    vis_ConnectDef(connect, 0, 0);
    /* generate refined mesh */
    msh_TriMeshRefine(trimesh, state, connect);
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43e_ref.bdf");

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("Generated refined mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);

    /* print node information */
    /* check for nodes generated on tagged line segment */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3], ix[4], ix[5]);
    }

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    vis_StateEnd(state);
    return 0;
}

6.13. Example 43f, Generate a 2D Planar Tri Mesh with Sizing

This example illustrates using node, line and region sizing to control mesh size. The domain consists of two side-by-side unit squares. An overall target edge length is set to .2. An unconnected point is placed in the center of the left region (region 1) and a size of .01 is specified there using msh_TriMeshSetNodeSizing(). A line sizing of .05 is placed on the left edge of region 1 using msh_TriMeshSetLineSizing() using an entity type of SYS_EDGE. A region sizing of .05 is defined for region 2 using the left edge of region 2. The function msh_TriMeshSetLineSizing() with an entity type of SYS_ELEM is used to set the region sizing. Note that the size for a region is set to the minimum of all region (SYS_ELEM) sizes set on lines which bound the region. The default growth rate of 1.5 is used to transition between regions of various sizes.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

#define MAX_NODE 7
#define MAX_LINE 8

/* 2 x 1 rectangle divided in the middle into two regions */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                                      {2., 0., 0.}, {2., 1., 0.}, {.5, .5, 0.}};

/* boundary lines of two regions */
static Vint lines[MAX_LINE][2] = {{1, 2}, {2, 4}, {4, 3}, {3, 1}, {2, 5}, {5, 6}, {6, 4}, {4, 2}};

/*----------------------------------------------------------------------
                      Generate a 2D Planar Tri Mesh with Sizing
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;
    Vint i;
    Vint numnp, numel;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* tag the two regions */
    msh_TriMeshSetLineAssoc(trimesh, VIS_GEOBODY, 1, SYS_ELEM, 1);
    msh_TriMeshSetLineAssoc(trimesh, VIS_GEOBODY, 5, SYS_ELEM, 2);

    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, 1);
    /* allow unconnected node */
    msh_TriMeshSetParami(trimesh, TRIMESH_UNCONNECT, SYS_ON);

    /* set mesh size */
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .2);

    /* refine about node 7 in region 1 */
    msh_TriMeshSetNodeSizing(trimesh, 7, .01);

    /* refine edge 4 in region 1 */
    msh_TriMeshSetLineSizing(trimesh, 4, SYS_EDGE, .05);
    /* refine region 2 */
    msh_TriMeshSetLineSizing(trimesh, 5, SYS_ELEM, .05);

    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43f.tri");
    /* generate */
    msh_TriMeshGenerate(trimesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43f.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    return 0;
}

6.14. Example 43g, Generate a 2D Planar Tri Mesh with Interior Sizing

This example is an alteration of Example 43f and illustrates specifying interior sizing. A boundary depth size is specified using msh_TriMeshSetLineDepthSizing(). This function allows a depth and size to be set for each boundary line. This feature is useful for generating isotropic boundary layers. In addition, a user specified element size is defined with a simple geometry shape in the interior of the domain using msh_TriMeshSetGeomSizing(). In this case a .05 by .15 by .25 ellipsoid. Note that only the x and y extents of the geometry are used for the sizing region.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"
#define MAX_NODE 6
#define MAX_LINE 8

/* 2 x 1 rectangle divided in the middle into two regions */
static Vdouble coords[MAX_NODE][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.}, {2., 0., 0.}, {2., 1., 0.}};

/* boundary lines of two regions */
static Vint lines[MAX_LINE][2] = {{1, 2}, {2, 4}, {4, 3}, {3, 1}, {2, 5}, {5, 6}, {6, 4}, {4, 2}};

/*----------------------------------------------------------------------
                      Generate a 2D Planar Tri Mesh with Interior Sizing
----------------------------------------------------------------------*/
int
main()
{
    msh_TriMesh* trimesh;
    vis_Connect* connect;
    Vint i;
    Vint numnp, numel;
    Vdouble xo[3], xa[3], xb[3], d[3];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create trimesh object */
    trimesh = msh_TriMeshBegin();
    msh_TriMeshDef(trimesh, MAX_NODE, MAX_LINE, 2);

    /* create connect object to receive mesh */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < MAX_NODE; i++) {
        msh_TriMeshSetNode(trimesh, i + 1, coords[i]);
    }
    /* define lines */
    for (i = 0; i < MAX_LINE; i++) {
        msh_TriMeshSetLine(trimesh, i + 1, lines[i]);
    }
    /* define line depth sizing in region 2 */
    for (i = 4; i < MAX_LINE; i++) {
        msh_TriMeshSetLineDepthSizing(trimesh, i + 1, .1, .02);
    }
    /* tag the two regions */
    msh_TriMeshSetLineAssoc(trimesh, VIS_GEOBODY, 1, SYS_ELEM, 1);
    msh_TriMeshSetLineAssoc(trimesh, VIS_GEOBODY, 5, SYS_ELEM, 2);

    /* generate quadratic elements */
    msh_TriMeshSetParami(trimesh, VIS_MESH_MAXI, 3);
    msh_TriMeshSetParami(trimesh, TRIMESH_SMOOTH, 1);

    /* set mesh size */
    msh_TriMeshSetParamd(trimesh, VIS_MESH_EDGELENGTH, .2);

    /* interior ellipsoidal shape for sizing */
    /* positioned in center of cube */
    xo[0] = .5;
    xo[1] = .5;
    xo[2] = 0.;
    /* oriented along x axis and at 45 degrees in y-z plane */
    xa[0] = 1.0;
    xa[1] = .5;
    xa[2] = .5;
    xb[0] = .5;
    xb[1] = 1.0;
    xb[2] = .5;
    /* lengths of ellipsoid along each axis */
    d[0] = .05;
    d[1] = .15;
    d[2] = .25;
    msh_TriMeshSetGeomSizing(trimesh, VIS_MESH_GEOMSPHERE, xo, xa, xb, d, 0.02);

    msh_TriMeshWrite(trimesh, SYS_ASCII, "exam43g.tri");
    /* generate */
    msh_TriMeshGenerate(trimesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* write out NASTRAN bulk data file */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam43g.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TriMeshEnd(trimesh);
    return 0;
}
../_images/vistools-exam43g.gif

Figure 28-48f, Generated Planar Tri Mesh with Interor Sizing

6.15. Example 44, Generate a 3D Volume Tet Mesh

This example illustrates generating a simple tetrahedral mesh in the interior of a cube using the TetMesh module. The input surface geometry consists of the 8 corner nodes of the cube and 12 triangular surface facets, 2 per face of the cube. The connectivity of the surface triangles is such that the normal to the triangles, determined by the right hand rule, points into the volume to be meshed. The number of input nodes and triangles are specified by calling msh_TetMeshDef(), then each node and triangle is input using msh_TetMeshSetNode() and msh_TetMeshSetTri() respectively.

The nodes and finite elements are actually generated by calling msh_TetMeshGenerate(). In this simple case no internal nodes are generated. The generated nodes and elements are entered into the Connect object which is passed in the argument list. Lastly, the generated nodes and elements are printed.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                               {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};

/*----------------------------------------------------------------------
                      Generate a 3D Volume Tet Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect* connect;

    Vint i, j;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[4];
    Vint num, aid, aids[10];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, 8, 12, 2);

    /* create connect object */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < 8; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
    }
    /* define a node association on node 2 */
    msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, 2, 1000);
    /* define an element edge association on tri 4, edge 2 */
    msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, 4, SYS_EDGE, 2, 10);
    /* define an element face association on tri 3 */
    msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOFACE, 3, SYS_FACE, 1, 100);
    /* define an element association on tri 3 */
    msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOBODY, 3, SYS_ELEM, 0, 200);
    msh_TetMeshSetTriAssoc(tetmesh, VIS_PROPID, 3, SYS_ELEM, 0, 1);
    /* set edgelength */
    msh_TetMeshSetParamd(tetmesh, VIS_MESH_EDGELENGTH, 1.);

    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam44.tet");
    /* generate */
    msh_TetMeshGenerate(tetmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %f, y= %f, z= %f\n", i, x[0], x[1], x[2]);
    }
    /* node associations */
    for (i = 1; i <= numnp; i++) {
        vis_ConnectAllNodeAssoc(connect, VIS_GEOVERT, i, &num, aids);
        for (j = 0; j < num; j++) {
            if (aids[j] == 0)
                continue;
            printf("id= %d  aid[%d]= %d\n", i, j, aids[j]);
        }
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3]);
    }
    /* element entity associations */
    for (i = 1; i <= numel; i++) {
        /* tet edges */
        for (j = 1; j <= 6; j++) {
            vis_ConnectElemEntAssoc(connect, VIS_GEOEDGE, SYS_EDGE, i, j, &aid);
            if (aid == 0)
                continue;
            printf("id= %d edge= %d, aid= %d\n", i, j, aid);
        }
        /* tet faces */
        for (j = 1; j <= 4; j++) {
            vis_ConnectElemEntAssoc(connect, VIS_GEOFACE, SYS_FACE, i, j, &aid);
            if (aid == 0)
                continue;
            printf("id= %d face= %d, aid= %d\n", i, j, aid);
        }
    }
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam44.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TetMeshEnd(tetmesh);
    return 0;
}

6.16. Example 44a, Generate a 3D Volume Tet Mesh with Interior Features

This example is an extension of Example 44 with the addition of an interior void, baffle face, floating face, floating edge and node. By default, a node point which been introduced in the interior which is not connected to any surface triangle, interior face or edge is reported as an error. The parameter TETMESH_UNCONNECT must be enabled to allow unconnected nodes to be inserted in the final mesh.

A monitor function has been installed to follow the progress of the meshing process. The function msh_TetMeshGetInteger() is called with type VIS_MESH_PROGRESS within the monitor function to query and print the current state of the meshing process.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/* fixed node points */
static Vdouble coords[25][3] = {{0., 0., 0.},
                                {5., 0., 0.},
                                {0., 3., 0.},
                                {5., 3., 0.},
                                {0., 0., 3.},
                                {5., 0., 3.},
                                {0., 3., 3.},
                                {5., 3., 3.},
                                /* unit cube void */
                                {1., 1., 1.},
                                {2., 1., 1.},
                                {1., 2., 1.},
                                {2., 2., 1.},
                                {1., 1., 2.},
                                {2., 1., 2.},
                                {1., 2., 2.},
                                {2., 2., 2.},
                                /* interior baffle face */
                                {3., 1., 1.},
                                {3., 1., 2.},
                                /* interior edge */
                                {4., 2., 1.},
                                {4., 2., 2.},
                                /* interior point */
                                {3., 2., 1.},
                                /* interior floating face */
                                {3.5, 1., 1.},
                                {3.5, 1., 2.},
                                {4.5, 1., 1.},
                                {4.5, 1., 2.}};

/* surface triangles */
static Vint tris[32][3] = {{1, 5, 3},
                           {3, 5, 7},
                           {4, 6, 2},
                           {8, 6, 4},
                           {1, 2, 5},
                           {5, 2, 6},
                           {7, 4, 3},
                           {8, 4, 7},
                           {1, 4, 2},
                           {4, 1, 3},
                           {5, 8, 7},
                           {5, 6, 8},
                           /* unit cube void */
                           {9, 11, 13},
                           {11, 15, 13},
                           {12, 10, 14},
                           {16, 12, 14},
                           {9, 13, 10},
                           {13, 14, 10},
                           {15, 11, 12},
                           {16, 15, 12},
                           {9, 10, 12},
                           {12, 11, 9},
                           {13, 15, 16},
                           {13, 16, 14},
                           /* interior double backed triangles */
                           {14, 18, 10},
                           {14, 10, 18},
                           {10, 18, 17},
                           {10, 17, 18},
                           /* interior double backed floating triangles */
                           {22, 25, 24},
                           {22, 24, 25},
                           {23, 22, 25},
                           {23, 25, 22}};

/* interior edge */
static Vint inedge[2] = {19, 20};

/*----------------------------------------------------------------------
                      monitor function
----------------------------------------------------------------------*/
static void
monitor(msh_TetMesh* tetmesh, Vobject* usrobj)
{
    int iparam[4];

    printf("monitor called\n");
    /* query and print progress */
    msh_TetMeshGetInteger(tetmesh, VIS_MESH_PROGRESS, iparam);
    printf("phase=   %d\n", iparam[0]);
    printf("numnp=   %d\n", iparam[1]);
    printf("numel=   %d\n", iparam[2]);
    printf("percent= %d\n", iparam[3]);
}

/*----------------------------------------------------------------------
                      Generate a 3D Volume Tet Mesh with Interior Features
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect* connect;

    Vint i;
    Vint numnode, numtris;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[4];

    /* number of nodes */
    numnode = 25;
    /* number of triangles */
    numtris = 32;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, numnode, numtris, 2);

    /* create connect object */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < numnode; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
    }
    /* define tris */
    for (i = 0; i < numtris; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
    }
    /* define interior edge */
    msh_TetMeshSetEdge(tetmesh, 1, inedge);

    /* allow unconnected nodes */
    msh_TetMeshSetParami(tetmesh, TETMESH_UNCONNECT, VIS_ON);

    /* set monitor function */
    msh_TetMeshSetFunction(tetmesh, TETMESH_FUN_MONITOR, (Vfunc*)monitor, NULL);

    /* generate */
    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam44a.tet");
    msh_TetMeshGenerate(tetmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %f, y= %f, z= %f\n", i, x[0], x[1], x[2]);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3]);
    }
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam44a.bdf");
    /* end objects */
    vis_ConnectEnd(connect);
    msh_TetMeshEnd(tetmesh);
    return 0;
}

6.17. Example 44b, Generate a 3D Volume Tet Mesh of a Simple Cube

This example is similar to Example 44. However, it eliminates all geometric complexities, meshes a simple cube, and labels vertices, edges, and faces that can be later retrieved on the generated mesh.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                               {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};

/*----------------------------------------------------------------------
                      Generate a 3D Volume Tet Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect* connect;
    vis_State* state;

    Vint i, j;
    Vint numnp, numel;
    Vdouble x[3], v;
    Vint nix, ix[4];
    Vint num, aid, aids[10];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, 8, 12, 2);

    /* create connect object */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define nodes */
    for (i = 0; i < 8; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
    }
    /* define a node association on node 2 */
    msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, 2, 1000);
    /* define an element edge association on tri 4, edge 2 */
    msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, 4, SYS_EDGE, 2, 10);
    /* define an element face association on tri 3 */
    msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOFACE, 3, SYS_FACE, 1, 100);

    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam44b.tet");
    /* generate */
    msh_TetMeshGenerate(tetmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %f, y= %f, z= %f\n", i, x[0], x[1], x[2]);
    }
    /* node associations */
    for (i = 1; i <= numnp; i++) {
        vis_ConnectAllNodeAssoc(connect, VIS_GEOVERT, i, &num, aids);
        for (j = 0; j < num; j++) {
            if (aids[j] == 0)
                continue;
            printf("id= %d  aid[%d]= %d\n", i, j, aids[j]);
        }
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d %d\n", i, ix[0], ix[1], ix[2], ix[3]);
    }
    /* element entity associations */
    for (i = 1; i <= numel; i++) {
        /* tet edges */
        for (j = 1; j <= 6; j++) {
            vis_ConnectElemEntAssoc(connect, VIS_GEOEDGE, SYS_EDGE, i, j, &aid);
            if (aid == 0)
                continue;
            printf("id= %d edge= %d, aid= %d\n", i, j, aid);
        }
        /* tet faces */
        for (j = 1; j <= 4; j++) {
            vis_ConnectElemEntAssoc(connect, VIS_GEOFACE, SYS_FACE, i, j, &aid);
            if (aid == 0)
                continue;
            printf("id= %d face= %d, aid= %d\n", i, j, aid);
        }
    }
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam44b.bdf");

    state = vis_StateBegin();
    vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
    for (i = 1; i <= numel; ++i) {
        v = .5;
        vis_StateSetDatadv(state, i, &v);
    }
    vis_StateWrite(state, NULL, SYS_ASCII, "exam44b.tet.sta1");
    vis_ConnectDef(connect, 0, 0);
    msh_TetMeshRefine(tetmesh, state, connect);
    vis_StateEnd(state);
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "exam44b.ref.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_TetMeshEnd(tetmesh);
    return 0;
}

6.18. Example 44c, Refine a 3D Volume Tet Mesh as an Orphan Mesh

This example extends Example 44b. Once an original tetrahedron mesh is generated its origin is erased by re-Defing the TetMesh object. The original mesh is specified using msh_TetMeshSetElem(). All edges in the cube are assigned a VIS_GEOEDGE association to ensure that they are preserved during mesh refinement.

A uniform new mesh size of 0.25 is employed throughout the volume, as this constant value is set in the State object required for the refinement process.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                               {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};

static Vint eassoc[12][3] = {{1, 0, 2},  {0, 3, 4},  {0, 5, 6}, {7, 0, 8}, {9, 0, 0}, {0, 0, 10},
                             {0, 11, 0}, {0, 0, 12}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}};

static Vint maxi = 2;

/*----------------------------------------------------------------------
      Generate a 3D Volume Tet Mesh and Refine it as an Orphan Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect *connect1, *connect2;
    vis_State* state;

    Vint i, j;
    Vint numnp, numel;
    Vdouble x[3], size;
    Vint nix, ix[4];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, 8, 12, maxi);

    /* create connect object */
    connect1 = vis_ConnectBegin();
    vis_ConnectPre(connect1, SYS_DOUBLE);
    vis_ConnectDef(connect1, 0, 0);

    /* define nodes */
    for (i = 0; i < 8; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
        msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, i + 1, i + 1);
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
        for (j = 0; j < 3; ++j) {
            if (eassoc[i][j]) {
                msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, i + 1, SYS_EDGE, j + 1, eassoc[i][j]);
            }
        }
    }
    /* generate */
    msh_TetMeshGenerate(tetmesh, connect1);
    vis_ConnectNumber(connect1, SYS_NODE, &numnp);
    vis_ConnectNumber(connect1, SYS_ELEM, &numel);
    printf("Generated original mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);
    vis_ConnectWrite(connect1, SYS_NASTRAN_BULKDATA, "exam44c_orig.bdf");

    /* redefine tetmesh to setup as taking orphan mesh */
    msh_TetMeshDef(tetmesh, numnp, 12, maxi);

    /* define nodes */
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect1, 1, &i, &x);
        msh_TetMeshSetNode(tetmesh, i, x);
        msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, i, i);
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
        for (j = 0; j < 3; ++j) {
            if (eassoc[i][j]) {
                msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, i + 1, SYS_EDGE, j + 1, eassoc[i][j]);
            }
        }
    }
    /* set tetmesh as background mesh */
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect1, i, &nix, ix);
        msh_TetMeshSetElem(tetmesh, i, SYS_SHAPETET, maxi, ix);
    }
    /* create state object with new sizes */
    state = vis_StateBegin();
    vis_StatePre(state, SYS_DOUBLE);
    vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
    size = 0.25;
    for (i = 1; i <= numel; i++) {
        vis_StateSetDatadv(state, i, &size);
    }
    /* refine mesh */
    connect2 = vis_ConnectBegin();
    vis_ConnectPre(connect2, SYS_DOUBLE);
    vis_ConnectDef(connect2, 0, 0);
    msh_TetMeshRefine(tetmesh, state, connect2);
    vis_ConnectNumber(connect2, SYS_NODE, &numnp);
    vis_ConnectNumber(connect2, SYS_ELEM, &numel);
    printf("Generated refined mesh\n");
    printf("numnp= %d\n", numnp);
    printf("numel= %d\n", numel);
    vis_ConnectWrite(connect2, SYS_NASTRAN_BULKDATA, "exam44c_ref.bdf");

    /* end objects */
    vis_ConnectEnd(connect1);
    vis_ConnectEnd(connect2);
    msh_TetMeshEnd(tetmesh);
    vis_StateEnd(state);
    return 0;
}

6.19. Example 44d, Successively Refine a 3D Linear Tet Mesh

Unlike Example 44c, this examples assumes that the initial tetrahedron mesh was generated by the same TetMesh object to be used for mesh refinement. As a result, there is no need to re-enter the boundary data or to specify the original mesh with msh_TetMeshSetElem().

The example also illustrates the ability of the msh_TetMeshRefine() method to be called several times. For simplicity, the uniform mesh size is halved from one refinement step to another.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                               {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};

static Vint maxi = 2;

static Vint eassoc[12][3] = {{1, 0, 2},  {0, 3, 4},  {0, 5, 6}, {7, 0, 8}, {9, 0, 0}, {0, 0, 10},
                             {0, 11, 0}, {0, 0, 12}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}};

/*----------------------------------------------------------------------
    Generate a Linear 3D Volume Tet Mesh and Perform Successive Refinements
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect *connect1, *connect2, *connecttmp;
    vis_State* state;
    vis_GridFun* gf;

    Vint i, j, iref;
    Vint numel;
    Vdouble size;
    Vchar buffer[80];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, 8, 12, maxi);

    /* create connect object */
    connect1 = vis_ConnectBegin();
    vis_ConnectPre(connect1, SYS_DOUBLE);
    vis_ConnectDef(connect1, 0, 0);

    /* define nodes */
    for (i = 0; i < 8; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
        /* define a node association on every vertex */
        msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, i + 1, i + 1);
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
        for (j = 0; j < 3; ++j) {
            if (eassoc[i][j]) {
                msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, i + 1, SYS_EDGE, j + 1, eassoc[i][j]);
            }
        }
    }
    /* define an element edge association on tri 4, edge 2 */

    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam44.tet");
    /* generate */
    msh_TetMeshGenerate(tetmesh, connect1);
    vis_ConnectNumber(connect1, SYS_ELEM, &numel);
    printf("Generated original mesh\n");
    vis_ConnectWrite(connect1, SYS_NASTRAN_BULKDATA, "exam44d_orig.bdf");

    /* create state object with new sizes and refine */
    gf = vis_GridFunBegin();
    state = vis_StateBegin();
    vis_StateSetObject(state, VIS_GRIDFUN, (Vobject*)gf);

    connect2 = vis_ConnectBegin();
    vis_ConnectPre(connect2, SYS_DOUBLE);
    size = 1.;
    for (iref = 0; iref < 3; ++iref) {
        size *= 0.5;
        vis_ConnectGridFun(connect1, gf);
        vis_StatePre(state, SYS_DOUBLE);
        vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
        for (i = 1; i <= numel; i++) {
            vis_StateSetDatadv(state, i, &size);
        }
        /* refine mesh */
        vis_ConnectDef(connect2, 0, 0);
        msh_TetMeshRefine(tetmesh, state, connect2);
        printf("Generated refined mesh %d\n", iref + 1);
        sprintf(buffer, "exam44d_ref%d.bdf", iref + 1);
        vis_ConnectWrite(connect2, SYS_NASTRAN_BULKDATA, buffer);
        vis_ConnectNumber(connect2, SYS_ELEM, &numel);

        /* swap connect objects */
        connecttmp = connect1;
        connect1 = connect2;
        connect2 = connecttmp;
    }

    /* end objects */
    vis_ConnectEnd(connect1);
    vis_ConnectEnd(connect2);
    msh_TetMeshEnd(tetmesh);
    vis_StateEnd(state);
    vis_GridFunEnd(gf);
    return 0;
}

6.20. Example 44e, Successively Refine a 3D Parabolic Tet Mesh

This example mimics Example 44d. However, it assumes that the mesh is parabolic.

#include <stdio.h>
#include <math.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble coords[26][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.}, {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.},
                                {1., 1., 1.}, {.5, 0., 0.}, {1., .5, 0.}, {.5, 1., 0.}, {0., .5, 0.}, {.5, 0., 1.}, {1., .5, 1.},
                                {.5, 1., 1.}, {0., .5, 1.}, {0., 0., .5}, {1., 0., .5}, {0., 1., .5}, {1., 1., .5}, {.5, .5, 0.},
                                {.5, .5, 1.}, {0., .5, .5}, {1., .5, .5}, {.5, 0., .5}, {.5, 1., .5}};

static Vint tris[12][6] = {{1, 5, 3, 17, 23, 12}, {3, 5, 7, 23, 16, 19}, {4, 6, 2, 24, 18, 10}, {8, 6, 4, 14, 24, 20},
                           {1, 2, 5, 9, 25, 17},  {5, 2, 6, 25, 18, 13}, {7, 4, 3, 26, 11, 19}, {8, 4, 7, 20, 26, 15},
                           {1, 4, 2, 21, 10, 9},  {4, 1, 3, 21, 12, 11}, {5, 8, 7, 22, 15, 16}, {5, 6, 8, 13, 14, 22}};

static Vint maxi = 3;

static Vint eassoc[12][3] = {{1, 0, 2},  {0, 3, 4},  {0, 5, 6}, {7, 0, 8}, {9, 0, 0}, {0, 0, 10},
                             {0, 11, 0}, {0, 0, 12}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}};

/*----------------------------------------------------------------------
    Generate a Quadratic 3D Volume Tet Mesh and Perform Successive Refinements
----------------------------------------------------------------------*/
int
main()
{
    msh_TetMesh* tetmesh;
    vis_Connect *connect1, *connect2, *connecttmp;
    vis_State* state;
    vis_GridFun* gf;

    Vint i, j, iref;
    Vint numel;
    Vdouble size;
    Vchar buffer[80];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create tetmesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, 26, 12, maxi);

    /* create connect object */
    connect1 = vis_ConnectBegin();
    vis_ConnectPre(connect1, SYS_DOUBLE);
    vis_ConnectDef(connect1, 0, 0);

    /* define nodes */
    for (i = 0; i < 26; i++) {
        msh_TetMeshSetNode(tetmesh, i + 1, coords[i]);
        /* define a node association on every corner vertex */
        if (i < 8) {
            msh_TetMeshSetNodeAssoc(tetmesh, VIS_GEOVERT, i + 1, i + 1);
        }
    }
    /* define tris */
    for (i = 0; i < 12; i++) {
        msh_TetMeshSetTri(tetmesh, i + 1, tris[i]);
        for (j = 0; j < 3; ++j) {
            if (eassoc[i][j]) {
                msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOEDGE, i + 1, SYS_EDGE, j + 1, eassoc[i][j]);
            }
        }
    }
    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam44e.tet");
    /* generate */
    msh_TetMeshGenerate(tetmesh, connect1);
    vis_ConnectNumber(connect1, SYS_ELEM, &numel);
    printf("Generated original mesh\n");
    vis_ConnectWrite(connect1, SYS_NASTRAN_BULKDATA, "exam44e_orig.bdf");

    /* create state object with new sizes and refine */
    gf = vis_GridFunBegin();
    state = vis_StateBegin();
    vis_StateSetObject(state, VIS_GRIDFUN, (Vobject*)gf);

    connect2 = vis_ConnectBegin();
    vis_ConnectPre(connect2, SYS_DOUBLE);
    size = 1.;
    for (iref = 0; iref < 3; ++iref) {
        size *= 0.5;
        vis_ConnectGridFun(connect1, gf);
        vis_StatePre(state, SYS_DOUBLE);
        vis_StateDef(state, numel, SYS_ELEM, SYS_NONE, VIS_SCALAR);
        for (i = 1; i <= numel; i++) {
            vis_StateSetDatadv(state, i, &size);
        }
        /* refine mesh */
        vis_ConnectDef(connect2, 0, 0);
        msh_TetMeshRefine(tetmesh, state, connect2);
        printf("Generated refined mesh %d\n", iref + 1);
        sprintf(buffer, "exam44e_ref%d.bdf", iref + 1);
        vis_ConnectWrite(connect2, SYS_NASTRAN_BULKDATA, buffer);
        vis_ConnectNumber(connect2, SYS_ELEM, &numel);

        /* swap connect objects */
        connecttmp = connect1;
        connect1 = connect2;
        connect2 = connecttmp;
    }

    /* end objects */
    vis_ConnectEnd(connect1);
    vis_ConnectEnd(connect2);
    msh_TetMeshEnd(tetmesh);
    vis_StateEnd(state);
    vis_GridFunEnd(gf);
    return 0;
}

6.21. Example 45, Generate Surface Mesh and 3D Volume Tet Mesh

This example illustrates generating tetrahedral meshes from a simple geometry. The geometry is a straight forward unit cube but it illustrates the complete process of generating a surface triangulation with a SurfMesh object using the geometry and associated geometry tags, passing this boundary triangulation to a TetMesh object and generating a final tetrahedral mesh with all the generated nodes and tetrahedra associated to the original geometry features.

The geometry points are defined and assigned the VIS_GEOVERT point association type. The value of the association is set to the point id. The nodes in the surface mesh generated by SurfMesh and the nodes in the final volume mesh generated by TetMesh which lie on these geometry points will have the VIS_GEOVERT node association type set to the original defined point id value.

The triangle edges are assigned VIS_GEOEDGE element entity association values set to the edge ids in the range 1-12. The triangle faces are assigned the VIS_GEOFACE element entity association values set to the face ids in the range 1-6. Finally the triangles are assigned VIS_GEOBODY element association values of 1, since there is only a single body. All tetrahedra edges lying on geometry edges will be given VIS_GEOEDGE element entity associations, all tetrahedra faces lying on geometry faces will be given VIS_GEOFACE element entity associations, and all tetrahedra will be given VIS_GEOBODY element associations.

The parameter, maxi is set to 3 to generate parabolic elements throughout. Setting the parameter to 2 will generate linear elements. The final output of the TetMesh object which is placed in a Connect object is both printed and written to a NASTRAN bulk data file.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"
/*
                 3----3-----4
                /.         /|
               / 4        / 2
             12  .      11  |
             /   1....1./...2
            7----7-----8   /
            |  9       |  10
            8 .        6 /
            |.         |/
            5----5-----6

           faces are ordered 1(-z), 2(+z) 3(-y) 4(+y) 5(-x) 6(+x)

           y
           |
           +-- x
          /
         z
*/

/* a simple unit cube */
static Vdouble xc[8][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.},
                           {0., 0., 1.}, {1., 0., 1.}, {0., 1., 1.}, {1., 1., 1.}};

static Vint tris[12][3] = {{1, 5, 3}, {3, 5, 7}, {4, 6, 2}, {8, 6, 4}, {1, 2, 5}, {5, 2, 6},
                           {7, 4, 3}, {8, 4, 7}, {1, 4, 2}, {4, 1, 3}, {5, 8, 7}, {5, 6, 8}};
static Vint tefl[12][3] = {{1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1}, {1, 0, 1}, {0, 1, 1},
                           {0, 1, 1}, {1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {0, 1, 1}, {1, 1, 0}};
/* edge associations */
static Vint ted[12][3] = {{9, 0, 4},  {0, 8, 12}, {0, 10, 2}, {1, 0, 11}, {1, 0, 9}, {0, 10, 5},
                          {0, 3, 12}, {11, 0, 7}, {0, 2, 1},  {0, 4, 3},  {0, 7, 8}, {5, 6, 0}};
/* face associations */
static Vint tfa[12] = {5, 5, 6, 6, 3, 3, 4, 4, 1, 1, 2, 2};

/* triangle normals */
static Vdouble norms[12][3][3] = {{{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}}, {{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}},
                                  {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},    {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},
                                  {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}}, {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}},
                                  {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},    {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},
                                  {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}}, {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}},
                                  {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}},    {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}}};

/*----------------------------------------------------------------------
                      Generate a Surface Mesh and Tet Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_SurfMesh* surfmesh;
    msh_TetMesh* tetmesh;
    vis_Connect *connectsurf, *connecttet;

    Vint i, j, k;
    Vint maxi;
    Vint numsurfpnts, numsurftris;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[10];
    Vint num, aid, aids[10];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create SurfMesh object */
    surfmesh = msh_SurfMeshBegin();
    msh_SurfMeshDef(surfmesh, 8, 12);

    /* create Connect object */
    connectsurf = vis_ConnectBegin();
    vis_ConnectPre(connectsurf, SYS_DOUBLE);
    vis_ConnectDef(connectsurf, 0, 0);

    /* generate parabolic elements */
    maxi = 3;
    /* define points */
    for (i = 0; i < 8; i++) {
        msh_SurfMeshSetPoint(surfmesh, i + 1, xc[i], 1);
        /* set VIS_GEOVERT association to the point id */
        msh_SurfMeshSetPointAssoc(surfmesh, VIS_GEOVERT, i + 1, i + 1);
    }
    /* define triangle connectivity */
    for (i = 0; i < 12; i++) {
        msh_SurfMeshSetTri(surfmesh, i + 1, tris[i], tefl[i]);
        msh_SurfMeshSetTriNorm(surfmesh, i + 1, norms[i]);
        /* set VIS_GEOEDGE and VIS_GEOFACE */
        for (j = 0; j < 3; j++) {
            msh_SurfMeshSetTriAssoc(surfmesh, VIS_GEOEDGE, i + 1, SYS_EDGE, j + 1, ted[i][j]);
        }
        msh_SurfMeshSetTriAssoc(surfmesh, VIS_GEOFACE, i + 1, SYS_FACE, 0, tfa[i]);
        /* set VIS_GEOBODY to 1, only a single body */
        msh_SurfMeshSetTriAssoc(surfmesh, VIS_GEOBODY, i + 1, SYS_ELEM, 0, 1);
    }

    /* set mesh parameters */
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_MAXI, maxi);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, .5);

    msh_SurfMeshWrite(surfmesh, SYS_ASCII, "exam45.srf");
    /* generate */
    msh_SurfMeshGenerate(surfmesh, connectsurf);
    if (msh_SurfMeshError(surfmesh)) {
        printf("surf mesh generation error\n");
        exit(1);
    }
    else {
        printf("surf mesh generation complete\n");
    }
    /* print number of surface nodes and elements */
    vis_ConnectNumber(connectsurf, SYS_NODE, &numsurfpnts);
    vis_ConnectNumber(connectsurf, SYS_ELEM, &numsurftris);
    printf(" Number of nodes= %d\n", numsurfpnts);
    printf(" Number of tris = %d\n", numsurftris);

    /* write generated tet mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connectsurf, SYS_NASTRAN_BULKDATA, "exam45_srf.bdf");

    /* create TetMesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, numsurfpnts, numsurftris, maxi);

    /* Fill TetMesh object with generated points on surface */
    msh_TetMeshConnect(tetmesh, connectsurf);

    /* set mesh parameters */
    msh_TetMeshSetParami(tetmesh, VIS_MESH_MAXI, maxi);
    msh_TetMeshSetParamd(tetmesh, VIS_MESH_EDGELENGTH, .5);

    /* create Connect object to hold generated tet mesh */
    connecttet = vis_ConnectBegin();
    vis_ConnectPre(connecttet, SYS_DOUBLE);
    vis_ConnectDef(connecttet, 0, 0);

    /* generate mesh */
    msh_TetMeshGenerate(tetmesh, connecttet);
    if (msh_TetMeshError(tetmesh)) {
        printf("tet mesh generation error\n");
        exit(1);
    }
    else {
        printf("tet mesh generation complete\n");
    }
    /* print number of surface nodes and elements */
    vis_ConnectNumber(connecttet, SYS_NODE, &numnp);
    vis_ConnectNumber(connecttet, SYS_ELEM, &numel);
    printf(" Number of nodes= %d\n", numnp);
    printf(" Number of tets = %d\n", numel);

    /* write generated tet mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connecttet, SYS_NASTRAN_BULKDATA, "exam45_tet.bdf");

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connecttet, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %e, y= %e, z= %e\n", i, x[0], x[1], x[2]);
    }
    /* node associations */
    for (i = 1; i <= numnp; i++) {
        vis_ConnectAllNodeAssoc(connecttet, VIS_GEOVERT, i, &num, aids);
        for (j = 0; j < num; j++) {
            if (aids[j] == 0)
                continue;
            printf("id= %d, VIS_GEOVERT[%d]= %d\n", i, j, aids[j]);
        }
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connecttet, i, &nix, ix);
        printf("id= %d  ix=", i);
        for (j = 0; j < nix; j++) {
            printf(" %d", ix[j]);
        }
        printf("\n");
    }
    /* element associations */
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemAssoc(connecttet, VIS_GEOBODY, 1, &i, &aid);
        if (aid == 0)
            continue;
        printf("id= %d, VIS_GEOBODY= %d\n", i, aid);
    }
    /* element entity associations */
    for (i = 1; i <= numel; i++) {
        /* tet edges */
        for (j = 1; j <= 6; j++) {
            vis_ConnectAllElemEntAssoc(connecttet, VIS_GEOEDGE, SYS_EDGE, i, j, &num, aids);
            for (k = 0; k < num; k++) {
                if (aids[k] == 0)
                    continue;
                printf("id= %d edge= %d, VIS_GEOEDGE[%d]= %d\n", i, j, k, aids[k]);
            }
        }
        /* tet faces */
        for (j = 1; j <= 4; j++) {
            vis_ConnectElemEntAssoc(connecttet, VIS_GEOFACE, SYS_FACE, i, j, &aid);
            if (aid == 0)
                continue;
            printf("id= %d face= %d, VIS_GEOFACE= %d\n", i, j, aid);
        }
    }
    /* end objects */
    vis_ConnectEnd(connectsurf);
    vis_ConnectEnd(connecttet);
    msh_SurfMeshEnd(surfmesh);
    msh_TetMeshEnd(tetmesh);
    return 0;
}

6.22. Example 45a, Generate Surface and 3D Volume Tet Mesh with Sizing

This example illustrates generating a tetrahedral mesh from a non-manifold, rectilinear solid containing an internal face. Node based and region based sizing is applied. A default edge length of .2 is applied for the surface mesh. In region 1, an unconnected node is inserted in the center to define a size of .02 using msh_TetMeshSetNodeSizing(). In region 2 a region size of .05 is applied using msh_TetMeshSetTriSizing().

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/* a cube with interior surface */
/* point coordinates */
static Vdouble xc[12][3] = {{0., 0., 0.}, {1., 0., 0.}, {0., 1., 0.}, {1., 1., 0.}, {0., 0., 1.}, {1., 0., 1.},
                            {0., 1., 1.}, {1., 1., 1.}, {0., 0., 2.}, {1., 0., 2.}, {0., 1., 2.}, {1., 1., 2.}};

/* triangle connectivity */
static Vint tris[22][3] = {{1, 5, 3}, {3, 5, 7},  {4, 6, 2},  {8, 6, 4},   {1, 2, 5},   {5, 2, 6},  {7, 4, 3},  {8, 4, 7},
                           {1, 4, 2}, {4, 1, 3},  {5, 8, 7},  {5, 6, 8},   {5, 9, 7},   {7, 9, 11}, {8, 10, 6}, {12, 10, 8},
                           {5, 6, 9}, {9, 6, 10}, {11, 8, 7}, {12, 8, 11}, {9, 12, 11}, {9, 10, 12}};
/* triangle right and left element associations */
static Vint ta[22][2] = {
/* outer boundary and non-manifold internal surface */
{1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 0}, {1, 2},
{1, 2}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}, {2, 0}};

/* triangle edge flags */
static Vint tefl[22][3] = {{1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1}, {1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1},
                           {0, 1, 1}, {0, 1, 1}, {0, 1, 1}, {1, 1, 0}, {1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1},
                           {1, 0, 1}, {0, 1, 1}, {0, 1, 1}, {1, 0, 1}, {0, 1, 1}, {1, 1, 0}};

/* triangle normals */
static Vdouble norms[22][3][3] = {{{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}}, {{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}},
                                  {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},    {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},
                                  {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}}, {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}},
                                  {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},    {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},
                                  {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}}, {{0., 0., -1.}, {0., 0., -1.}, {0., 0., -1.}},
                                  {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}},    {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}},

                                  {{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}}, {{-1., 0., 0.}, {-1., 0., 0.}, {-1., 0., 0.}},
                                  {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},    {{1., 0., 0.}, {1., 0., 0.}, {1., 0., 0.}},
                                  {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}}, {{0., -1., 0.}, {0., -1., 0.}, {0., -1., 0.}},
                                  {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},    {{0., 1., 0.}, {0., 1., 0.}, {0., 1., 0.}},
                                  {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}},    {{0., 0., 1.}, {0., 0., 1.}, {0., 0., 1.}}};

/*----------------------------------------------------------------------
                      Generate Surface and 3D Volume Tet Mesh with Sizing
----------------------------------------------------------------------*/
int
main()
{
    msh_SurfMesh* surfmesh;
    msh_TetMesh* tetmesh;
    vis_Connect *connectsurf, *connecttet;

    Vint i;
    Vint numsurfpnts, numsurftris;
    Vint numnp, numel;
    Vdouble x[3];
    Vint nix, ix[3], aid;

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create SurfMesh object */
    surfmesh = msh_SurfMeshBegin();
    msh_SurfMeshDef(surfmesh, 12, 22);

    /* create Connect object */
    connectsurf = vis_ConnectBegin();
    vis_ConnectPre(connectsurf, SYS_DOUBLE);
    vis_ConnectDef(connectsurf, 0, 0);

    /* define points */
    for (i = 0; i < 12; i++) {
        msh_SurfMeshSetPoint(surfmesh, i + 1, xc[i], 1);
    }
    /* define triangle connectivity */
    for (i = 0; i < 22; i++) {
        msh_SurfMeshSetTri(surfmesh, i + 1, tris[i], tefl[i]);
        msh_SurfMeshSetTriNorm(surfmesh, i + 1, norms[i]);
        msh_SurfMeshSetTriBack(surfmesh, i + 1, ta[i][0], ta[i][1]);
        /* set property associations */
        if (ta[i][0]) {
            msh_SurfMeshSetTriAssoc(surfmesh, VIS_GEOBODY, i + 1, SYS_ELEM, 0, ta[i][0]);
            msh_SurfMeshSetTriAssoc(surfmesh, VIS_PROPID, i + 1, SYS_ELEM, 0, ta[i][0]);
        }
        if (ta[i][1]) {
            msh_SurfMeshSetTriAssoc(surfmesh, VIS_GEOBODY, i + 1, SYS_ELEM, -1, ta[i][1]);
            msh_SurfMeshSetTriAssoc(surfmesh, VIS_PROPID, i + 1, SYS_ELEM, -1, ta[i][1]);
        }
    }
    /* set mesh parameters */
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_MAXI, 2);
    msh_SurfMeshSetParami(surfmesh, SURFMESH_NONMANIFOLD, VIS_ON);
    msh_SurfMeshSetParami(surfmesh, SURFMESH_INTSURFBACK, VIS_ON);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, .2);

    /* write out file encapsulating meshing specifications */
    msh_SurfMeshWrite(surfmesh, SYS_ASCII, "exam45a.srf");
    /* generate */
    msh_SurfMeshGenerate(surfmesh, connectsurf);
    if (msh_SurfMeshError(surfmesh)) {
        printf("surf mesh generation error\n");
        exit(1);
    }
    else {
        printf("surf mesh generation complete\n");
    }
    msh_SurfMeshEnd(surfmesh);

    /* print number of surface nodes and elements */
    vis_ConnectNumber(connectsurf, SYS_NODE, &numsurfpnts);
    vis_ConnectNumber(connectsurf, SYS_ELEM, &numsurftris);
    printf(" Number of nodes= %d\n", numsurfpnts);
    printf(" Number of tris = %d\n", numsurftris);

    /* write generated surf mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connectsurf, SYS_NASTRAN_BULKDATA, "exam45a-srf.bdf");

    /* create TetMesh object */
    tetmesh = msh_TetMeshBegin();
    msh_TetMeshDef(tetmesh, numsurfpnts + 1, numsurftris, 2);

    /* Fill TetMesh object with generated points on surface */
    for (i = 1; i <= numsurfpnts; i++) {
        vis_ConnectCoordsdv(connectsurf, 1, &i, (Vdouble(*)[3])x);
        msh_TetMeshSetNode(tetmesh, i, x);
    }

    /* Fill TetMesh object with generated triangles */
    for (i = 1; i <= numsurftris; i++) {
        vis_ConnectElemNode(connectsurf, i, &nix, ix);
        msh_TetMeshSetTri(tetmesh, i, ix);
        vis_ConnectElemAssoc(connectsurf, VIS_PROPID, 1, &i, &aid);
        if (aid) {
            msh_TetMeshSetTriAssoc(tetmesh, VIS_PROPID, i, SYS_ELEM, 0, aid);
        }
        vis_ConnectElemAssoc(connectsurf, VIS_GEOBODY, 1, &i, &aid);
        if (aid) {
            msh_TetMeshSetTriAssoc(tetmesh, VIS_GEOBODY, i, SYS_ELEM, 0, aid);
        }
    }
    /* set sizing */
    msh_TetMeshSetParamd(tetmesh, VIS_MESH_EDGELENGTH, .2);
    /* add sizing at unconnected node in center of region 1 */
    x[0] = .5;
    x[1] = .5;
    x[2] = .5;
    msh_TetMeshSetNode(tetmesh, numsurfpnts + 1, x);
    msh_TetMeshSetNodeSizing(tetmesh, numsurfpnts + 1, .02);
    msh_TetMeshSetParami(tetmesh, TETMESH_UNCONNECT, SYS_ON);
    /* set region sizing */
    for (i = 1; i <= numsurftris; i++) {
        vis_ConnectElemAssoc(connectsurf, VIS_GEOBODY, 1, &i, &aid);
        if (aid == 2) {
            msh_TetMeshSetTriSizing(tetmesh, i, .05);
        }
    }

    vis_ConnectEnd(connectsurf);

    /* create Connect object to hold generated tet mesh */
    connecttet = vis_ConnectBegin();
    vis_ConnectPre(connecttet, SYS_DOUBLE);
    vis_ConnectDef(connecttet, 0, 0);

    /* write out file encapsulating meshing specifications */
    msh_TetMeshWrite(tetmesh, SYS_ASCII, "exam45a.tet");
    /* generate mesh */
    msh_TetMeshGenerate(tetmesh, connecttet);
    if (msh_TetMeshError(tetmesh)) {
        printf("tet mesh generation error\n");
        exit(1);
    }
    else {
        printf("tet mesh generation complete\n");
    }
    msh_TetMeshEnd(tetmesh);
    /* print number of surface nodes and elements */
    vis_ConnectNumber(connecttet, SYS_NODE, &numnp);
    vis_ConnectNumber(connecttet, SYS_ELEM, &numel);
    printf(" Number of nodes= %d\n", numnp);
    printf(" Number of elems= %d\n", numel);

    /* write generated tet mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connecttet, SYS_NASTRAN_BULKDATA, "exam45a-tet.bdf");

    vis_ConnectEnd(connecttet);
    return 0;
}

6.23. Example 46, Generate 3D Curve Mesh Using CurvMesh

This example illustrates generating parabolic line elements along a curve defined by 4 input line segments. The input curve consists of a straight line segment .75 long along the positive y-axis followed by a 90 degree curve of radius .25 followed by a straight line segment 1.75 long along the positive x-axis. The number of input points and the number and order of input line segments are specified by calling msh_CurvMeshDef(), then each point and line is input using msh_CurvMeshSetPoint() and msh_CurvMeshSetLine() respectively. The input point numbers 2 and 4 are marked as preserved. This ensures that a node is generated at this location in the output mesh. All free points (those points connected to only one input line segment) are automatically preserved. The function msh_CurvMeshSetLineTang() sets the line edgepoint tangent vectors. These tangents are used to construct a cubic geometry along the line.

A line association is defined along input line 1. As a result all nodes generated along input line 1 will be given the the specified association. In a similar manner, a point association is defined at point 5. The node generated under the point will be given the specified association.

The density of the generated mesh and the order of line elements to generate are specified using msh_CurvMeshSetParamd() and msh_CurvMeshSetParami(). In this case a target edge length of .4 is specified with no line element spanning more than 30. degrees. The size of the generated line elements is the minimum of these two size criteria.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble xc[5][3] = {{0., 0., 0.}, {0., .75, 0.}, {.0732233, .9267767, 0.}, {.25, 1., 0.}, {2., 1., 0.}};

static Vint pc[5] = {0, 1, 0, 1, 0};

static Vint lix[4][2] = {{1, 2}, {2, 3}, {3, 4}, {4, 5}};
static Vdouble vex[4][2][3] = {
{{0., 1., 0.}, {0., 1., 0.}}, {{0., 1., 0.}, {.707, .707, 0.}}, {{.707, .707, 0.}, {1., 0., 0.}}, {{1., 0., 0.}, {1., 0., 0.}}};

/*----------------------------------------------------------------------
                      Generate a Boundary Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_CurvMesh* curvmesh;
    vis_Connect* connect;

    Vint i;
    Vint numnp, numel;
    Vdouble x[3];
    Vint aid, aid1;
    Vint nix, ix[3];

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create CurvMesh object */
    curvmesh = msh_CurvMeshBegin();
    msh_CurvMeshDef(curvmesh, 5, 4);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectDef(connect, 0, 0);

    /* define points */
    for (i = 0; i < 5; i++) {
        msh_CurvMeshSetPoint(curvmesh, i + 1, xc[i], pc[i]);
    }
    /* define parabolic line segment connectivity */
    for (i = 0; i < 4; i++) {
        msh_CurvMeshSetLine(curvmesh, i + 1, lix[i]);
        msh_CurvMeshSetLineTang(curvmesh, i + 1, vex[i]);
    }
    /* set assoc of 100 at line 1 */
    msh_CurvMeshSetLineAssoc(curvmesh, VIS_MISCID, 1, 100);

    /* set assoc of 200 at point 5 */
    msh_CurvMeshSetPointAssoc(curvmesh, VIS_MISCID1, 5, 200);

    /* set mesh parameters */
    msh_CurvMeshSetParami(curvmesh, VIS_MESH_MAXI, 3);
    msh_CurvMeshSetParamd(curvmesh, VIS_MESH_SPANANGLE, 30.);
    msh_CurvMeshSetParamd(curvmesh, VIS_MESH_EDGELENGTH, .4);

    /* generate */
    msh_CurvMeshGenerate(curvmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    vis_ConnectWrite(connect, SYS_SDRC_UNIVERSAL, "exam46.unv");

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        vis_ConnectNodeAssoc(connect, VIS_MISCID, 1, &i, &aid);
        vis_ConnectNodeAssoc(connect, VIS_MISCID1, 1, &i, &aid1);
        printf("id= %d  x= %f, y= %f, z= %f,  aid= %d, aid1= %d\n", i, x[0], x[1], x[2], aid, aid1);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d\n", i, ix[0], ix[1], ix[2]);
    }
    /* end objects */
    vis_ConnectEnd(connect);
    msh_CurvMeshEnd(curvmesh);
    return 0;
}

6.24. Example 47, Generate 3D Surface Mesh Using SurfMesh

This example illustrates generating 3 node triangular elements on a curved surface defined by 12 points and 10 triangular facets. The input surface represents a simple folded plate. This geometry illustrates the ability of the SurfMesh module to automatically refine the generated mesh in areas of high curvature and transition the element size to a target size in flat areas.

The number of input points and the number and order of input triangles are specified by calling msh_SurfMeshDef(), then each point and triangle is input using msh_SurfMeshSetPoint() and msh_SurfMeshSetTri() respectively. In general new points are placed on the surface defined by the triangular facets and any defined point may be moved unless it is explicitly preserved. Note that the four corner input points are marked as preserved in the function msh_SurfMeshSetPoint() so that they do not move.

The density of the generated mesh and the order of elements to generate are specified using msh_SurfMeshSetParamd() and msh_SurfMeshSetParami(). The dimensions of the plate are approximately 1. by 2. and it has been folded in the middle through 90 degrees with a radius of curvature of .1 . A target edge length of .2 is specified with no element spanning more than 30. degrees. The size of the generated elements is the minimum of these two size criteria. A minimum edge length of .02 is specified. This ensures that no unnecessarily small triangles are produced. The call msh_SurfMeshGenerate() generates nodes and triangular elements and enters them into the Connect object.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

static Vdouble xc[12][3] = {
{0.000e+00, 0.000e+00, 0.000e+00}, {0.000e+00, 9.000e-01, 0.000e+00}, {1.333e-02, 9.500e-01, 0.000e+00},
{5.000e-02, 9.866e-01, 0.000e+00}, {1.000e-01, 1.000e+00, 0.000e+00}, {1.000e+00, 1.000e+00, 0.000e+00},
{0.000e+00, 0.000e+00, 1.000e+00}, {0.000e+00, 9.000e-01, 1.000e+00}, {1.333e-02, 9.500e-01, 1.000e+00},
{5.000e-02, 9.866e-01, 1.000e+00}, {1.000e-01, 1.000e+00, 1.000e+00}, {1.000e+00, 1.000e+00, 1.000e+00}};

static Vint pc[12] = {1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1};

static Vint tix[10][3] = {{1, 2, 7},  {2, 8, 7},  {2, 3, 8},   {3, 9, 8},  {3, 4, 9},
                          {4, 10, 9}, {4, 5, 10}, {5, 11, 10}, {5, 6, 11}, {6, 12, 11}};

static Vint efl[3] = {0, 0, 0};

/*----------------------------------------------------------------------
                      Generate a Surface Mesh
----------------------------------------------------------------------*/
int
main()
{
    msh_SurfMesh* surfmesh = NULL;
    vis_Connect* connect = NULL;

    Vint i = 0;
    Vint numnp = 0;
    Vint numel = 0;
    Vdouble x[3] = {0};
    Vint nix = 0;
    Vint ix[3] = {0};

    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    /* create SurfMesh object */
    surfmesh = msh_SurfMeshBegin();
    msh_SurfMeshDef(surfmesh, 12, 10);

    /* create Connect object */
    connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* define points */
    for (i = 0; i < 12; i++) {
        msh_SurfMeshSetPoint(surfmesh, i + 1, xc[i], pc[i]);
    }
    /* define triangle connectivity */
    for (i = 0; i < 10; i++) {
        msh_SurfMeshSetTri(surfmesh, i + 1, tix[i], efl);
    }

    /* set mesh parameters */
    msh_SurfMeshSetParami(surfmesh, VIS_MESH_MAXI, 2);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_SPANANGLE, 30.);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_MINEDGELENGTH, .02);
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, .2);

    /* write out file encapsulating meshing specifications */
    msh_SurfMeshWrite(surfmesh, SYS_ASCII, "exam47.srf");
    /* generate */
    msh_SurfMeshGenerate(surfmesh, connect);

    /* print generated nodes and elements */
    vis_ConnectNumber(connect, SYS_NODE, &numnp);
    vis_ConnectNumber(connect, SYS_ELEM, &numel);
    printf("numnp= %d, numel= %d\n", numnp, numel);

    /* print node information */
    printf("Node information\n");
    for (i = 1; i <= numnp; i++) {
        vis_ConnectCoordsdv(connect, 1, &i, (Vdouble(*)[3])x);
        printf("id= %d  x= %f, y= %f, z= %f\n", i, x[0], x[1], x[2]);
    }
    /* print element information */
    printf("Element information\n");
    for (i = 1; i <= numel; i++) {
        vis_ConnectElemNode(connect, i, &nix, ix);
        printf("id= %d  ix= %d %d %d\n", i, ix[0], ix[1], ix[2]);
    }
    /* write generated mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "example47_SurfMeshGenerate.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_SurfMeshEnd(surfmesh);
    return 0;
}

6.25. Example 47 [Parallel], Generate 3D Surface Mesh Using Parallel SurfMesh

This example illustrates generating surface triangle elements in a complex geometry. The input file is a .srf file read by the SurfMesh

The target is to perform this in parallel. To define the size of the thread pool to be used, please use msh_SurfMeshSetNumThreads(). To generate and execute please use msh_SurfMeshGenerate().

The algorithm for the parallel strategy consists of the following steps:

Step 1: Gather all elements that are not complying with the surfmesh parameters. Step 2: Call ParMetis to generate a given number of partitions. Step 3: Per Partition, call surfMesh blocking the boundary edges. Step 4: Per partition, identify two element groups. Compliant elements and Non-Compliant elements. Step 5: Use the non-compliant elements as input for the next iteration.

In this example, the algorithm will:

ITERATION 1 Split the domain into 16 partitions, we consume a threadpool using 4 threads to work dynamically. We generated 16 partition meshes. Those meshes will have compliant and non-compliant elements.

ITERATION 2 Taking only as input the non-compliant elements, we partition that element group to generate newly 16 partitions. We block the new partition-boundaries and we generate surface meshes. From that second iteration, we split two element groups: the compliant and non-compliant elements for next iteration.

ITERATION 3 Finally, we take the the non-compliant elements group from the previous iteration. Since this is the last iteration, we call surfMesh in sequential.

ITERATION 1,2,3 + FINAL MERGE After 3 Iterations, we gather the element group of compliant elements of each iteration and we obtain a perfect puzzle!

The partition strategy can be changed using msh_SurfMeshSetPartitionStrategy(). Refer to that method for further details.

The call msh_SurfMesh_EXPERIMENTAL_Generate() generates nodes and triangular elements and enters them into the Connect object.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"
#include <vector>
#include <string>

/*----------------------------------------------------------------------
        Generate a Surface Mesh in Parallel
----------------------------------------------------------------------*/
int
main(int argc, char** argv)
{
    vsy_LicenseValidate(CEETRON_SAM_LICENSE);

    std::string filename = argv[1]; // Expecting "PATH_to/walek.srf";
    constexpr Vint threadCount = 4;
    constexpr Vdouble targetEdgeLength = 0.6;
    std::vector<Vint> partitionStrategy = {16, 16};

    /* create Connect object */
    vis_Connect* connect = vis_ConnectBegin();
    vis_ConnectPre(connect, SYS_DOUBLE);
    vis_ConnectDef(connect, 0, 0);

    /* create SurfMesh object */
    msh_SurfMesh* surfmesh = msh_SurfMeshBegin();

    msh_SurfMeshRead(surfmesh, SYS_ASCII, filename.c_str());
    Vint errorFlag = msh_SurfMeshError(surfmesh);
    if (errorFlag) {
        printf("A problem occurred during reading of input file\n");
        exit(errorFlag);
    }

    /* set mesh parameters */
    msh_SurfMeshSetParamd(surfmesh, VIS_MESH_EDGELENGTH, targetEdgeLength);

    /* Define the parallel strategy. Set the number of partitions per iteration of the algorithm.*/
    msh_SurfMeshSetPartitionStrategy(surfmesh, static_cast<Vint>(partitionStrategy.size()), partitionStrategy.data());

    /* Define the number of threads to use during the parallel computation and call Experimental Generate */
    msh_SurfMeshSetNumThreads(surfmesh, threadCount);
    msh_SurfMeshGenerate(surfmesh, connect);

    errorFlag = msh_SurfMeshError(surfmesh);

    if (errorFlag) {
        printf("A problem occurred during execution of parallel surfmesher\n");
    }

    /*Optionally, print some recorded data available in a dictionary report*/
    vsy_DblDict* keyValueDict = nullptr;
    msh_SurfMeshGetExecutionReport(surfmesh, &keyValueDict);
    // vsy_DblDictPrint(keyValueDict);

    /* print output nodes and elements */
    Vint nodeCount = 0;
    Vint elementCount = 0;
    vis_ConnectNumber(connect, SYS_NODE, &nodeCount);
    vis_ConnectNumber(connect, SYS_ELEM, &elementCount);
    printf("nodeCount= %d, elementCount= %d\n", nodeCount, elementCount);

    // EXPECTED output
    /*
     * Percentage of Compliant Elements [%], v= 98.278854
     * Percentage of Non-Compliant Elements [%], v= 1.721146
     * Total Surfmesher Time [s], v= 58.630623
     * Iteration 3 Time[s], v= 14.753183
     * Iteration 2 Time[s], v= 4.554247
     * Iteration 1 Time[s], v= 38.485378
     * nodeCount= 923153, elementCount= 1846298
     */

    /* write generated mesh in NASTRAN bulk data format */
    vis_ConnectWrite(connect, SYS_NASTRAN_BULKDATA, "example47_SurfMeshParallel.bdf");

    /* end objects */
    vis_ConnectEnd(connect);
    msh_SurfMeshEnd(surfmesh);
    return 0;
}
../_images/vistools_exam47parallel.png

6.26. Example 47a, Generate 3D Surface Mesh with Associations

This example is an extension of Example 47 in which node and element associations are assigned to the input points, triangles and triangle edges to be propagated into the final mesh. Triangle vertex normals are also assigned. In addition, alterations are made to the meshing parameters to output 6 node parabolic triangles rather than 3 node linear triangles. Node association types VIS_MISCID and VIS_MISCID1 are used in this example, however any association type may be assigned. Node associations are assigned to points 1 and 7 using msh_SurfMeshSetPointAssoc(). Besides specifying the association, this function also marks the node to be preserved. After the triangles are defined, node associations are defined on triangle edges and faces. The triangle edge association on triangle 2 edge 1 creates a preserved edge. This preserved edge along with the “free” edges of triangle 1 and 2 (free edges are always preserved) bounds triangles 1 and 2. This means that node and element associations can be safely applied to the set of triangles 1 and 2 with predictable results. The end points of the preserved edge (triangle 1 edge 2) are automatically preserved as “topological” points. Any point connected to 1 preserved edge or 3 or more preserved edges is a topologically preserved point. The VIS_PROPID association is used to propagate property ids 1 and 2 to the two sets of triangles bounded by preserved edges.

The function msh_SurfMeshSetParami() is used to indicate that parabolic triangles are to be output. The mesh density and spanning angle parameters are doubled from those used in Example 47 where linear elements were generated. The call msh_SurfMeshGenerate() generates nodes and triangular elements and enters them into the Connect object. The resulting mesh is printed and also written to file exam47a.bdf in NASTRAN bulk data file format. Note that vis/exam/exam30ddev.cpp may be used to display the generated mesh saved in the file.

In the figure below, triangle and point numbers are displayed along with triangle vertex normals in orange. Preserved edges are black, free edges are red and interior, unpreserved edges are in white.

#include <stdio.h>
#include "sam/base/base.h"
#include "sam/msh/msh.h"
#include "sam/base/license.h"
#include "sam/CEETRON_SAM_license.h"

/* point coordinates */
static Vdouble xc[12][3] = {
{0.000e+00, 0.000e+00, 0.000e+00}, {0.000e+00, 9.000e-01, 0.000e+00}, {1.333e-02, 9.500e-01, 0.000e+00},
{5.000e-02, 9.866e-01, 0.000e+00}, {1.000e-01, 1.000e+00, 0.000e+00}, {1.000e+00, 1.000e+00, 0.000e+00},
{0.000e+00, 0.000e+00, 1.000e+00}, {0.000e+00, 9.000e-01, 1.000e+00}, {1.333e-02, 9.500e-01, 1.000e+00},
{5.000e-02, 9.866e-01, 1.000e+00}, {1.000e-01, 1.000e+00, 1.000e+00}, {1.000e+00, 1.000e+00, 1.000e+00}};
/* preserved point flags */
static Vint pc[12] = {1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1};
/* triangle connectivity */
static Vint tix[10][3] = {{1, 2, 7},  {2, 8, 7},  {2, 3, 8},   {3, 9, 8},  {3, 4, 9},
                          {4, 10, 9}, {4, 5, 10}, {5, 11, 10}, {5, 6, 11}, {6, 12, 11}};
/* triangle preserved edge flags */
static Vint efl[10][3] = {{0, 0, 0},<