Example a_backup/meshing_api: Tetrahedral meshing including boundary extraction and element creation with the Cmgui C API.

This example program demonstrates how to integrate third party meshing algorithms into software using the Cmgui library. It calls Cmgui APIs to extract boundary geometry for input to a meshing library (Netgen in this example), and invokes Cmgui element APIs to create the resulting mesh and coordinate field.

The program first creates a linear cube using the API, then draws an isosurface of a closed sphere inside it. (In a a real situation this could be generated from a segmented medical image.) This isosurface is exported to an STL file and Netgen APIs are used to read the STL file and generate a tetrahedral mesh within its boundaries. Finally, a Cmgui mesh consisting of nodes, elements and coordinate field is created using the Cmgui API from the Netgen output, and this is exported to file.

Note that the Cmgui API supports alternative approaches for extracting geometry without writing to an external file. These include iteration over elements and field evaluation to obtain exact surface boundaries. Surface graphics can also be converted to a mesh for processing in the same manner.

Screenshot of example a_backup/meshing_api


The comfile run by this example is as follows:

/* ***** BEGIN LICENSE BLOCK *****
 * Version: MPL 1.1/GPL 2.0/LGPL 2.1
 *
 * The contents of this file are subject to the Mozilla Public License Version
 * 1.1 (the "License"); you may not use this file except in compliance with
 * the License. You may obtain a copy of the License at
 * http://www.mozilla.org/MPL/
 *
 * Software distributed under the License is distributed on an "AS IS" basis,
 * WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License
 * for the specific language governing rights and limitations under the
 * License.
 *
 * The Original Code is cmgui meshing_api example.
 *
 * The Initial Developer of the Original Code is
 * Auckland Uniservices Ltd, Auckland, New Zealand.
 * Portions created by the Initial Developer are Copyright (C) 2010
 * the Initial Developer. All Rights Reserved.
 *
 * Contributor(s):
 *
 * Alternatively, the contents of this file may be used under the terms of
 * either the GNU General Public License Version 2 or later (the "GPL"), or
 * the GNU Lesser General Public License Version 2.1 or later (the "LGPL"),
 * in which case the provisions of the GPL or the LGPL are applicable instead
 * of those above. If you wish to allow use of your version of this file only
 * under the terms of either the GPL or the LGPL, and not to allow others to
 * use your version of this file under the terms of the MPL, indicate your
 * decision by deleting the provisions above and replace them with the notice
 * and other provisions required by the GPL or the LGPL. If you do not delete
 * the provisions above, a recipient may use your version of this file under
 * the terms of any one of the MPL, the GPL or the LGPL.
 *
 * ***** END LICENSE BLOCK ***** */

extern "C"
{
#include "api/cmiss_context.h"
#include "api/cmiss_field.h"
#include "api/cmiss_field_finite_element.h"
#include "api/cmiss_field_module.h"
#include "api/cmiss_element.h"
#include "api/cmiss_graphic.h"
#include "api/cmiss_graphics_module.h"
#include "api/cmiss_node.h"
#include "api/cmiss_region.h"
#include "api/cmiss_rendition.h"
#include "api/cmiss_status.h"
#include "api/cmiss_tessellation.h"
}
#include <iostream>

namespace nglib {
#include <nglib.h>
}
using namespace nglib;

using namespace std;


/*****************************************************************************//**
 * The following function creates a single element linear lagrange cube with 8 nodes
 * using the finite element APIs provided in cmgui.
 */
