How to use set geometry from 10 node tetrahedrons ?

[2109...@gmail.com]

How to generate higher-order geometry and use it in MoFEM ?

[Lukasz Kaczmraczyk]

If you are reading mesh with 10-node tetrahedrons, you have nodes in the middle of each edge. Since you have four corner nodes and six edges that give 10.

MoFEM does not work directly on 10-node elements but applies some general strategy to represent geometry using arbitrary order and heterogeneous approximation base. That allows for greater flexibility. 

1) To represent geometry you can define field,

ierr = m_field.add_field("MESH_NODE_POSITIONS",H1, AINSWORTH_LEGENDRE_BASE ,3); CHKERRQ(ierr);

• 1st argument is name of field which you like to have
  
• 2nd is approximation space
  
• 3rd approximation base  
  
• 4th option indicate that field has three coefficients
  

See for details this link add_field.

2) Add tetrahedrons to that field

ierr = m_field.add_ents_to_field_by_TETs(0,"MESH_NODE_POSITIONS"); CHKERRQ(ierr);

• 1st argument tells to add all tetrahedrons (and lower dimension entities adjacent to tetrahedrons) in root mesh to the field,
• 2nd name of the field to which add tetrahedrons

3) Set approximation order, for simplicity 2nd order to all edges

ierr = m_field.set_field_order(0,MBEDGE,"MESH_NODE_POSITIONS",2); CHKERRQ(ierr);
ierr = m_field.set_field_order(0,MBVERTEX,"MESH_NODE_POSITIONS",1); CHKERRQ(ierr);

• 1st tells that all base functions from edges of the field have polynomial of order 2
• 2nd name is name of the filed to which you set approximation order
• 3rd is order of polynomial in that case second

The same command apply for vertices, but in that case base functions from degrees of freedom from vertices have always polynomial  1 at least for that space. Setting other value will generate error.

4) Build field

 ierr = m_field.build_fields(); CHKERRQ(ierr);

That command construct data structures and allocate memory on the mesh. 

5) The last step is to project information from 10 node tetrahedrons on hierarchical approximation field, as follows

Projection10NodeCoordsOnField ent_method_material(m_field,"MESH_NODE_POSITIONS");
ierr = m_field.loop_dofs("MESH_NODE_POSITIONS",ent_method_material); CHKERRQ(ierr);

You using implementation of entity operator to that, i.e.  Projection10NodeCoordsOnField.

Another question is how you can generate 10-node tetrahedrons for example in cubit, you do that by setting block properties, for example for volume 1

set duplicate block elements on
block 2 volume 1
block 2 element type TETRA10