integrules

Contents Functions

integrules

Purpose

Command function for FEM integration rule support.

Description

This function groups integration rule manipulation utilities used by various elements. In terms of notations, a field u is interpolated within an element by shapes functions N_i and values of the field at nodes u_i

u(x,y,z)=

∑

N_i(r,s,t)u_i (9.16)

The relation between physical coordinates x,y,z and element coordinates r,s,t is itself described by a mapping associated with shape functions. When computing an integral, one selects a number of Gauss points r_g,s_g,t_g and associated weights w_g leading to an approximation of the integral as

(9.17)

where J is the determinant of the Jacobian of the transform from reference to physical coordinates. The field .wjdet is used to denote the local value of the product J w_g. The following calls generate the reference EltConst data structure, see section 7.15.4.

Gauss

This command supports the definition of Gauss points and associated weights. It is called with integrules('Gauss Topology',RuleNumber). Supported topologies are 1d (line), q2d (2D quadrangle), t2d (2D triangle), t3d (3D tetrahedron), p3d (3D prism), h3d (3D hexahedron). integrules('Gauss q2d') will list available 2D quadrangle rules.

Integ -3 is always the default rule for the order of the element.
-2 a rule at nodes.
-1 the rule at center.

    [ -3]    [ 0x1 double]    'element dependent default'        
    [ -2]    [ 0x1 double]    'node'        
    [ -1]    [ 1x4 double]    'center'      
    [102]    [ 4x4 double]    'gefdyn 2x2'  
    [  2]    [ 4x4 double]    'standard 2x2'
    [109]    [ 9x4 double]    'Q4WT'        
    [103]    [ 9x4 double]    'gefdyn 3x3'  
    [104]    [16x4 double]    'gefdyn 4x4'  
    [  9]    [ 9x4 double]    '9 point'     
    [  3]    [ 9x4 double]    'standard 3x3'
    [  2]    [ 4x4 double]    'standard 2x2'
    [ 13]    [13x4 double]    '2x2 and 3x3'

bar1,beam1,beam3

For integration rule selection, these elements use the 1D rules which list you can find using integrules('Gauss1d').

Geometric orientation convention for segment is • (1) → (2)

One can show the edge using elt_name edge (e.g. beaml edge).

t3p,t6p

Vertex coordinates of the reference element can be found using r1=integrules('tria3');r1.xi.

Figure 9.1: t3p reference element.

Vertex coordinates of the reference element can be found using r1=integrules('tria6');r1.xi.

Figure 9.2: t6p reference element.

For integration rule selection, these elements use the 2D triangle rules which list you can find using integrules('Gausst2d').

Geometric orientation convention for triangle is to number anti-clockwise in the two-dimensional case (in the three-dimensional case, there is no orientation).
• edge [1]: (1) → (2) (nodes 4, 5,... if there are supplementary nodes) • edge [2]: (2) → (3) (...) • edge [3]: (3) → (1) (...)

One can show the edges or faces using elt_name edge or elt_name face (e.g. t3p edge).

q4p, q5p, q8p

Vertex coordinates of the reference element can be found using r1=integrules('quad4');r1.xi.

Figure 9.3: q4p reference element.

Figure 9.4: q5p reference element.

Vertex coordinates of the reference element can be found using the r1=integrules('quadb');r1.xi.

Figure 9.5: q8p reference element.

For integration rule selection, these elements use the 2D quadrangle rules which list you can find using integrules('Gaussq2d').

Geometric orientation convention for quadrilateral is to number anti-clockwise (same remark as for the triangle)
• edge [1]: (1) → (2) (nodes 5, 6, ...) • edge [2]: (2) → (3) (...) • edge [3]: (3) → (4) • edge [4]: (4) → (1)

One can show the edges or faces using elt_name edge or elt_name face (e.g. q4p edge).

tetra4,tetra10

3D tetrahedron geometries with linear and quadratic shape functions. Vertex coordinates of the reference element can be found using r1=integrules('tetra4');r1.xi (or command 'tetra10').

Figure 9.6: tetra4 reference element.

Figure 9.7: tetra10 reference element.

For integration rule selection, these elements use the 3D pentahedron rules which list you can find using integrules('Gausst3d').

Geometric orientation convention for tetrahedron is to have trihedral (12^→,13^→,14^→) direct (ij^→ designates the vector from point i to point j).
• edge [1]: (1) → (2) (nodes 5, ...) • edge [2]: (2) → (3) (...) • edge [3]: (3) → (1)
• edge [4]: (1) → (4) • edge [5]: (2) → (4) • edge [6]: (3) → (4) (nodes ..., p)

All faces, seen from the exterior, are described anti-clockwise:
• face [1]: (1) (3) (2) (nodes p+1, ...) • face [2]: (1) (4) (3) (...)
• face [3]: (1) (2) (4) • face [4]: (2) (3) (4)