void create_cube_element(Cmiss_region_id root_region)
{
	Cmiss_region_id region = Cmiss_region_create_child(root_region, "cube");
	Cmiss_field_module_id field_module = Cmiss_region_get_field_module(region);
	Cmiss_field_module_begin_change(field_module);

	Cmiss_field_id coordinates_field =
		Cmiss_field_module_create_finite_element(field_module, /*number_of_components*/3);
	Cmiss_field_set_name(coordinates_field, "coordinates");
	Cmiss_field_set_attribute_integer(coordinates_field, CMISS_FIELD_ATTRIBUTE_IS_COORDINATE, 1);
	Cmiss_field_set_attribute_integer(coordinates_field, CMISS_FIELD_ATTRIBUTE_IS_MANAGED, 1);

	/* create nodes & set field parameters */

	Cmiss_nodeset_id nodeset = Cmiss_field_module_find_nodeset_by_name(field_module, "cmiss_nodes");
	Cmiss_node_template_id node_template1 = Cmiss_nodeset_create_node_template(nodeset);
	Cmiss_node_template_define_field(node_template1, coordinates_field);
	double node_coordinates[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 }
	};
	Cmiss_field_cache_id field_cache = Cmiss_field_module_create_cache(field_module);
	for (int i = 0; i < 8; i++)
	{
		Cmiss_node_id node = Cmiss_nodeset_create_node(nodeset, i+1, node_template1);
		Cmiss_field_cache_set_node(field_cache, node);
		Cmiss_field_assign_real(coordinates_field, field_cache, /*number_of_values*/3, node_coordinates[i]);
		Cmiss_node_destroy(&node);
	}
	Cmiss_field_cache_destroy(&field_cache);
	Cmiss_node_template_destroy(&node_template1);

	/* define element field interpolation */

	Cmiss_mesh_id mesh = Cmiss_field_module_find_mesh_by_dimension(field_module, /*dimension*/3);
	Cmiss_element_template_id element_template = Cmiss_mesh_create_element_template(mesh);
	Cmiss_element_template_set_shape_type(element_template, CMISS_ELEMENT_SHAPE_CUBE);
	Cmiss_element_template_set_number_of_nodes(element_template, 8);
	Cmiss_element_basis_id trilinear_basis = Cmiss_field_module_create_element_basis(
		field_module, /*dimension*/3, CMISS_BASIS_FUNCTION_LINEAR_LAGRANGE);
	const int cube_local_node_indexes[8] = { 1, 2, 3, 4, 5, 6, 7, 8};
	Cmiss_element_template_define_field_simple_nodal(element_template, coordinates_field,
		/*component_number*/-1, trilinear_basis, /*number_of_nodes*/8, cube_local_node_indexes);
	Cmiss_element_basis_destroy(&trilinear_basis);

	/* create element */
	for (int i = 1; i <= 8; i++)
	{
		Cmiss_node_id node = Cmiss_nodeset_find_node_by_identifier(nodeset, i);
		Cmiss_element_template_set_node(element_template, i, node);
		Cmiss_node_destroy(&node);
	}
	Cmiss_mesh_define_element(mesh, -1, element_template);
	Cmiss_element_template_destroy(&element_template);

	Cmiss_mesh_destroy(&mesh);
	Cmiss_nodeset_destroy(&nodeset);
	Cmiss_field_destroy(&coordinates_field);
	Cmiss_field_module_end_change(field_module);
	Cmiss_field_module_destroy(&field_module);
	Cmiss_region_destroy(®ion);
}

/*****************************************************************************//**
 * Create an spherical isosurface inside the generated cubic element.
 * This function uses a few execute commands. These commands may be replaced in the
 * future.
 */
void create_isosurface(Cmiss_context_id context)
{
	Cmiss_region_id root_region = Cmiss_context_get_default_region(context);
	Cmiss_region_id child = Cmiss_region_find_child_by_name(root_region, "cube");
	Cmiss_context_execute_command(context, "gfx define field cube/two constant 2");
	Cmiss_context_execute_command(context, "gfx define field cube/sum multiply fields two coordinates");
	Cmiss_context_execute_command(context, "gfx define field cube/one constant -1");
	Cmiss_context_execute_command(context, "gfx define field cube/final add fields sum one");
	Cmiss_context_execute_command(context, "gfx define field cube/sphere magnitude field final");
	Cmiss_graphics_module_id graphics_module = Cmiss_context_get_default_graphics_module(context);
	Cmiss_graphics_module_enable_renditions(graphics_module, root_region);
	Cmiss_rendition_id rendition = Cmiss_graphics_module_get_rendition(graphics_module, child);
	Cmiss_rendition_execute_command(rendition, "general clear");
	Cmiss_rendition_execute_command(rendition,
		"iso_surfaces iso_scalar sphere iso_values 0.5 coordinate coordinates use_elements select_on material default selected_material default_selected");
	Cmiss_graphic_id graphic = Cmiss_rendition_get_first_graphic(rendition);
	Cmiss_tessellation_id tessellation = Cmiss_graphics_module_create_tessellation(graphics_module);
	const int minimum_divisions[3] = {16, 16, 16};
	Cmiss_tessellation_set_minimum_divisions(tessellation, 3, &(minimum_divisions[0]));
	Cmiss_graphic_set_tessellation(graphic, tessellation);
	Cmiss_tessellation_destroy(&tessellation);
	Cmiss_graphic_destroy(&graphic);
	Cmiss_rendition_destroy(&rendition);
	Cmiss_graphics_module_destroy(&graphics_module);
	Cmiss_region_destroy(&child);
	Cmiss_region_destroy(&root_region);
}

