The source file miscellaneous_ex6.C with comments:
#include "libmesh/elem.h"
#include "libmesh/face_tri3.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/mesh_tetgen_interface.h"
#include "libmesh/mesh_triangle_holes.h"
#include "libmesh/mesh_triangle_interface.h"
#include "libmesh/node.h"
#include "libmesh/serial_mesh.h"
Bring in everything from the libMesh namespace
using namespace libMesh;
Major functions called by main
void triangulate_domain();
void tetrahedralize_domain();
Helper routine for tetrahedralize_domain(). Adds the points and elements
of a convex hull generated by TetGen to the input mesh
void add_cube_convex_hull_to_mesh(MeshBase& mesh, Point lower_limit, Point upper_limit);
Begin the main program.
int main (int argc, char** argv)
{
Initialize libMesh and any dependent libaries, like in example 2.
LibMeshInit init (argc, argv);
libmesh_example_assert(2 <= LIBMESH_DIM, "2D support");
std::cout << "Triangulating an L-shaped domain with holes" << std::endl;
1.) 2D triangulation of L-shaped domain with three holes of different shape
triangulate_domain();
libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
std::cout << "Tetrahedralizing a prismatic domain with a hole" << std::endl;
2.) 3D tetrahedralization of rectangular domain with hole.
tetrahedralize_domain();
return 0;
}
void triangulate_domain()
{
#ifdef LIBMESH_HAVE_TRIANGLE
Use typedefs for slightly less typing.
typedef TriangleInterface::Hole Hole;
typedef TriangleInterface::PolygonHole PolygonHole;
typedef TriangleInterface::ArbitraryHole ArbitraryHole;
Libmesh mesh that will eventually be created.
Mesh mesh(2);
The points which make up the L-shape:
mesh.add_point(Point( 0. , 0.));
mesh.add_point(Point( 0. , -1.));
mesh.add_point(Point(-1. , -1.));
mesh.add_point(Point(-1. , 1.));
mesh.add_point(Point( 1. , 1.));
mesh.add_point(Point( 1. , 0.));
Declare the TriangleInterface object. This is where
we can set parameters of the triangulation and where the
actual triangulate function lives.
TriangleInterface t(mesh);
Customize the variables for the triangulation
t.desired_area() = .01;
A Planar Straight Line Graph (PSLG) is essentially a list
of segments which have to exist in the final triangulation.
For an L-shaped domain, Triangle will compute the convex
hull of boundary points if we do not specify the PSLG.
The PSLG algorithm is also required for triangulating domains
containing holes
t.triangulation_type() = TriangleInterface::PSLG;
Turn on/off Laplacian mesh smoothing after generation.
By default this is on.
t.smooth_after_generating() = true;
Define holes...
hole_1 is a circle (discretized by 50 points)
hole_1 is a circle (discretized by 50 points)
PolygonHole hole_1(Point(-0.5, 0.5), // center
0.25, // radius
50); // n. points
hole_2 is itself a triangle
PolygonHole hole_2(Point(0.5, 0.5), // center
0.1, // radius
3); // n. points
hole_3 is an ellipse of 100 points which we define here
Point ellipse_center(-0.5, -0.5);
const unsigned int n_ellipse_points=100;
std::vector<Point> ellipse_points(n_ellipse_points);
const Real
dtheta = 2*libMesh::pi / static_cast<Real>(n_ellipse_points),
a = .1,
b = .2;
for (unsigned int i=0; i<n_ellipse_points; ++i)
ellipse_points[i]= Point(ellipse_center(0)+a*cos(i*dtheta),
ellipse_center(1)+b*sin(i*dtheta));
ArbitraryHole hole_3(ellipse_center, ellipse_points);
Create the vector of Hole*'s ...
std::vector<Hole*> holes;
holes.push_back(&hole_1);
holes.push_back(&hole_2);
holes.push_back(&hole_3);
... and attach it to the triangulator object
t.attach_hole_list(&holes);
Triangulate!
t.triangulate();
Write the result to file
mesh.write("delaunay_l_shaped_hole.e");
#endif // LIBMESH_HAVE_TRIANGLE
}
void tetrahedralize_domain()
{
#ifdef LIBMESH_HAVE_TETGEN
The algorithm is broken up into several steps:
1.) A convex hull is constructed for a rectangular hole.
2.) A convex hull is constructed for the domain exterior.
3.) Neighbor information is updated so TetGen knows there is a convex hull
4.) A vector of hole points is created.
5.) The domain is tetrahedralized, the mesh is written out, etc.
The mesh we will eventually generate
The mesh we will eventually generate
SerialMesh mesh(3);
Lower and Upper bounding box limits for a rectangular hole within the unit cube.
Point hole_lower_limit(0.2, 0.2, 0.4);
Point hole_upper_limit(0.8, 0.8, 0.6);
1.) Construct a convex hull for the hole
add_cube_convex_hull_to_mesh(mesh, hole_lower_limit, hole_upper_limit);
2.) Generate elements comprising the outer boundary of the domain.
add_cube_convex_hull_to_mesh(mesh, Point(0.,0.,0.), Point(1., 1., 1.));
3.) Update neighbor information so that TetGen can verify there is a convex hull.
mesh.find_neighbors();
4.) Set up vector of hole points
std::vector<Point> hole(1);
hole[0] = Point( 0.5*(hole_lower_limit + hole_upper_limit) );
5.) Set parameters and tetrahedralize the domain
0 means "use TetGen default value"
0 means "use TetGen default value"
Real quality_constraint = 2.0;
The volume constraint determines the max-allowed tetrahedral
volume in the Mesh. TetGen will split cells which are larger than
this size
Real volume_constraint = 0.001;
Construct the Delaunay tetrahedralization
TetGenMeshInterface t(mesh);
t.triangulate_conformingDelaunayMesh_carvehole(hole,
quality_constraint,
volume_constraint);
Find neighbors, etc in preparation for writing out the Mesh
mesh.prepare_for_use();
Finally, write out the result
mesh.write("hole_3D.e");
#endif // LIBMESH_HAVE_TETGEN
}
void add_cube_convex_hull_to_mesh(MeshBase& mesh, Point lower_limit, Point upper_limit)
{
#ifdef LIBMESH_HAVE_TETGEN
SerialMesh cube_mesh(3);
unsigned n_elem = 1;
MeshTools::Generation::build_cube(cube_mesh,
n_elem,n_elem,n_elem, // n. elements in each direction
lower_limit(0), upper_limit(0),
lower_limit(1), upper_limit(1),
lower_limit(2), upper_limit(2),
HEX8);
The pointset_convexhull() algorithm will ignore the Hex8s
in the Mesh, and just construct the triangulation
of the convex hull.
TetGenMeshInterface t(cube_mesh);
t.pointset_convexhull();
Now add all nodes from the boundary of the cube_mesh to the input mesh.
Map from "node id in cube_mesh" -> "node id in mesh". Initially inserted with a dummy value, later to be assigned a value by the input mesh.
Map from "node id in cube_mesh" -> "node id in mesh". Initially inserted with a dummy value, later to be assigned a value by the input mesh.
std::map<unsigned,unsigned> node_id_map;
typedef std::map<unsigned,unsigned>::iterator iterator;
{
MeshBase::element_iterator it = cube_mesh.elements_begin();
const MeshBase::element_iterator end = cube_mesh.elements_end();
for ( ; it != end; ++it)
{
Elem* elem = *it;
for (unsigned s=0; s<elem->n_sides(); ++s)
if (elem->neighbor(s) == NULL)
{
Add the node IDs of this side to the set
AutoPtr<Elem> side = elem->side(s);
for (unsigned n=0; n<side->n_nodes(); ++n)
node_id_map.insert( std::make_pair(side->node(n), /*dummy_value=*/0) );
}
}
}
For each node in the map, insert it into the input mesh and keep
track of the ID assigned.
for (iterator it=node_id_map.begin(); it != node_id_map.end(); ++it)
{
Id of the node in the cube mesh
unsigned id = (*it).first;
Pointer to node in the cube mesh
Node* old_node = cube_mesh.node_ptr(id);
Add geometric point to input mesh
Node* new_node = mesh.add_point ( *old_node );
Track ID value of new_node in map
(*it).second = new_node->id();
}
With the points added and the map data structure in place, we are
ready to add each TRI3 element of the cube_mesh to the input Mesh
with proper node assignments
{
MeshBase::element_iterator el = cube_mesh.elements_begin();
const MeshBase::element_iterator end_el = cube_mesh.elements_end();
for (; el != end_el; ++el)
{
Elem* old_elem = *el;
if (old_elem->type() == TRI3)
{
Elem* new_elem = mesh.add_elem(new Tri3);
Assign nodes in new elements. Since this is an example,
we'll do it in several steps.
for (unsigned i=0; i<old_elem->n_nodes(); ++i)
{
Locate old node ID in the map
iterator it = node_id_map.find(old_elem->node(i));
Check for not found
if (it == node_id_map.end())
{
libMesh::err << "Node id " << old_elem->node(i) << " not found in map!" << std::endl;
libmesh_error();
}
Mapping to node ID in input mesh
unsigned new_node_id = (*it).second;
Node pointer assigned from input mesh
new_elem->set_node(i) = mesh.node_ptr(new_node_id);
}
}
}
}
#endif // LIBMESH_HAVE_TETGEN
}
The source file miscellaneous_ex6.C without comments:
#include "libmesh/elem.h"
#include "libmesh/face_tri3.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/mesh_tetgen_interface.h"
#include "libmesh/mesh_triangle_holes.h"
#include "libmesh/mesh_triangle_interface.h"
#include "libmesh/node.h"
#include "libmesh/serial_mesh.h"
using namespace libMesh;
void triangulate_domain();
void tetrahedralize_domain();
void add_cube_convex_hull_to_mesh(MeshBase& mesh, Point lower_limit, Point upper_limit);
int main (int argc, char** argv)
{
LibMeshInit init (argc, argv);
libmesh_example_assert(2 <= LIBMESH_DIM, "2D support");
std::cout << "Triangulating an L-shaped domain with holes" << std::endl;
triangulate_domain();
libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
std::cout << "Tetrahedralizing a prismatic domain with a hole" << std::endl;
tetrahedralize_domain();
return 0;
}
void triangulate_domain()
{
#ifdef LIBMESH_HAVE_TRIANGLE
typedef TriangleInterface::Hole Hole;
typedef TriangleInterface::PolygonHole PolygonHole;
typedef TriangleInterface::ArbitraryHole ArbitraryHole;
Mesh mesh(2);
mesh.add_point(Point( 0. , 0.));
mesh.add_point(Point( 0. , -1.));
mesh.add_point(Point(-1. , -1.));
mesh.add_point(Point(-1. , 1.));
mesh.add_point(Point( 1. , 1.));
mesh.add_point(Point( 1. , 0.));
TriangleInterface t(mesh);
t.desired_area() = .01;
t.triangulation_type() = TriangleInterface::PSLG;
t.smooth_after_generating() = true;
PolygonHole hole_1(Point(-0.5, 0.5), // center
0.25, // radius
50); // n. points
PolygonHole hole_2(Point(0.5, 0.5), // center
0.1, // radius
3); // n. points
Point ellipse_center(-0.5, -0.5);
const unsigned int n_ellipse_points=100;
std::vector<Point> ellipse_points(n_ellipse_points);
const Real
dtheta = 2*libMesh::pi / static_cast<Real>(n_ellipse_points),
a = .1,
b = .2;
for (unsigned int i=0; i<n_ellipse_points; ++i)
ellipse_points[i]= Point(ellipse_center(0)+a*cos(i*dtheta),
ellipse_center(1)+b*sin(i*dtheta));
ArbitraryHole hole_3(ellipse_center, ellipse_points);
std::vector<Hole*> holes;
holes.push_back(&hole_1);
holes.push_back(&hole_2);
holes.push_back(&hole_3);
t.attach_hole_list(&holes);
t.triangulate();
mesh.write("delaunay_l_shaped_hole.e");
#endif // LIBMESH_HAVE_TRIANGLE
}
void tetrahedralize_domain()
{
#ifdef LIBMESH_HAVE_TETGEN
SerialMesh mesh(3);
Point hole_lower_limit(0.2, 0.2, 0.4);
Point hole_upper_limit(0.8, 0.8, 0.6);
add_cube_convex_hull_to_mesh(mesh, hole_lower_limit, hole_upper_limit);
add_cube_convex_hull_to_mesh(mesh, Point(0.,0.,0.), Point(1., 1., 1.));
mesh.find_neighbors();
std::vector<Point> hole(1);
hole[0] = Point( 0.5*(hole_lower_limit + hole_upper_limit) );
Real quality_constraint = 2.0;
Real volume_constraint = 0.001;
TetGenMeshInterface t(mesh);
t.triangulate_conformingDelaunayMesh_carvehole(hole,
quality_constraint,
volume_constraint);
mesh.prepare_for_use();
mesh.write("hole_3D.e");
#endif // LIBMESH_HAVE_TETGEN
}
void add_cube_convex_hull_to_mesh(MeshBase& mesh, Point lower_limit, Point upper_limit)
{
#ifdef LIBMESH_HAVE_TETGEN
SerialMesh cube_mesh(3);
unsigned n_elem = 1;
MeshTools::Generation::build_cube(cube_mesh,
n_elem,n_elem,n_elem, // n. elements in each direction
lower_limit(0), upper_limit(0),
lower_limit(1), upper_limit(1),
lower_limit(2), upper_limit(2),
HEX8);
TetGenMeshInterface t(cube_mesh);
t.pointset_convexhull();
std::map<unsigned,unsigned> node_id_map;
typedef std::map<unsigned,unsigned>::iterator iterator;
{
MeshBase::element_iterator it = cube_mesh.elements_begin();
const MeshBase::element_iterator end = cube_mesh.elements_end();
for ( ; it != end; ++it)
{
Elem* elem = *it;
for (unsigned s=0; s<elem->n_sides(); ++s)
if (elem->neighbor(s) == NULL)
{
AutoPtr<Elem> side = elem->side(s);
for (unsigned n=0; n<side->n_nodes(); ++n)
node_id_map.insert( std::make_pair(side->node(n), /*dummy_value=*/0) );
}
}
}
for (iterator it=node_id_map.begin(); it != node_id_map.end(); ++it)
{
unsigned id = (*it).first;
Node* old_node = cube_mesh.node_ptr(id);
Node* new_node = mesh.add_point ( *old_node );
(*it).second = new_node->id();
}
{
MeshBase::element_iterator el = cube_mesh.elements_begin();
const MeshBase::element_iterator end_el = cube_mesh.elements_end();
for (; el != end_el; ++el)
{
Elem* old_elem = *el;
if (old_elem->type() == TRI3)
{
Elem* new_elem = mesh.add_elem(new Tri3);
for (unsigned i=0; i<old_elem->n_nodes(); ++i)
{
iterator it = node_id_map.find(old_elem->node(i));
if (it == node_id_map.end())
{
libMesh::err << "Node id " << old_elem->node(i) << " not found in map!" << std::endl;
libmesh_error();
}
unsigned new_node_id = (*it).second;
new_elem->set_node(i) = mesh.node_ptr(new_node_id);
}
}
}
}
#endif // LIBMESH_HAVE_TETGEN
}
The console output of the program:
***************************************************************
* Running Example miscellaneous_ex6:
* mpirun -np 2 example-devel -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
Triangulating an L-shaped domain with holes
Tetrahedralizing a prismatic domain with a hole
************************************************************************************************************************
*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***
************************************************************************************************************************
---------------------------------------------- PETSc Performance Summary: ----------------------------------------------
/workspace/libmesh/examples/miscellaneous/miscellaneous_ex6/.libs/lt-example-devel on a intel-12. named hbar.ices.utexas.edu with 2 processors, by benkirk Tue Feb 5 13:41:23 2013
Using Petsc Release Version 3.3.0, Patch 2, Fri Jul 13 15:42:00 CDT 2012
Max Max/Min Avg Total
Time (sec): 1.481e-01 1.00001 1.481e-01
Objects: 1.000e+00 1.00000 1.000e+00
Flops: 0.000e+00 0.00000 0.000e+00 0.000e+00
Flops/sec: 0.000e+00 0.00000 0.000e+00 0.000e+00
MPI Messages: 0.000e+00 0.00000 0.000e+00 0.000e+00
MPI Message Lengths: 0.000e+00 0.00000 0.000e+00 0.000e+00
MPI Reductions: 1.000e+00 1.00000
Flop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)
e.g., VecAXPY() for real vectors of length N --> 2N flops
and VecAXPY() for complex vectors of length N --> 8N flops
Summary of Stages: ----- Time ------ ----- Flops ----- --- Messages --- -- Message Lengths -- -- Reductions --
Avg %Total Avg %Total counts %Total Avg %Total counts %Total
0: Main Stage: 1.4803e-01 99.9% 0.0000e+00 0.0% 0.000e+00 0.0% 0.000e+00 0.0% 0.000e+00 0.0%
------------------------------------------------------------------------------------------------------------------------
See the 'Profiling' chapter of the users' manual for details on interpreting output.
Phase summary info:
Count: number of times phase was executed
Time and Flops: Max - maximum over all processors
Ratio - ratio of maximum to minimum over all processors
Mess: number of messages sent
Avg. len: average message length
Reduct: number of global reductions
Global: entire computation
Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().
%T - percent time in this phase %f - percent flops in this phase
%M - percent messages in this phase %L - percent message lengths in this phase
%R - percent reductions in this phase
Total Mflop/s: 10e-6 * (sum of flops over all processors)/(max time over all processors)
------------------------------------------------------------------------------------------------------------------------
Event Count Time (sec) Flops --- Global --- --- Stage --- Total
Max Ratio Max Ratio Max Ratio Mess Avg len Reduct %T %f %M %L %R %T %f %M %L %R Mflop/s
------------------------------------------------------------------------------------------------------------------------
--- Event Stage 0: Main Stage
------------------------------------------------------------------------------------------------------------------------
Memory usage is given in bytes:
Object Type Creations Destructions Memory Descendants' Mem.
Reports information only for process 0.
--- Event Stage 0: Main Stage
Viewer 1 0 0 0
========================================================================================================================
Average time to get PetscTime(): 9.53674e-08
Average time for MPI_Barrier(): 8.10623e-07
Average time for zero size MPI_Send(): 2.0504e-05
#PETSc Option Table entries:
-ksp_right_pc
-log_summary
-pc_type bjacobi
-sub_pc_factor_levels 4
-sub_pc_factor_zeropivot 0
-sub_pc_type ilu
#End of PETSc Option Table entries
Compiled without FORTRAN kernels
Compiled with full precision matrices (default)
sizeof(short) 2 sizeof(int) 4 sizeof(long) 8 sizeof(void*) 8 sizeof(PetscScalar) 8 sizeof(PetscInt) 4
Configure run at: Thu Nov 8 11:21:02 2012
Configure options: --with-debugging=false --COPTFLAGS=-O3 --CXXOPTFLAGS=-O3 --FOPTFLAGS=-O3 --with-clanguage=C++ --with-shared-libraries=1 --with-mpi-dir=/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1 --with-mumps=true --download-mumps=1 --with-metis=true --download-metis=1 --with-parmetis=true --download-parmetis=1 --with-superlu=true --download-superlu=1 --with-superludir=true --download-superlu_dist=1 --with-blacs=true --download-blacs=1 --with-scalapack=true --download-scalapack=1 --with-hypre=true --download-hypre=1 --with-blas-lib="[/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_intel_lp64.so,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_sequential.so,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_core.so]" --with-lapack-lib="[/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64/libmkl_lapack95_lp64.a]"
-----------------------------------------
Libraries compiled on Thu Nov 8 11:21:02 2012 on daedalus.ices.utexas.edu
Machine characteristics: Linux-2.6.32-279.1.1.el6.x86_64-x86_64-with-redhat-6.3-Carbon
Using PETSc directory: /opt/apps/ossw/libraries/petsc/petsc-3.3-p2
Using PETSc arch: intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt
-----------------------------------------
Using C compiler: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpicxx -wd1572 -O3 -fPIC ${COPTFLAGS} ${CFLAGS}
Using Fortran compiler: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpif90 -fPIC -O3 ${FOPTFLAGS} ${FFLAGS}
-----------------------------------------
Using include paths: -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/include -I/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/include -I/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/include
-----------------------------------------
Using C linker: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpicxx
Using Fortran linker: /opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/bin/mpif90
Using libraries: -Wl,-rpath,/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -L/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -lpetsc -lX11 -Wl,-rpath,/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -L/opt/apps/ossw/libraries/petsc/petsc-3.3-p2/intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt/lib -lcmumps -ldmumps -lsmumps -lzmumps -lmumps_common -lpord -lHYPRE -lpthread -lsuperlu_dist_3.0 -lparmetis -lmetis -lscalapack -lblacs -lsuperlu_4.3 -Wl,-rpath,/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64 -L/opt/apps/sysnet/intel/12.1/mkl/10.3.12.361/lib/intel64 -lmkl_lapack95_lp64 -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -Wl,-rpath,/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/lib -L/opt/apps/ossw/libraries/mpich2/mpich2-1.4.1p1/sl6/intel-12.1/lib -Wl,-rpath,/opt/apps/sysnet/intel/12.1/composer_xe_2011_sp1.7.256/compiler/lib/intel64 -L/opt/apps/sysnet/intel/12.1/composer_xe_2011_sp1.7.256/compiler/lib/intel64 -Wl,-rpath,/usr/lib/gcc/x86_64-redhat-linux/4.4.6 -L/usr/lib/gcc/x86_64-redhat-linux/4.4.6 -lmpichf90 -lifport -lifcore -lm -lm -lmpichcxx -ldl -lmpich -lopa -lmpl -lrt -lpthread -limf -lsvml -lipgo -ldecimal -lcilkrts -lstdc++ -lgcc_s -lirc -lirc_s -ldl
-----------------------------------------
----------------------------------------------------------------------------------------------------------------------
| Processor id: 0 |
| Num Processors: 2 |
| Time: Tue Feb 5 13:41:23 2013 |
| OS: Linux |
| HostName: hbar.ices.utexas.edu |
| OS Release: 2.6.32-279.1.1.el6.x86_64 |
| OS Version: #1 SMP Tue Jul 10 11:24:23 CDT 2012 |
| Machine: x86_64 |
| Username: benkirk |
| Configuration: ./configure '--enable-everything' |
| '--prefix=/workspace/libmesh/install' |
| 'CXX=mpicxx' |
| 'CC=mpicc' |
| 'F77=mpif77' |
| 'FC=mpif90' |
| 'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2' |
| 'PETSC_ARCH=intel-12.1-mkl-intel-10.3.12.361-mpich2-1.4.1p1-cxx-opt' |
| 'SLEPC_DIR=/opt/apps/ossw/libraries/slepc/slepc-3.3-p2-petsc-3.3-p2-cxx-opt' |
| 'TRILINOS_DIR=/opt/apps/ossw/libraries/trilinos/trilinos-10.12.2/sl6/intel-12.1/mpich2-1.4.1p1/mkl-intel-10.3.12.361'|
| 'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/intel-12.1' |
----------------------------------------------------------------------------------------------------------------------
-----------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=0.160684, Active time=0.090826 |
-----------------------------------------------------------------------------------------------------------
| Event nCalls Total Time Avg Time Total Time Avg Time % of Active Time |
| w/o Sub w/o Sub With Sub With Sub w/o S With S |
|-----------------------------------------------------------------------------------------------------------|
| |
| |
| Mesh |
| find_neighbors() 5 0.0260 0.005197 0.0264 0.005280 28.61 29.07 |
| renumber_nodes_and_elem() 6 0.0009 0.000142 0.0009 0.000142 0.94 0.94 |
| write() 2 0.0084 0.004214 0.0086 0.004276 9.28 9.42 |
| |
| MeshCommunication |
| compute_hilbert_indices() 2 0.0171 0.008563 0.0171 0.008563 18.86 18.86 |
| find_global_indices() 2 0.0055 0.002770 0.0251 0.012529 6.10 27.59 |
| parallel_sort() 2 0.0016 0.000787 0.0017 0.000828 1.73 1.82 |
| |
| MeshTools::Generation |
| build_cube() 2 0.0002 0.000087 0.0002 0.000087 0.19 0.19 |
| |
| MetisPartitioner |
| partition() 1 0.0245 0.024492 0.0367 0.036732 26.97 40.44 |
| |
| Parallel |
| allgather() 6 0.0003 0.000047 0.0003 0.000057 0.31 0.38 |
| broadcast() 2 0.0000 0.000006 0.0000 0.000006 0.01 0.01 |
| max(scalar) 129 0.0004 0.000003 0.0004 0.000003 0.49 0.49 |
| max(vector) 31 0.0002 0.000008 0.0005 0.000017 0.26 0.58 |
| min(bool) 153 0.0004 0.000003 0.0004 0.000003 0.49 0.49 |
| min(scalar) 124 0.0013 0.000011 0.0013 0.000011 1.49 1.49 |
| min(vector) 31 0.0003 0.000009 0.0007 0.000022 0.32 0.74 |
| probe() 10 0.0006 0.000061 0.0006 0.000061 0.67 0.67 |
| receive() 10 0.0001 0.000009 0.0007 0.000071 0.10 0.78 |
| send() 10 0.0000 0.000005 0.0000 0.000005 0.05 0.05 |
| send_receive() 14 0.0001 0.000007 0.0009 0.000065 0.12 1.00 |
| sum() 6 0.0000 0.000006 0.0001 0.000014 0.04 0.09 |
| |
| Parallel::Request |
| wait() 10 0.0000 0.000004 0.0000 0.000004 0.05 0.05 |
| |
| Partitioner |
| set_node_processor_ids() 1 0.0013 0.001311 0.0015 0.001485 1.44 1.63 |
| set_parent_processor_ids() 1 0.0013 0.001325 0.0013 0.001325 1.46 1.46 |
| single_partition() 2 0.0000 0.000017 0.0000 0.000017 0.04 0.04 |
-----------------------------------------------------------------------------------------------------------
| Totals: 562 0.0908 100.00 |
-----------------------------------------------------------------------------------------------------------
***************************************************************
* Done Running Example miscellaneous_ex6:
* mpirun -np 2 example-devel -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
Miscellaneous Example 6 - Meshing with LibMesh's TetGen and Triangle Interfaces
LibMesh provides interfaces to both Triangle and TetGen for generating Delaunay triangulations and tetrahedralizations in two and three dimensions (respectively).
Local header files