One can show the edges or faces using elt_name edge or elt_name face (e.g. tetra10 face).

penta6, penta15

3D prism geometries with linear and quadratic shape functions. Vertex coordinates of the reference element can be found using r1=integrules('penta6');r1.xi (or command 'penta15').

Figure 9.8: penta6 reference element.

Figure 9.9: penta15 reference element.

For integration rule selection, these elements use the 3D pentahedron rules which list you can find using integrules('Gaussp3d').

Geometric orientation convention for pentahedron is to have trihedral (12^→,13^→,14^→) direct
• edge [1]: (1) → (2) (nodes 7, ...) • edge [2]: (2) → (3) (...) • edge [3]: (3) → (1)
• edge [4]: (1) → (4) • edge [5]: (2) → (5) • edge [6]: (3) → (6)
• edge [7]: (4) → (5) • edge [8]: (5) → (6) • edge [9]: (6) → (4) (nodes ..., p)

All faces, seen from the exterior, are described anti-clockwise.
• face [1] : (1) (3) (2) (nodes p+1, ...) • face [2] : (1) (4) (6) (3) • face [3] : (1) (2) (5) (4)
• face [4] : (4) (5) (6) • face [5] : (2) (3) (5) (6)

One can show the edges or faces using elt_name edge or elt_name face (e.g. penta15 face).

hexa8, hexa20, hexa21, hexa27

3D brick geometries, using linear hexa8, and quadratic shape functions. Vertex coordinates of the reference element can be found using r1=integrules('hexa8');r1.xi (or command 'hexa20', 'hexa27').

Figure 9.10: hexa8 reference topology.

Figure 9.11: hexa20 reference topology.

For integration rule selection, these elements use the 3D hexahedron rules which list you can find using integrules('Gaussh3d').

Geometric orientation convention for hexahedron is to have trihedral (12^→,14^→,15^→) direct
• edge [1]: (1) → (2) (nodes 9, ...) • edge [2]: (2) → (3) (...) • edge [3]: (3) → (4)
• edge [4]: (4) → (1) • edge [5]: (1) → (5) • edge [6]: (2) → (6)
• edge [7]: (3) → (7) • edge [8]: (4) → (8) • edge [9]: (5) → (6)
• edge [10]: (6) → (7) • edge [11]: (7) → (8) • edge [12]: (8) → (5) (nodes ..., p)

All faces, seen from the exterior, are described anti-clockwise.
• face [1] : (1) (4) (3) (2) (nodes p+1, ...) • face [2] : (1) (5) (8) (4)
• face [3] : (1) (2) (6) (5) • face [4] : (5) (6) (7) (8)
• face [5] : (2) (3) (7) (6) • face [6] : (3) (4) (8) (7)

One can show the edges or faces using elt_name edge or elt_name face (e.g. hexa8 face).

BuildNDN

The commands are extremely low level utilities to fill the .NDN field for a given set of nodes. The calling format is of_mk('BuildNDN',type,rule,nodeE) where type is an int32 that specifies the rule to be used : 2 for 2D, 3 for 3D, 31 for 3D with xyz sorting of NDN columns, 23 for surface in a 3D model, 13 for a 3D line. A negative value can be used to switch to the .m file implementation in integrules.

The 23 rule generates a transformation with the first axis along N,r, the second axis orthogonal in the plane tangent to N,r, N,s and the third axis locally normal to the element surface. If a local material orientation is provided in columns 5 to 7 of nodeE then the material x axis is defined by projection on the surface. One recalls that columns of nodeE are field based on the InfoAtNode field and the first three labels should be 'v1x','v1y','v1z'.

With the 32 rule if a local material orientation is provided in columns 5 to 7 for x and 8 to 10 for y the spatial derivatives of the shape functions are given in this local frame.

The rule structure is described earlier in this section and node has three columns that give the positions in the nodes of the current element. The rule.NDN and rule.jdet fields are modified. They must have the correct size before the call is made or severe crashes can be experienced.

If a rule.bas field is defined (9× Nw), each column is filled to contain the local basis at the integration point for 23 and 13 types. If a rule.J field with (4× Nw), each column is filled to contain the jacobian at the integration point for 23.

model=femesh('testhexa8');  nodeE=model.Node(:,5:7);
opt=integrules('hexa8',-1);
nodeE(:,5:10)=0; nodeE(:,7)=1;  nodeE(:,8)=1; % xe=z and ye=y
integrules('buildndn',32,opt,nodeE)

model=femesh('testquad4'); nodeE=model.Node(:,5:7);
opt=integrules('q4p',-1);opt.bas=zeros(9,opt.Nw);opt.J=zeros(4,opt.Nw);
nodeE(:,5:10)=0; nodeE(:,5:6)=1;  % xe= along [1,1,0]
integrules('buildndn',23,opt,nodeE)