/*****************************************************************************//**
 * Export the isosurface created previously into a stl file. This will create a
 * file in the local hard disk. But it can probably be stored in the memory instead.
 */
void export_stl(Cmiss_context_id context, const char *file_name)
{
	char command[200];
	sprintf(command, "gfx export stl file %s scene default", file_name);
	Cmiss_context_execute_command(context, command);
}

/*****************************************************************************//**
 * Convert a mesh from netgen into linear simplex model in cmgui using finite
 * element APIs. This function will export the element and node into an exregion
 * file for testing purpose.
 */
void mesh_to_simplex_conversion(Cmiss_region_id root_region, Ng_Mesh *ng_mesh)
{
	const int number_of_nodes = Ng_GetNP(ng_mesh);

	Cmiss_region_id region = Cmiss_region_create_child(root_region, "final_mesh");
	Cmiss_field_module_id field_module = Cmiss_region_get_field_module(region);
	Cmiss_field_module_begin_change(field_module);

	Cmiss_field_id coordinates_field =
		Cmiss_field_module_create_finite_element(field_module, /*number_of_components*/3);
	Cmiss_field_set_name(coordinates_field, "coordinates");
	Cmiss_field_set_attribute_integer(coordinates_field, CMISS_FIELD_ATTRIBUTE_IS_COORDINATE, 1);
	Cmiss_field_set_attribute_integer(coordinates_field, CMISS_FIELD_ATTRIBUTE_IS_MANAGED, 1);

	/* create nodes & set field parameters */
	Cmiss_nodeset_id nodeset = Cmiss_field_module_find_nodeset_by_name(field_module, "cmiss_nodes");
	Cmiss_node_template_id node_template1 = Cmiss_nodeset_create_node_template(nodeset);
	Cmiss_node_template_define_field(node_template1, coordinates_field);
	double coor_temp[3];
	Cmiss_field_cache_id field_cache = Cmiss_field_module_create_cache(field_module);
	for (int i = 0; i < number_of_nodes; i++)
	{
		Ng_GetPoint(ng_mesh, i+1/*index in netgen starts from 1*/, coor_temp);
		Cmiss_node_id node = Cmiss_nodeset_create_node(nodeset, i+1, node_template1);
		Cmiss_field_cache_set_node(field_cache, node);
		Cmiss_field_assign_real(coordinates_field, field_cache, /*number_of_values*/3, coor_temp);
		Cmiss_node_destroy(&node);
	}
	Cmiss_field_cache_destroy(&field_cache);
	Cmiss_node_template_destroy(&node_template1);

	/* define element field interpolation */
	Cmiss_mesh_id mesh = Cmiss_field_module_find_mesh_by_dimension(field_module, /*dimension*/3);
	Cmiss_element_template_id element_template = Cmiss_mesh_create_element_template(mesh);
	Cmiss_element_template_set_shape_type(element_template, CMISS_ELEMENT_SHAPE_TETRAHEDRON);
	Cmiss_element_template_set_number_of_nodes(element_template, 4);
	Cmiss_element_basis_id linear_tet_basis = Cmiss_field_module_create_element_basis(
		field_module, /*dimension*/3, CMISS_BASIS_FUNCTION_LINEAR_SIMPLEX);
	const int linear_tet_local_node_indexes[4] = { 1, 2, 3, 4 };
	Cmiss_element_template_define_field_simple_nodal(element_template, coordinates_field,
		/*component_number*/-1, linear_tet_basis, /*number_of_nodes*/4, linear_tet_local_node_indexes);
	Cmiss_element_basis_destroy(&linear_tet_basis);

	/* create element */
	const int number_of_elements = Ng_GetNE(ng_mesh);
	int nodal_idx[4];
	for (int i = 0 ; i < number_of_elements; i++)
	{
		Ng_GetVolumeElement(ng_mesh, i+1, nodal_idx);
		for (int k = 0; k < 4; k++)
		{
			Cmiss_node_id node = Cmiss_nodeset_find_node_by_identifier(nodeset, nodal_idx[k]);
			Cmiss_element_template_set_node(element_template, k+1, node);
			Cmiss_node_destroy(&node);
		}
		Cmiss_mesh_define_element(mesh, -1, element_template);
	}

	Cmiss_element_template_destroy(&element_template);
	Cmiss_mesh_destroy(&mesh);
	Cmiss_nodeset_destroy(&nodeset);
	Cmiss_field_destroy(&coordinates_field);
	Cmiss_field_module_end_change(field_module);
	Cmiss_field_module_destroy(&field_module);
	Cmiss_region_destroy(®ion);
}

