![]() |
FrontISTR
5.7.0
Large-scale structural analysis program with finit element method
|
This module encapsulate the basic functions of all elements provide by this software. More...
Functions/Subroutines | |
integer(kind=kind(2)) function | getspacedimension (etype) |
Obtain the space dimension of the element. More... | |
integer(kind=kind(2)) function | getnumberofnodes (etype) |
Obtain number of nodes of the element. More... | |
integer(kind=kind(2)) function | getnumberofsubface (etype) |
Obtain number of sub-surface. More... | |
subroutine | getsubface (intype, innumber, outtype, nodes) |
Find the definition of surface of the element. More... | |
integer function | numofquadpoints (fetype) |
Obtains the number of quadrature points of the element. More... | |
subroutine | getquadpoint (fetype, np, pos) |
Fetch the coordinate of gauss point. More... | |
real(kind=kreal) function | getweight (fetype, np) |
Fetch the weight value in given gauss point. More... | |
subroutine | getshapederiv (fetype, localcoord, shapederiv) |
Calculate derivatives of shape function in natural coordinate system. More... | |
subroutine | getshape2ndderiv (fetype, localcoord, shapederiv) |
Calculate the 2nd derivative of shape function in natural coordinate system. More... | |
subroutine | getshapefunc (fetype, localcoord, func) |
Calculate the shape function in natural coordinate system. More... | |
subroutine | getnodalnaturalcoord (fetype, nncoord) |
subroutine | getglobalderiv (fetype, nn, localcoord, elecoord, det, gderiv) |
Calculate shape derivative in global coordinate system. More... | |
real(kind=kreal) function | getdeterminant (fetype, nn, localcoord, elecoord) |
Calculate shape derivative in global coordinate system. More... | |
subroutine | getjacobian (fetype, nn, localcoord, elecoord, det, jacobian, inverse) |
calculate Jacobian matrix, its determinant and inverse More... | |
real(kind=kreal) function, dimension(3) | surfacenormal (fetype, nn, localcoord, elecoord) |
Calculate normal of 3d-surface. More... | |
real(kind=kreal) function, dimension(2) | edgenormal (fetype, nn, localcoord, elecoord) |
Calculate normal of 2d-edge. More... | |
subroutine | tangentbase (fetype, nn, localcoord, elecoord, tangent) |
Calculate base vector of tangent space of 3d surface. More... | |
subroutine | curvature (fetype, nn, localcoord, elecoord, l2ndderiv, normal, curv) |
Calculate curvature tensor at a point along 3d surface. More... | |
subroutine | getelementcenter (fetype, localcoord) |
Return natural coordinate of the center of surface element. More... | |
integer function | isinsideelement (fetype, localcoord, clearance) |
if a point is inside a surface element -1: No; 0: Yes; >0: Node's (vertex) number More... | |
integer function | isinside3delement (fetype, localcoord, clearance) |
if a point is inside a surface element -1: No; 0: Yes; >0: Node's (vertex) number More... | |
subroutine | getvertexcoord (fetype, cnode, localcoord) |
Get the natural coord of a vertex node. More... | |
subroutine | extrapolatevalue (lpos, fetype, nnode, pvalue, ndvalue) |
This subroutine extrapolate a point value into elemental nodes. More... | |
subroutine | interapolatevalue (lpos, fetype, nnode, pvalue, ndvalue) |
This subroutine interapolate element nodes value into a point value. More... | |
subroutine | gauss2node (fetype, gaussv, nodev) |
This subroutine extroplate value in quadrature point to element nodes. More... | |
real(kind=kreal) function | getreferencelength (fetype, nn, localcoord, elecoord) |
This function calculates reference length at a point in surface. More... | |
Variables | |
integer, parameter | fe_unknown = -1 |
integer, parameter | fe_line2n = 111 |
integer, parameter | fe_line3n = 112 |
integer, parameter | fe_tri3n = 231 |
integer, parameter | fe_tri6n = 232 |
integer, parameter | fe_tri6nc = 2322 |
integer, parameter | fe_quad4n = 241 |
integer, parameter | fe_quad8n = 242 |
integer, parameter | fe_truss = 301 |
integer, parameter | fe_tet4n = 341 |
integer, parameter | fe_tet4n_pipi = 3414 |
integer, parameter | fe_tet10n = 342 |
integer, parameter | fe_tet10nc = 3422 |
integer, parameter | fe_prism6n = 351 |
integer, parameter | fe_prism15n = 352 |
integer, parameter | fe_hex8n = 361 |
integer, parameter | fe_hex20n = 362 |
integer, parameter | fe_hex27n = 363 |
integer, parameter | fe_if_line2n = 511 |
integer, parameter | fe_beam2n = 611 |
integer, parameter | fe_beam3n = 612 |
integer, parameter | fe_beam341 = 641 |
integer, parameter | fe_tri6n_shell = 732 |
integer, parameter | fe_dsg3_shell = 733 |
integer, parameter | fe_mitc3_shell = 731 |
integer, parameter | fe_mitc4_shell = 741 |
integer, parameter | fe_mitc8_shell = 742 |
integer, parameter | fe_mitc9_shell = 743 |
integer, parameter | fe_mitc3_shell361 = 761 |
integer, parameter | fe_mitc4_shell361 = 781 |
integer, parameter | fe_nodesmooth_tet4n = 881 |
integer, parameter | fe_edgesmooth_tet4n = 891 |
integer, parameter | fe_tri3n_patch = 1031 |
integer, parameter | fe_tri6n_patch = 1032 |
integer, parameter | fe_quad4n_patch = 1041 |
integer, parameter | fe_quad8n_patch = 1042 |
This module encapsulate the basic functions of all elements provide by this software.
After providing a specified elemental type, this module provide functions to fetch
Elemental surfaces (Nodes' number)
Those would be supposed to adopt this module include
And something else ...
subroutine elementinfo::curvature | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(2), intent(in) | localcoord, | ||
real(kind=kreal), dimension(3,nn), intent(in) | elecoord, | ||
real(kind=kreal), dimension(3,2,2), intent(out) | l2ndderiv, | ||
real(kind=kreal), dimension(3), intent(in), optional | normal, | ||
real(kind=kreal), dimension(2,2), intent(out), optional | curv | ||
) |
Calculate curvature tensor at a point along 3d surface.
[in] | fetype | type of surface element |
[in] | nn | number of elemental nodes |
[in] | localcoord | position |
[in] | elecoord | nodes coordinates of element |
[out] | l2ndderiv | 2nd derivative of shape function |
[in] | normal | normal direction of surface |
[out] | curv | curvature tensor |
Definition at line 967 of file element.f90.
real(kind=kreal) function, dimension(2) elementinfo::edgenormal | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(1), intent(in) | localcoord, | ||
real(kind=kreal), dimension(2,nn), intent(in) | elecoord | ||
) |
Calculate normal of 2d-edge.
[in] | fetype | type of surface element |
[in] | nn | number of elemental nodes |
[in] | localcoord | position |
[in] | elecoord | nodes coordinates of element |
Definition at line 905 of file element.f90.
subroutine elementinfo::extrapolatevalue | ( | real(kind=kreal), dimension(:), intent(in) | lpos, |
integer, intent(in) | fetype, | ||
integer, intent(in) | nnode, | ||
real(kind=kreal), dimension(:), intent(in) | pvalue, | ||
real(kind=kreal), dimension(:,:), intent(out) | ndvalue | ||
) |
This subroutine extrapolate a point value into elemental nodes.
[in] | lpos | position of value given |
[in] | fetype | element type |
[in] | nnode | number of element node |
[in] | pvalue | value to be extropolated |
[out] | ndvalue | equivalent nodal value |
Definition at line 1171 of file element.f90.
subroutine elementinfo::gauss2node | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(:,:), intent(in) | gaussv, | ||
real(kind=kreal), dimension(:,:), intent(out) | nodev | ||
) |
This subroutine extroplate value in quadrature point to element nodes.
[in] | fetype | element type |
[in] | gaussv | values in quadrature points |
[out] | nodev | values in nodes |
Definition at line 1204 of file element.f90.
real(kind=kreal) function elementinfo::getdeterminant | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:,:), intent(in) | elecoord | ||
) |
Calculate shape derivative in global coordinate system.
[in] | fetype | element type |
[in] | nn | number of elemental nodes |
[in] | localcoord | curr position with natural coord |
[in] | elecoord | nodal coord of curr element |
Definition at line 792 of file element.f90.
subroutine elementinfo::getelementcenter | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(2), intent(out) | localcoord | ||
) |
Return natural coordinate of the center of surface element.
[in] | fetype | type of surface element |
[out] | localcoord | center coordinate |
Definition at line 1012 of file element.f90.
subroutine elementinfo::getglobalderiv | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:,:), intent(in) | elecoord, | ||
real(kind=kreal), intent(out) | det, | ||
real(kind=kreal), dimension(:,:), intent(out) | gderiv | ||
) |
Calculate shape derivative in global coordinate system.
[in] | fetype | element type |
[in] | nn | number of elemental nodes |
[in] | localcoord | curr position with natural coord |
[in] | elecoord | nodal coord of curr element |
[out] | det | nodal coord of curr element |
[out] | gderiv | shape derivative in global coordinate system |
Definition at line 741 of file element.f90.
subroutine elementinfo::getjacobian | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:,:), intent(in) | elecoord, | ||
real(kind=kreal), intent(out) | det, | ||
real(kind=kreal), dimension(:,:), intent(out) | jacobian, | ||
real(kind=kreal), dimension(:,:), intent(out) | inverse | ||
) |
calculate Jacobian matrix, its determinant and inverse
[in] | fetype | element type |
[in] | nn | number of element nodes |
[in] | localcoord | curr position with natural coord |
[in] | elecoord | nodal coord of curr element |
[out] | det | nodal coord of curr element |
[out] | inverse | inverse of jacobian |
Definition at line 820 of file element.f90.
subroutine elementinfo::getnodalnaturalcoord | ( | integer, intent(in) | fetype, |
real(kind = kreal), dimension(:, :), intent(out) | nncoord | ||
) |
Definition at line 699 of file element.f90.
integer(kind=kind(2)) function elementinfo::getnumberofnodes | ( | integer, intent(in) | etype | ) |
Obtain number of nodes of the element.
[in] | etype | element type |
Definition at line 131 of file element.f90.
integer(kind=kind(2)) function elementinfo::getnumberofsubface | ( | integer, intent(in) | etype | ) |
Obtain number of sub-surface.
[in] | etype | element type |
Definition at line 168 of file element.f90.
subroutine elementinfo::getquadpoint | ( | integer, intent(in) | fetype, |
integer, intent(in) | np, | ||
real(kind=kreal), dimension(:), intent(out) | pos | ||
) |
Fetch the coordinate of gauss point.
[in] | fetype | element type |
[in] | np | number of curr quadrature point |
[out] | pos | natural coord of curr quadrature point |
Definition at line 489 of file element.f90.
real(kind=kreal) function elementinfo::getreferencelength | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(2), intent(in) | localcoord, | ||
real(kind=kreal), dimension(3,nn), intent(in) | elecoord | ||
) |
This function calculates reference length at a point in surface.
[in] | fetype | surface element type |
[in] | nn | number of elemental nodes |
[in] | localcoord | natural coordinates |
[in] | elecoord | nodes coordinates of surface element |
Definition at line 1259 of file element.f90.
subroutine elementinfo::getshape2ndderiv | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:,:,:), intent(out) | shapederiv | ||
) |
Calculate the 2nd derivative of shape function in natural coordinate system.
[in] | fetype | elemental type |
[in] | localcoord | natural points |
[out] | shapederiv | 2nd order shape derivatives |
Definition at line 622 of file element.f90.
subroutine elementinfo::getshapederiv | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:,:), intent(out) | shapederiv | ||
) |
Calculate derivatives of shape function in natural coordinate system.
[in] | fetype | input element type |
[in] | localcoord | natural points |
[out] | shapederiv | derivative of shape function |
Definition at line 578 of file element.f90.
subroutine elementinfo::getshapefunc | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(:), intent(in) | localcoord, | ||
real(kind=kreal), dimension(:), intent(out) | func | ||
) |
Calculate the shape function in natural coordinate system.
[in] | fetype | input element type |
[in] | localcoord | natural points |
[out] | func | shape function |
Definition at line 647 of file element.f90.
integer(kind=kind(2)) function elementinfo::getspacedimension | ( | integer, intent(in) | etype | ) |
Obtain the space dimension of the element.
[in] | etype | element type |
Definition at line 117 of file element.f90.
subroutine elementinfo::getsubface | ( | integer, intent(in) | intype, |
integer, intent(in) | innumber, | ||
integer, intent(out) | outtype, | ||
integer, dimension(:), intent(out) | nodes | ||
) |
Find the definition of surface of the element.
[in] | intype | input element type, with space dimension n |
[in] | innumber | number of sub-surface |
[out] | outtype | output element type, with space dimension n-1 |
[out] | nodes | output face nodes' ID |
Definition at line 193 of file element.f90.
subroutine elementinfo::getvertexcoord | ( | integer, intent(in) | fetype, |
integer, intent(in) | cnode, | ||
real(kind=kreal), dimension(2), intent(out) | localcoord | ||
) |
Get the natural coord of a vertex node.
[in] | fetype | type of surface element |
[in] | cnode | current node |
[out] | localcoord | natural coord |
Definition at line 1136 of file element.f90.
real(kind=kreal) function elementinfo::getweight | ( | integer, intent(in) | fetype, |
integer, intent(in) | np | ||
) |
Fetch the weight value in given gauss point.
[in] | fetype | element type |
[in] | np | number of curr quadrature point |
Definition at line 535 of file element.f90.
subroutine elementinfo::interapolatevalue | ( | real(kind=kreal), dimension(:), intent(in) | lpos, |
integer, intent(in) | fetype, | ||
integer, intent(in) | nnode, | ||
real(kind=kreal), dimension(:), intent(out) | pvalue, | ||
real(kind=kreal), dimension(:,:), intent(in) | ndvalue | ||
) |
This subroutine interapolate element nodes value into a point value.
[in] | lpos | position of value given |
[in] | fetype | element type |
[in] | nnode | number of element node |
[out] | pvalue | value to be extropolated |
[in] | ndvalue | equivalent nodal value |
Definition at line 1187 of file element.f90.
integer function elementinfo::isinside3delement | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(3), intent(inout) | localcoord, | ||
real(kind=kreal), optional | clearance | ||
) |
if a point is inside a surface element -1: No; 0: Yes; >0: Node's (vertex) number
[in] | fetype | type of surface element |
[in,out] | localcoord | natural coord |
clearance | clearance used for judgement |
Definition at line 1094 of file element.f90.
integer function elementinfo::isinsideelement | ( | integer, intent(in) | fetype, |
real(kind=kreal), dimension(2), intent(inout) | localcoord, | ||
real(kind=kreal), optional | clearance | ||
) |
if a point is inside a surface element -1: No; 0: Yes; >0: Node's (vertex) number
[in] | fetype | type of surface element |
[in,out] | localcoord | natural coord |
clearance | clearance used for judgement |
Definition at line 1029 of file element.f90.
integer function elementinfo::numofquadpoints | ( | integer, intent(in) | fetype | ) |
Obtains the number of quadrature points of the element.
[in] | fetype | element type |
Definition at line 450 of file element.f90.
real(kind=kreal) function, dimension(3) elementinfo::surfacenormal | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(2), intent(in) | localcoord, | ||
real(kind=kreal), dimension(3,nn), intent(in) | elecoord | ||
) |
Calculate normal of 3d-surface.
[in] | fetype | type of surface element |
[in] | nn | number of elemental nodes |
[in] | localcoord | position |
[in] | elecoord | nodes coordinates of element |
Definition at line 870 of file element.f90.
subroutine elementinfo::tangentbase | ( | integer, intent(in) | fetype, |
integer, intent(in) | nn, | ||
real(kind=kreal), dimension(2), intent(in) | localcoord, | ||
real(kind=kreal), dimension(3,nn), intent(in) | elecoord, | ||
real(kind=kreal), dimension(3,2), intent(out) | tangent | ||
) |
Calculate base vector of tangent space of 3d surface.
[in] | fetype | type of surface element |
[in] | nn | number of elemental nodes |
[in] | localcoord | position |
[in] | elecoord | nodes coordinates of element |
[out] | tangent | two tangent vectors |
Definition at line 933 of file element.f90.
integer, parameter elementinfo::fe_beam2n = 611 |
Definition at line 87 of file element.f90.
integer, parameter elementinfo::fe_beam341 = 641 |
Definition at line 89 of file element.f90.
integer, parameter elementinfo::fe_beam3n = 612 |
Definition at line 88 of file element.f90.
integer, parameter elementinfo::fe_dsg3_shell = 733 |
Definition at line 92 of file element.f90.
integer, parameter elementinfo::fe_edgesmooth_tet4n = 891 |
Definition at line 102 of file element.f90.
integer, parameter elementinfo::fe_hex20n = 362 |
Definition at line 82 of file element.f90.
integer, parameter elementinfo::fe_hex27n = 363 |
Definition at line 83 of file element.f90.
integer, parameter elementinfo::fe_hex8n = 361 |
Definition at line 81 of file element.f90.
integer, parameter elementinfo::fe_if_line2n = 511 |
Definition at line 85 of file element.f90.
integer, parameter elementinfo::fe_line2n = 111 |
Definition at line 67 of file element.f90.
integer, parameter elementinfo::fe_line3n = 112 |
Definition at line 68 of file element.f90.
integer, parameter elementinfo::fe_mitc3_shell = 731 |
Definition at line 93 of file element.f90.
integer, parameter elementinfo::fe_mitc3_shell361 = 761 |
Definition at line 98 of file element.f90.
integer, parameter elementinfo::fe_mitc4_shell = 741 |
Definition at line 94 of file element.f90.
integer, parameter elementinfo::fe_mitc4_shell361 = 781 |
Definition at line 99 of file element.f90.
integer, parameter elementinfo::fe_mitc8_shell = 742 |
Definition at line 95 of file element.f90.
integer, parameter elementinfo::fe_mitc9_shell = 743 |
Definition at line 96 of file element.f90.
integer, parameter elementinfo::fe_nodesmooth_tet4n = 881 |
Definition at line 101 of file element.f90.
integer, parameter elementinfo::fe_prism15n = 352 |
Definition at line 80 of file element.f90.
integer, parameter elementinfo::fe_prism6n = 351 |
Definition at line 79 of file element.f90.
integer, parameter elementinfo::fe_quad4n = 241 |
Definition at line 72 of file element.f90.
integer, parameter elementinfo::fe_quad4n_patch = 1041 |
Definition at line 106 of file element.f90.
integer, parameter elementinfo::fe_quad8n = 242 |
Definition at line 73 of file element.f90.
integer, parameter elementinfo::fe_quad8n_patch = 1042 |
Definition at line 107 of file element.f90.
integer, parameter elementinfo::fe_tet10n = 342 |
Definition at line 77 of file element.f90.
integer, parameter elementinfo::fe_tet10nc = 3422 |
Definition at line 78 of file element.f90.
integer, parameter elementinfo::fe_tet4n = 341 |
Definition at line 75 of file element.f90.
integer, parameter elementinfo::fe_tet4n_pipi = 3414 |
Definition at line 76 of file element.f90.
integer, parameter elementinfo::fe_tri3n = 231 |
Definition at line 69 of file element.f90.
integer, parameter elementinfo::fe_tri3n_patch = 1031 |
Definition at line 104 of file element.f90.
integer, parameter elementinfo::fe_tri6n = 232 |
Definition at line 70 of file element.f90.
integer, parameter elementinfo::fe_tri6n_patch = 1032 |
Definition at line 105 of file element.f90.
integer, parameter elementinfo::fe_tri6n_shell = 732 |
Definition at line 91 of file element.f90.
integer, parameter elementinfo::fe_tri6nc = 2322 |
Definition at line 71 of file element.f90.
integer, parameter elementinfo::fe_truss = 301 |
Definition at line 74 of file element.f90.
integer, parameter elementinfo::fe_unknown = -1 |
Definition at line 65 of file element.f90.