/*****************************************************************************//**
 * This function read in a stl file and will be loaded into Netgen and meshed.
 * This funciton will also calls the mesh_to_simplex_conversion function which
 * will convert the mesh into simplex element using cmgui APIs.
 */
void netgen_mesh_generation(Cmiss_region_id root_region, const char *file_name)
{
	Ng_Init();
	Ng_Meshing_Parameters *mp=new Ng_Meshing_Parameters();
	Ng_STL_Geometry *geom=Ng_STL_LoadGeometry(file_name);
	Ng_Mesh *mesh = Ng_NewMesh();
	int return_code=Ng_STL_InitSTLGeometry(geom);
	return_code=Ng_STL_MakeEdges(geom, mesh, mp);
	if (return_code!=NG_OK)
	{
		Ng_Exit();
	}

	return_code=Ng_STL_GenerateSurfaceMesh(geom, mesh, mp);
	if (return_code!=NG_OK)
	{
		Ng_Exit();
	}

	return_code=Ng_GenerateVolumeMesh(mesh,mp);
	if (return_code!=NG_OK)
	{
		Ng_Exit();
	}
	mesh_to_simplex_conversion(root_region, mesh);
	Ng_DeleteMesh(mesh);
	delete mp;
	Ng_Exit ();
}

/*****************************************************************************//**
 * Using combinations of execute command and APIs function to create visualisation
 * of the final mesh(simplex elements).
 */
void show_final_mesh(Cmiss_context_id context)
{
	Cmiss_region_id root_region = Cmiss_context_get_default_region(context);
	Cmiss_context_execute_command(context, "gfx define faces egroup final_mesh");
	Cmiss_region_id child = Cmiss_region_find_child_by_name(root_region, "final_mesh");
	Cmiss_graphics_module_id graphics_module = Cmiss_context_get_default_graphics_module(context);
	Cmiss_rendition_id rendition = Cmiss_graphics_module_get_rendition(graphics_module, child);
	//Cmiss_rendition_execute_command(rendition, "surface coordinate coordinates");
	Cmiss_rendition_execute_command(rendition,"surfaces coordinate coordinates");
	Cmiss_context_execute_command(context, "gfx create window 1");
	Cmiss_region_id cube = Cmiss_region_find_child_by_name(root_region, "cube");
	Cmiss_rendition_id cube_rendition = Cmiss_graphics_module_get_rendition(graphics_module, cube);
	Cmiss_graphic_id graphic = Cmiss_rendition_get_first_graphic(cube_rendition);
	Cmiss_graphic_set_visibility_flag(graphic, 0);
	Cmiss_graphic_destroy(&graphic);
	Cmiss_rendition_destroy(&cube_rendition);
	Cmiss_rendition_destroy(&rendition);
	Cmiss_graphics_module_destroy(&graphics_module);
	Cmiss_region_destroy(&child);
	Cmiss_region_destroy(&root_region);
}

int main()
{
	int return_code = 0;
	Cmiss_context_id context = Cmiss_context_create("mesh_api");

	if (context && (CMISS_OK == Cmiss_context_enable_user_interface(context, NULL)))
	{
		Cmiss_region *root_region = Cmiss_context_get_default_region(context);
		create_cube_element(root_region);
		create_isosurface(context);
		export_stl(context, "cube.stl");
		netgen_mesh_generation(root_region, "cube.stl");
		show_final_mesh(context);
		Cmiss_context_execute_command(context, "gfx write elements nodes group final_mesh final_mesh.exelem");
		/* Uncomment the following 2line to visualise the final mesh*/
		//show_final_mesh(context);
		// Cmiss_context_run_main_loop(context);
		Cmiss_region_destroy(&root_region);
		Cmiss_context_destroy(&context);
	}

	return return_code;
}


Files used by this example are:

Name             Modified     Size

meshing_api.cpp 20-Apr-2012 14k Makefile 20-Apr-2012 2.8k

Download the entire example:

Name                                  Modified     Size

examples_a_backup_meshing_api.tar.gz 12-Aug-2014 22k

Testing status by version:

StatusTestedReal time (s)
x86-linux
cmgui-gtk-debugFailureSun Mar 6 00:06:35 20160
last breakThu Aug 7 16:04:00 20140
last successFri Jan 11 00:13:00 20133

Testing status by file:


Html last generated: Sun Mar 6 05:50:59 2016

Input last modified: Fri Apr 20 16:00:14 2012


CMISS Help / Examples / a_backup / meshing_api