Elmer FEM solver
Elmer is an open source finite element software for multiphysical problems
 All Classes Files Functions Variables Typedefs Macros Groups Pages
elementdescription Module Reference
integer, parameter max_element_nodes = 256
 
subroutine addelementdescription (element, BasisTerms)
 
subroutine initializeelementdescriptions ()
 
type(elementtype_t) function,
pointer 
getelementtype (code, CompStabFlag)
 
subroutine stabparam (Element, Nodes, n, mK, hK)
 
real(kind=dp) function interpolateinelement1d (element, x, u)
 
subroutine nodalbasisfunctions1d (y, element, u)
 
real(kind=dp) function firstderivative1d (element, x, u)
 
subroutine nodalfirstderivatives1d (y, element, u)
 
real(kind=dp) function secondderivatives1d (element, x, u)
 
real(kind=dp) function interpolateinelement2d (element, x, u, v)
 
subroutine nodalbasisfunctions2d (y, element, u, v)
 
real(kind=dp) function firstderivativeinu2d (element, x, u, v)
 
real(kind=dp) function firstderivativeinv2d (element, x, u, v)
 
subroutine nodalfirstderivatives2d (y, element, u, v)
 
real(kind=dp) function,
dimension(2, 2) 
secondderivatives2d (element, x, u, v)
 
real(kind=dp) function interpolateinelement3d (element, x, u, v, w)
 
subroutine nodalbasisfunctions3d (y, element, u, v, w)
 
real(kind=dp) function firstderivativeinu3d (element, x, u, v, w)
 
real(kind=dp) function firstderivativeinv3d (element, x, u, v, w)
 
real(kind=dp) function firstderivativeinw3d (element, x, u, v, w)
 
subroutine nodalfirstderivatives3d (y, element, u, v, w)
 
real(kind=dp) function,
dimension(3, 3) 
secondderivatives3d (element, x, u, v, w)
 
subroutine nodalbasisfunctions (n, Basis, element, u, v, w)
 
subroutine nodalfirstderivatives (n, dLBasisdx, element, u, v, w)
 
recursive function elementinfo (Element, Nodes, u, v, w, detJ, Basis, dBasisdx, ddBasisddx, SecondDerivatives, Bubbles, BasisDegree)
 
real(kind=dp) function elementsize (Element, Nodes)
 
recursive function rtandbdmelementinfo (Element, Nodes, u, v, w, F, detF, Basis, RTBasis, DivRTBasis, BDM, BasisDegree)
 
logical function piolatransformationdata (nn, Element, Nodes, F, DetF, dLBasisdx)
 
real(kind=dp) function,
dimension(3) 
crossproduct (v1, v2)
 
logical function whitneyelementinfo (Element, Basis, dBasisdx, nedges, WhitneyBasis, dWhitneyBasisdx)
 
logical function elementmetric (nDOFs, Elm, Nodes, Metric, DetG, dLBasisdx, LtoGMap)
 
subroutine globalfirstderivativesinternal (elm, nodes, df, gx, gy, gz, Metric, dLBasisdx)
 
subroutine globalfirstderivatives (Elm, Nodes, df, gx, gy, gz, Metric, dLBasisdx)
 
real(kind=dp) function interpolateinelement (elm, f, u, v, w, Basis)
 
subroutine globalsecondderivatives (elm, nodes, f, values, u, v, w, Metric, dBasisdx)
 
integer function, dimension(:,:),
pointer 
lgetedgemap (ElementFamily)
 
real(kind=dp) function elementdiameter (elm, nodes)
 
logical function triangleinside (nx, ny, nz, x, y, z)
 
logical function quadinside (nx, ny, nz, x, y, z)
 
logical function tetrainside (nx, ny, nz, x, y, z)
 
logical function brickinside (nx, ny, nz, x, y, z)
 
subroutine checknormaldirection (Boundary, Normal, x, y, z, turn)
 
real(kind=dp) function,
dimension(3) 
normalvector (Boundary, BoundaryNodes, u0, v0, Check)
 
real(kind=dp) function,
dimension(3) 
surfacevector (Boundary, BoundaryNodes, u, v)
 
real(kind=dp) function linefaceintersection (FaceElement, FaceNodes, Rinit, Rfin, u, v)
 
real(kind=dp) function linefaceintersection2 (FaceElement, FaceNodes, Rinit, Rfin, Intersect)
 
real(kind=dp) function pointfacedistance (BoundaryElement, BoundaryNodes, Coord, Normal, u0, v0)
 
subroutine globaltolocal (u, v, w, x, y, z, Element, ElementNodes)
 
subroutine invertmatrix3x3 (G, GI, detG)
 
integer function, dimension(3) gettrianglefacedirection (Element, FaceMap)
 
integer function, dimension(4) getsquarefacedirection (Element, FaceMap)
 
logical function wedgeordering (ordering)
 
logical function whitney2elementinfo (Element, Basis, dBasisdx, nfaces, WhitneyBasis, dWhitneyBasisdx)
 

Member Function/Subroutine Documentation

subroutine elementdescription::addelementdescription ( type(elementtype_t), target  element,
integer, dimension(:)  BasisTerms 
)

Add an element description to global list of element types.

Parameters
basistermsList of terms in the basis function that should be included for this element type.
elementStructure holding element type description

References compute1dpbasis(), and linearalgebra::invertmatrix().

Referenced by initializeelementdescriptions().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function elementdescription::brickinside ( real(kind=dp), dimension(:)  nx,
real(kind=dp), dimension(:)  ny,
real(kind=dp), dimension(:)  nz,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z 
)

Figure out if given point x,y,z is inside a brick, whose node coordinates are given in nx,ny,nz. Method: Divide to tetrahedrons.

References nx, ny, and tetrainside().

Here is the call graph for this function:

subroutine elementdescription::checknormaldirection ( type(element_t), pointer  Boundary,
real(kind=dp), dimension(3)  Normal,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z,
logical, optional  turn 
)

Normal will point from more dense material to less dense or outwards, if no elements on the other side.

References messages::fatal(), interpolateinelement(), normal(), nx, and ny.

Referenced by freesurface::meancurvature(), freesurface::moveboundary(), and normalvector().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) elementdescription::crossproduct ( real(kind=dp), dimension(3)  v1,
real(kind=dp), dimension(3)  v2 
)

Perform the cross product of two vectors.

real(kind=dp) function elementdescription::elementdiameter ( type(element_t)  elm,
type(nodes_t)  nodes 
)

Figure out element diameter parameter for stablization.

References lgetedgemap().

Referenced by meshutils::displacemesh(), materialmodels::effectiveviscosity(), localmatrix(), meshutils::meshstabparams(), adaptive::refinemesh(), and stabparam().

Here is the call graph for this function:

Here is the caller graph for this function:

recursive function elementdescription::elementinfo ( type(element_t), target  Element,
type(nodes_t)  Nodes,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w,
real(kind=dp)  detJ,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:,:), optional  dBasisdx,
real(kind=dp), dimension(:,:,:), optional  ddBasisddx,
logical, optional  SecondDerivatives,
logical, optional  Bubbles,
integer, dimension(:), optional  BasisDegree 
)

Return basis function values, basis function global first and, if requested, second derivatives at given point in local coordinates.Also return square root of element coordinate system metrics determinant (=sqrt(det(J^TJ))).

Parameters
elementElement structure
nodesElement nodal coordinates.
u1st Local coordinate at which to calculate the basis function.
v2nd local coordinate.
w3rd local coordinate.
detjSquare root of determinant of element coordinate system metric
basisBasis function values at (u,v,w)
dbasisdxGlobal first derivatives of basis functions at (u,v,w)
ddbasisddxGlobal second derivatives of basis functions at (u,v,w) if requested
basisdegreeDegree of each basis function in Basis(:) vector. May be used with P element basis functions
secondderivativesAre the second derivatives needed? (still present for historical reasons)
bubblesAre the bubbles to be avaluated.

References pelementbase::brickbubblepbasis(), pelementbase::brickedgepbasis(), pelementbase::brickfacepbasis(), pelementbase::brickpyraedgepbasis(), coordinatesystems::coordinatesystemdimension(), pelementbase::dbrickbubblepbasis(), pelementbase::dbrickedgepbasis(), pelementbase::dbrickfacepbasis(), pelementbase::dbrickpyraedgepbasis(), pelementbase::dlinebubblepbasis(), pelementbase::dpyramidbubblepbasis(), pelementbase::dpyramidedgepbasis(), pelementbase::dpyramidfacepbasis(), pelementbase::dquadbubblepbasis(), pelementbase::dquadedgepbasis(), pelementbase::dquadpyraedgepbasis(), pelementbase::dtetrabubblepbasis(), pelementbase::dtetraedgepbasis(), pelementbase::dtetrafacepbasis(), pelementbase::dtrianglebubblepbasis(), pelementbase::dtriangleebubblepbasis(), pelementbase::dtriangleedgepbasis(), pelementbase::dwedgebubblepbasis(), pelementbase::dwedgeedgepbasis(), pelementbase::dwedgefacepbasis(), elementmetric(), messages::error(), messages::fatal(), pelementmaps::getbrickedgemap(), pelementmaps::getbrickfacemap(), pelementmaps::getbubbledofs(), getelementtype(), pelementmaps::getpyramidedgemap(), pelementmaps::getpyramidfacemap(), pelementmaps::getquadedgemap(), getsquarefacedirection(), pelementmaps::gettriangleedgemap(), gettrianglefacedirection(), pelementmaps::getwedgeedgemap(), pelementmaps::getwedgefacemap(), globalsecondderivatives(), pelementmaps::ispelement(), pelementbase::linebubblepbasis(), nodalbasisfunctions(), nodalfirstderivatives(), pelementbase::pyramidbubblepbasis(), pelementbase::pyramidedgepbasis(), pelementbase::pyramidfacepbasis(), pelementbase::quadbubblepbasis(), pelementbase::quadedgepbasis(), pelementbase::quadpyraedgepbasis(), pelementbase::tetrabubblepbasis(), pelementbase::tetraedgepbasis(), pelementbase::tetrafacepbasis(), pelementbase::trianglebubblepbasis(), pelementbase::triangleebubblepbasis(), pelementbase::triangleedgepbasis(), pelementbase::wedgebubblepbasis(), pelementbase::wedgeedgepbasis(), pelementbase::wedgefacepbasis(), and wedgeordering().

Referenced by asciipointstomesh(), differentials::axiscurl(), backstressgeneralcompose(), boundaryflux(), boundaryintegrals(), boundarylocalmatrix(), bulkassembly(), bulkintegrals(), solverutils::calculatenodalweights(), particleutils::characteristicelementtime(), compressibilityintegrate(), computeacousticimpedance(), computeaveragevelocity(), computefieldsatpoint(), computemasscenter(), computemissingintegrals(), computepotential(), computepressureerror(), computestress(), computestressandstrain(), computeveloatpoint(), computevelocityerror(), compvecpot(), coordinateintegrals(), coupledvelocitymatrix(), differentials::curl(), dcrboundaryresidual(), dcredgeresidual(), dcrinsideresidual(), diffuseconvectivebboundary(), diffuseconvectiveboundary(), diffuseconvective::diffuseconvectiveboundary(), diffuseconvectivecompose(), diffuseconvective::diffuseconvectivecompose(), diffuseconvectivegenbboundary(), diffuseconvectivegeneral::diffuseconvectivegenboundary(), diffuseconvectivegenboundary(), diffuseconvectivegeneral::diffuseconvectivegencompose(), diffuseconvectivegencompose(), materialmodels::effectiveconductivity(), materialmodels::effectiveviscosity(), elastboundaryresidual(), elastedgeresidual(), elastinsideresidual(), electricboundaryresidual(), electricedgeresidual(), electricforceintegrate(), electricinsideresidual(), elementutils::elementarea(), elementenergy(), elementsize(), epsilonwall(), explicitstabilisationmatrix(), flowboundaryresidual(), flowedgeresidual(), flowinsideresidual(), elementutils::fluxintegrate(), forceintegrate(), materialmodels::frictionheat(), fsiintegration(), generalcurrent(), generalelectricflux(), globaltolocalcoords(), harmmagaxiscompose(), heatboundaryresidual(), heatedgeresidual(), heatinsideresidual(), heavisideintegrate(), inertialmoment(), integconns(), integmassconsistent(), integovera(), integralconstraint(), integratematrix(), integratesource(), differentials::jouleheat(), jouleintegrate(), laplacecompose(), solverutils::laplacematrixassembly(), levelsettimestep(), elementutils::lineintegrate(), defutils::localbcbdofs(), defutils::localbcintegral(), localboundary(), localboundarymatrix(), localbulkmatrix(), localfluxbc(), localinterfacematrix(), localjumps(), localmatrix(), localmatrixbc(), localmatrixboundary(), localreleaserate(), localstress(), differentials::lorentzforce(), lorentzforceave(), lumpedcartesianmass(), lumpedfluidicforce(), lumpedsprings(), magnetodynamics2dharmonic(), magnetodynamicscalcfields(), defutils::mapgausspoints(), solverutils::massmatrixassembly(), maxwellaxis::maxwellaxisboundary(), maxwellaxis::maxwellaxiscompose(), maxwell::maxwellboundary(), maxwell::maxwellcompose(), maxwellgeneral::maxwellgeneralboundary(), maxwellgeneral::maxwellgeneralcompose(), maxwellstresstensorintegrate(), freesurface::meancurvature(), meshboundary(), navierstokes::navierstokesboundary(), navierstokesboundary(), navierstokes::navierstokescompose(), navierstokescylindrical::navierstokescylindricalboundary(), navierstokescylindrical::navierstokescylindricalcompose(), navierstokesgeneral::navierstokesgeneralboundary(), navierstokesgeneral::navierstokesgeneralcompose(), navierstokes::navierstokeswalllaw(), neohookeanlocalmatrix(), particleutils::particleelementinfo(), particledynamicsstuff::particlewallproc(), phasechangesolve(), multigrid::pmgsolve(), poissonboltzmanncompose(), freesurface::poissonsolve(), polylineintegrals(), potential(), pressureintegrate(), pressurelaplacematrix(), pressuremassmatrix(), projecttoplane(), particleutils::randompointinelement(), schurcomplementimpedancematrix(), schurcomplementmatrix(), schurcomplementslipmatrix(), slipmatrix(), stabparam(), statcurrentboundary(), statcurrentcompose(), statelecboundary(), stateleccompose(), statmagaxiscompose(), statmagcartesiancompose(), steadyphasechange(), stressboundary(), stresslocal::stressboundary(), stressboundaryresidual(), stresslocal::stresscompose(), stressedgeresidual(), stresslocal::stressforcecompose(), stressgeneralboundary(), stressgeneral::stressgeneralboundary(), stressgeneral::stressgeneralcompose(), stressgeneralcompose(), stressinsideresidual(), surfacecenterpoints(), surfaceforceintegration(), surfaceimpedanceintegration(), elementutils::surfaceintegrate(), surfacelocalmatrix(), torque(), totalchargebc(), v2f(), velocityimpedancematrix(), velocitylocalmatrix(), velocityslipmatrix(), navierstokes::vmswalls(), elementutils::volumeintegrate(), meshutils::weightedprojector(), meshutils::weightedprojector2(), writevtifile(), and writevtufile().

logical function elementdescription::elementmetric ( integer  nDOFs,
type(element_t)  Elm,
type(nodes_t)  Nodes,
real(kind=dp), dimension(:,:)  Metric,
real(kind=dp)  DetG,
real(kind=dp), dimension(:,:)  dLBasisdx,
real(kind=dp), dimension(3,3)  LtoGMap 
)

Compute contravariant metric tensor (=J^TJ)^-1 of element coordinate system, and square root of determinant of covariant metric tensor (=sqrt(det(J^TJ)))

Parameters
ndofsNumber of active nodes in element
elmElement structure
nodeselement nodal coordinates
metricContravariant metric tensor
dlbasisdxDerivatives of element basis function with respect to local coordinates
detgSQRT of determinant of element coordinate metric
ltogmapMapping between local and global coordinates

References coordinatesystems::coordinatesystemdimension(), messages::error(), g(), messages::info(), and invertmatrix3x3().

Referenced by elementinfo().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function elementdescription::elementsize ( type(element_t)  Element,
type(nodes_t)  Nodes 
)

Returns just the size of the element at its center. providing a more economical way than calling ElementInfo.

References elementinfo(), and messages::fatal().

Referenced by particleutils::characteristicelementtime(), computefieldsatpoint(), computeveloatpoint(), particleutils::initializeparticles(), and putelementsinchildquadrants().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivative1d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u 
)

Given element structure return value of the first partial derivative with respect to local coordinate of a quantity x given at element nodes at local coordinate point u inside the element. Element basis functions are used to compute the value.

Parameters
elementelement structure
uPoint at which to evaluate the partial derivative
xNodal values of the quantity whose partial derivative we want to know

Referenced by globaltolocal(), freesurface::meancurvature(), freesurface::moveboundary(), and normalvector().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivativeinu2d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v 
)

Given element structure return value of the first partial derivative with respect to local coordinate u of i quantity x given at element nodes at local coordinate point u,v inside the element. Element basis functions are used to compute the value.

Referenced by globaltolocal(), freesurface::meancurvature(), freesurface::moveboundary(), and normalvector().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivativeinu3d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Given element structure return value of the first partial derivative with respect to local coordinate u of a quantity x given at element nodes at local coordinate point u,v,w inside the element. Element basis functions are used to compute the value.

Referenced by globaltolocal(), and nodalfirstderivatives().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivativeinv2d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v 
)

Given element structure return value of the first partial derivative with respect to local coordinate v of i quantity x given at element nodes at local coordinate point u,v inside the element. Element basis functions are used to compute the value.

Referenced by globaltolocal(), freesurface::meancurvature(), freesurface::moveboundary(), and normalvector().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivativeinv3d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Given element structure return value of the first partial derivative with respect to local coordinate v of a quantity x given at element nodes at local coordinate point u,v,w inside the element. Element basis functions are used to compute the value.

Referenced by globaltolocal(), and nodalfirstderivatives().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::firstderivativeinw3d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Given element structure return value of the first partial derivatives with respect to local coordinate w of a quantity x given at element nodes at local coordinate point u,v,w inside the element. Element basis functions are used to compute the value.

Referenced by globaltolocal(), and nodalfirstderivatives().

Here is the caller graph for this function:

type(elementtype_t) function, pointer elementdescription::getelementtype ( integer  code,
logical, optional  CompStabFlag 
)

Given element type code return pointer to the corresponding element type structure.

References stabparam(), and messages::warn().

Referenced by compvecpot(), meshutils::createlinemesh(), createsurfaces(), elementinfo(), meshutils::findmeshedges2d(), meshutils::findmeshedges3d(), meshutils::findmeshfaces3d(), defutils::gausspointsboundary(), isosurfacesolver(), meshutils::loadmesh(), defutils::localbcintegral(), localmatrix(), meshutils::meshextrude(), navierstokes::navierstokescompose(), navierstokescylindrical::navierstokescylindricalcompose(), parallelprojecttoplane(), polylineintegrals(), projecttoplane(), resulttopost(), and resulttoresult().

Here is the call graph for this function:

Here is the caller graph for this function:

integer function, dimension(4) elementdescription::getsquarefacedirection ( type(element_t)  Element,
integer, dimension(4)  FaceMap 
)

Given element and its face map (for some square face of element ), this routine returns global direction of square face so that functions are continuous over element boundaries.

Referenced by elementinfo().

Here is the caller graph for this function:

integer function, dimension(3) elementdescription::gettrianglefacedirection ( type(element_t)  Element,
integer, dimension(3)  FaceMap 
)

Given element and its face map (for some triangular face of element ), this routine returns global direction of triangle face so that functions are continuous over element boundaries.

References sort().

Referenced by elementinfo().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::globalfirstderivatives ( type(element_t)  Elm,
type(nodes_t)  Nodes,
real(kind=dp), dimension(:)  df,
real(kind=dp)  gx,
real(kind=dp)  gy,
real(kind=dp)  gz,
real(kind=dp), dimension(:,:)  Metric,
real(kind=dp), dimension(:,:)  dLBasisdx 
)

Given element structure return value of the first partial derivative with respect to global coordinates of a quantity f given at element nodes at local coordinate point u,v,w inside the element. Element basis functions are used to compute the value.

References globalfirstderivativesinternal().

Here is the call graph for this function:

subroutine elementdescription::globalfirstderivativesinternal ( type(element_t)  elm,
type(nodes_t)  nodes,
real(kind=dp), dimension(:)  df,
real(kind=dp)  gx,
real(kind=dp)  gy,
real(kind=dp)  gz,
real(kind=dp), dimension(:,:)  Metric,
real(kind=dp), dimension(:,:)  dLBasisdx 
)

Given element structure return value of the first partial derivatives with respect to global coordinates of a quantity x given at element nodes at local coordinate point u,v,w inside the element. Element basis functions are used to compute the value. This is internal version,and shoudnt usually be called directly by the user, but trough the wrapper routine GlobalFirstDerivatives.

References coordinatesystems::coordinatesystemdimension().

Referenced by globalfirstderivatives().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::globalsecondderivatives ( type(element_t)  elm,
type(nodes_t)  nodes,
real(kind=dp), dimension(:)  f,
real(kind=dp), dimension(:,:)  values,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w,
real(kind=dp), dimension(:,:)  Metric,
real(kind=dp), dimension(:,:), optional  dBasisdx 
)

Compute elementwise matrix of second partial derivatives at given point u,v,w in global coordinates.

References coordinatesystems::coordinatesystemdimension(), secondderivatives1d(), secondderivatives2d(), and secondderivatives3d().

Referenced by elementinfo().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::globaltolocal ( real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z,
type(element_t), pointer  Element,
type(nodes_t)  ElementNodes 
)

Convert global coordinates x,y,z inside element to local coordinates u,v,w of the element.

Todo:
Change to support p elements

References coordinatesystems::coordinatesystemdimension(), firstderivative1d(), firstderivativeinu2d(), firstderivativeinu3d(), firstderivativeinv2d(), firstderivativeinv3d(), firstderivativeinw3d(), interpolateinelement(), generalutils::solvelinsys2x2(), generalutils::solvelinsys3x3(), and messages::warn().

Referenced by boundaryflux(), particledynamicsstuff::particlewallproc(), and interpolation::pointinelement().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::initializeelementdescriptions ( )

Read the element description input file and add the element types to a global list. The file is assumed to be found under the name $ELMER_HOME/lib/elements.def This is the first routine the user of the element utilities should call in his/her code.

References addelementdescription(), and generalutils::readandtrim().

Referenced by elmersolver(), gebhardtfactors(), resulttopost(), resulttoresult(), and viewfactors().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function elementdescription::interpolateinelement ( type(element_t)  elm,
real(kind=dp), dimension(:)  f,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w,
real(kind=dp), dimension(:), optional  Basis 
)

Given element structure return value of a quantity x given at element nodes at local coordinate point u inside the element. Element basis functions are used to compute the value. This is just a wrapper routine and will call the real function according to element dimension.

References interpolateinelement1d(), interpolateinelement2d(), and interpolateinelement3d().

Referenced by checknormaldirection(), globaltolocal(), interpolatemeshtomeshq(), elementutils::lineintegrate(), normalvector(), savescalars(), and surfacevector().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function elementdescription::interpolateinelement1d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u 
)

Given element structure return value of a quantity x given at element nodes at local coordinate point u inside the element. Element basis functions are used to compute the value. This is for 1D elements, and shouldnt propably be called directly by the user but trough the wrapper routine InterpolateInElement.

Parameters
elementelement structure
uPoint at which to evaluate the value
xNodal values of the quantity whose value we want to know

Referenced by interpolateinelement().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::interpolateinelement2d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v 
)

Given element structure return value of a quantity x given at element nodes at local coordinate point (u,vb) inside the element. Element basis functions are used to compute the value.This is for 2D elements, and shouldnt propably be called directly by the user but trough the wrapper routine InterpolateInElement.

Parameters
elementelement structure
uPoint at which to evaluate the partial derivative
vPoint at which to evaluate the partial derivative
xNodal values of the quantity whose partial derivative we want to know

Referenced by interpolateinelement(), and meshutils::splitmeshequal().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::interpolateinelement3d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Given element structure return value of a quantity x given at element nodes at local coordinate point (u,v,w) inside the element. Element basis functions are used to compute the value. This is for 3D elements, and shouldnt propably be called directly by the user but trough the wrapper routine InterpolateInElement.

Referenced by interpolateinelement(), nodalbasisfunctions(), and meshutils::splitmeshequal().

Here is the caller graph for this function:

subroutine elementdescription::invertmatrix3x3 ( real(kind=dp), dimension(3,3)  G,
real(kind=dp), dimension(3,3)  GI,
real(kind=dp)  detG 
)

References g().

Referenced by elementmetric().

Here is the call graph for this function:

Here is the caller graph for this function:

integer function, dimension(:,:), pointer elementdescription::lgetedgemap ( integer  ElementFamily)

Referenced by elementdiameter(), and rtandbdmelementinfo().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::linefaceintersection ( type(element_t), pointer  FaceElement,
type(nodes_t)  FaceNodes,
real(kind=dp), dimension(3)  Rinit,
real(kind=dp), dimension(3)  Rfin,
real(kind=dp), optional  u,
real(kind=dp), optional  v 
)

This subroutine tests where the intersection between the line defined by two points and a plane (or line) defined by a boundary element meet. There is an intersection if ( 0 < Lambda < 1 ). Of all intersections the first one is that with the smallest positive lambda.

References normal(), normalvector(), and surfacevector().

Referenced by particleutils::segmentelementintersection().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function elementdescription::linefaceintersection2 ( type(element_t), pointer  FaceElement,
type(nodes_t)  FaceNodes,
real(kind=dp), dimension(3)  Rinit,
real(kind=dp), dimension(3)  Rfin,
logical  Intersect 
)

This subroutine performs a similar test as above using slightly different strategy.

References linearalgebra::invertmatrix().

Referenced by particleutils::segmentelementinside(), and particleutils::segmentelementintersection2().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::nodalbasisfunctions ( integer  n,
real(kind=dp), dimension(:)  Basis,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

References pelementbase::bricknodalpbasis(), interpolateinelement3d(), pelementmaps::ispbrick(), pelementmaps::ispelement(), pelementmaps::isppyramid(), pelementmaps::ispquad(), pelementmaps::isptetra(), pelementmaps::isptriangle(), pelementmaps::ispwedge(), nodalbasisfunctions1d(), nodalbasisfunctions2d(), nodalbasisfunctions3d(), pelementbase::pyramidnodalpbasis(), pelementbase::quadnodalpbasis(), pelementbase::tetranodalpbasis(), pelementbase::trianglenodalpbasis(), and pelementbase::wedgenodalpbasis().

Referenced by elementinfo(), parallelprojecttoplane(), and projecttoplane().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::nodalbasisfunctions1d ( real(kind=dp), dimension(:)  y,
type(element_t)  element,
real(kind=dp)  u 
)
Parameters
elementelement structure
uPoint at which to evaluate the value
yvalue of the quantity y = x(u)

Referenced by nodalbasisfunctions().

Here is the caller graph for this function:

subroutine elementdescription::nodalbasisfunctions2d ( real(kind=dp), dimension(:)  y,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v 
)
Parameters
elementelement structure
uPoint at which to evaluate the value
vPoint at which to evaluate the value
yvalue of the quantity y = x(u,v)

Referenced by nodalbasisfunctions().

Here is the caller graph for this function:

subroutine elementdescription::nodalbasisfunctions3d ( real(kind=dp), dimension(:)  y,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Referenced by nodalbasisfunctions().

Here is the caller graph for this function:

subroutine elementdescription::nodalfirstderivatives ( integer  n,
real(kind=dp), dimension(:,:)  dLBasisdx,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

References pelementbase::dbricknodalpbasis(), pelementbase::dpyramidnodalpbasis(), pelementbase::dquadnodalpbasis(), pelementbase::dtetranodalpbasis(), pelementbase::dtrianglenodalpbasis(), pelementbase::dwedgenodalpbasis(), firstderivativeinu3d(), firstderivativeinv3d(), firstderivativeinw3d(), pelementmaps::ispbrick(), pelementmaps::ispelement(), pelementmaps::isppyramid(), pelementmaps::ispquad(), pelementmaps::isptetra(), pelementmaps::isptriangle(), pelementmaps::ispwedge(), nodalfirstderivatives1d(), nodalfirstderivatives2d(), and nodalfirstderivatives3d().

Referenced by elementinfo().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine elementdescription::nodalfirstderivatives1d ( real(kind=dp), dimension(:,:)  y,
type(element_t)  element,
real(kind=dp)  u 
)
Parameters
uPoint at which to evaluate the partial derivative
yvalue of the quantity y = /
elementelement structure

Referenced by nodalfirstderivatives().

Here is the caller graph for this function:

subroutine elementdescription::nodalfirstderivatives2d ( real(kind=dp), dimension(:,:)  y,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v 
)

Referenced by nodalfirstderivatives().

Here is the caller graph for this function:

subroutine elementdescription::nodalfirstderivatives3d ( real(kind=dp), dimension(:,:)  y,
type(element_t)  element,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Referenced by nodalfirstderivatives().

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) elementdescription::normalvector ( type(element_t), pointer  Boundary,
type(nodes_t)  BoundaryNodes,
real(kind=dp), optional  u0,
real(kind=dp), optional  v0,
logical, optional  Check 
)

Gives the normal vector of a boundary element. For noncurved elements the normal vector does not depend on the local coordinate while otherwise it does. There are different uses of the function where some do not have the luxury of knowing the local coordinates and hence the center point is used as default.

References checknormaldirection(), messages::fatal(), firstderivative1d(), firstderivativeinu2d(), firstderivativeinv2d(), interpolateinelement(), normal(), nx, and ny.

Referenced by addheatfluxbc(), advectiondiffusionsolver(), solverutils::averageboundarynormals(), boundaryflux(), boundaryintegrals(), boundarylocalmatrix(), bulkassembly(), meshutils::checkinterfacemeshangle(), computeacousticimpedance(), computeaveragevelocity(), computenormaldisplacement(), computepotential(), dcrboundaryresidual(), dcredgeresidual(), diffuseconvectivebboundary(), diffuseconvectivegenbboundary(), dirichletafromb(), elastboundaryresidual(), elastedgeresidual(), electricboundaryresidual(), electricedgeresidual(), electricforceintegrate(), epsilonwall(), flowboundaryresidual(), flowedgeresidual(), elementutils::fluxintegrate(), forceintegrate(), fsiintegration(), heatboundaryresidual(), heatedgeresidual(), helmholtz_smoluchowski_comp(), integmassconsistent(), integratematrix(), linefaceintersection(), elementutils::lineintegrate(), localboundary(), localboundarymatrix(), localfluxbc(), localinterfacematrix(), localjumps(), localmatrix(), localmatrixbc(), localmatrixboundary(), particleutils::locateparticleinmeshmarch(), lumpedfluidicforce(), maxwellaxis::maxwellaxisboundary(), maxwell::maxwellboundary(), maxwellgeneral::maxwellgeneralboundary(), maxwellstresstensorintegrate(), meltingheat(), meshboundary(), movingelstatsolver(), navierstokes::navierstokesboundary(), navierstokesboundary(), navierstokescylindrical::navierstokescylindricalboundary(), navierstokesgeneral::navierstokesgeneralboundary(), navierstokes::navierstokeswalllaw(), particledynamicsstuff::particlewallproc(), pointfacedistance(), polylineintegrals(), pressureintegrate(), setinitialconditions(), slipmatrix(), statelecboundary(), stresslocal::stressboundary(), stressboundaryresidual(), stressedgeresidual(), stressgeneral::stressgeneralboundary(), stressgeneralboundary(), surfaceforceintegration(), surfaceimpedanceintegration(), surfacelocalmatrix(), totalchargebc(), velocityimpedancematrix(), velocitylocalmatrix(), velocityslipmatrix(), viewfactors(), and navierstokes::vmswalls().

Here is the call graph for this function:

logical function elementdescription::piolatransformationdata ( integer  nn,
type(element_t)  Element,
type(nodes_t)  Nodes,
real(kind=dp), dimension(:,:)  F,
real(kind=dp)  DetF,
real(kind=dp), dimension(:,:)  dLBasisdx 
)

This function returns data for performing the Piola transformation.

Parameters
nnThe number of classic nodes used in the element mapping
elementElement structure
nodesData corresponding to the classic element nodes
fThe gradient of the element mapping
detfThe determinant of the gradient matrix (or the Jacobian matrix)
dlbasisdxDerivatives of nodal basis functions with respect to local coordinates

Referenced by rtandbdmelementinfo().

Here is the caller graph for this function:

real(kind=dp) function elementdescription::pointfacedistance ( type(element_t), pointer  BoundaryElement,
type(nodes_t)  BoundaryNodes,
real(kind=dp), dimension(3)  Coord,
real(kind=dp), dimension(3)  Normal,
real(kind=dp), optional  u0,
real(kind=dp), optional  v0 
)

This subroutine computes the signed distance of a point from a surface.

References normal(), normalvector(), and surfacevector().

Referenced by particleutils::getparticleelementintersection(), and particlewallcontact().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function elementdescription::quadinside ( real(kind=dp), dimension(:)  nx,
real(kind=dp), dimension(:)  ny,
real(kind=dp), dimension(:)  nz,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z 
)

Figure out if given point x,y,z is inside a quadrilateral, whose node coordinates are given in nx,ny,nz. Method: Invert the basis functions....

References nx, and ny.

recursive function elementdescription::rtandbdmelementinfo ( type(element_t), target  Element,
type(nodes_t)  Nodes,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w,
real(kind=dp), dimension(3,3)  F,
real(kind=dp)  detF,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:,:)  RTBasis,
real(kind=dp), dimension(:)  DivRTBasis,
logical, optional  BDM,
integer, dimension(:), optional  BasisDegree 
)

Return RT/BDM basis function values and their divergence at the integration point. The data for performing the Piola transformation is also returned. Note that the reference element is chosen as in the p-approximation so that the reference element faces have the same area. This choice simplifies the associated assembly procedure.

Parameters
elementElement structure
nodesData corresponding to the classic element nodes
u1st Local coordinate at which the basis functions are evaluated
v2nd local coordinate
w3rd local coordinate
fThe gradient of the element map K=f(k), with k the reference element
detfThe determinant of the gradient matrix (or the Jacobian matrix)
basisStandard nodal basis functions evaluated at (u,v,w)
rtbasisBasis functions for spanning RT/BDM approximation on reference element
divrtbasisThe divergence of RT/BDM basis functions with respect to the local coordinates
bdmIf .TRUE., a basis for BDM space is constructed
basisdegreeThis a dummy parameter at the moment

References crossproduct(), pelementbase::dtetranodalpbasis(), pelementbase::dtrianglenodalpbasis(), messages::fatal(), lgetedgemap(), piolatransformationdata(), pelementbase::tetranodalpbasis(), and pelementbase::trianglenodalpbasis().

Here is the call graph for this function:

real(kind=dp) function elementdescription::secondderivatives1d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u 
)

Given element structure return value of the second partial derivative with respect to local coordinate of a quantity x given at element nodes at local coordinate point u inside the element. Element basis functions are used to compute the value.

Parameters
elementelement structure
uPoint at which to evaluate the partial derivative
xNodal values of the quantity whose partial derivative we want to know

Referenced by globalsecondderivatives().

Here is the caller graph for this function:

real(kind=dp) function, dimension (2,2) elementdescription::secondderivatives2d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v 
)

Given element structure return value of the second partial derivatives with respect to local coordinates of a quantity x given at element nodes at local coordinate point u,v inside the element. Element basis functions are used to compute the value.

Referenced by globalsecondderivatives().

Here is the caller graph for this function:

real(kind=dp) function, dimension (3,3) elementdescription::secondderivatives3d ( type(element_t)  element,
real(kind=dp), dimension(:)  x,
real(kind=dp)  u,
real(kind=dp)  v,
real(kind=dp)  w 
)

Given element structure return value of the second partial derivatives with respect to local coordinates of i quantity x given at element nodes at local coordinate point u,v inside the element. Element basis functions are used to compute the value.

Referenced by globalsecondderivatives().

Here is the caller graph for this function:

subroutine elementdescription::stabparam ( type(element_t), pointer  Element,
type(nodes_t)  Nodes,
integer  n,
real(kind=dp)  mK,
real(kind=dp), optional  hK 
)

Compute convection diffusion equation stab. parameter for each and every element of the model by solving the largest eigenvalue of.

Lu = Gu, L = (^2 u,^ w), G = ( u, w)

References coordinatesystems::coordinatesystemdimension(), elementdiameter(), elementinfo(), g(), integration::gausspoints(), and messages::info().

Referenced by meshutils::displacemesh(), getelementtype(), meshutils::meshstabparams(), freesurface::moveboundary(), and shearcorrectionfactor().

Here is the call graph for this function:

Here is the caller graph for this function:

real(kind=dp) function, dimension(3) elementdescription::surfacevector ( type(element_t), pointer  Boundary,
type(nodes_t)  BoundaryNodes,
real(kind=dp), optional  u,
real(kind=dp), optional  v 
)

Returns a point that is most importantly supposed to be on the surface For noncurved elements this may simply be the mean while otherwise there may be a need to find the surface node using the local coordinates. Hence the optional parameters. Typically the NormalVector and SurfaceVector should be defined at the same position.

References interpolateinelement(), nx, and ny.

Referenced by linefaceintersection(), and pointfacedistance().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function elementdescription::tetrainside ( real(kind=dp), dimension(:)  nx,
real(kind=dp), dimension(:)  ny,
real(kind=dp), dimension(:)  nz,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z 
)

Figure out if given point x,y,z is inside a tetrahedron, whose node coordinates are given in nx,ny,nz. Method: Invert the basis functions....

References nx, and ny.

Referenced by brickinside().

Here is the caller graph for this function:

logical function elementdescription::triangleinside ( real(kind=dp), dimension(:)  nx,
real(kind=dp), dimension(:)  ny,
real(kind=dp), dimension(:)  nz,
real(kind=dp)  x,
real(kind=dp)  y,
real(kind=dp)  z 
)

Figure out if given point x,y,z is inside a triangle, whose node coordinates are given in nx,ny,nz. Method: Invert the basis functions....

References nx, and ny.

logical function elementdescription::wedgeordering ( integer, dimension(4), intent(in)  ordering)

Function checks if given local numbering of a square face is legal for wedge element.

Referenced by elementinfo().

Here is the caller graph for this function:

logical function elementdescription::whitney2elementinfo ( type(element_t)  Element,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:,:)  dBasisdx,
integer  nfaces,
real(kind=dp), dimension(:,:)  WhitneyBasis,
real(kind=dp), dimension(:,:,:)  dWhitneyBasisdx 
)

Returns Whitney basis vector functions (face elements W2 in 3D) values and Whitney basis function global first derivatives at given point in local coordinates. IMPORTANT NOTE: For performing approximation by using face elements consider utilizing the function RTElementInfo which is newer, offers better extendability and is expected to be the place where future developments are directed.

References coordinatesystems::coordinatesystemdimension(), messages::error(), and messages::fatal().

Referenced by compvecpot().

Here is the call graph for this function:

Here is the caller graph for this function:

logical function elementdescription::whitneyelementinfo ( type(element_t)  Element,
real(kind=dp), dimension(:)  Basis,
real(kind=dp), dimension(:,:)  dBasisdx,
integer  nedges,
real(kind=dp), dimension(:,:)  WhitneyBasis,
real(kind=dp), dimension(:,:,:)  dWhitneyBasisdx 
)

Returns Whitney basis vector functions (edge elements W1 in 3D) values and Whitney basis function global first derivatives at given point in local coordinates.

References coordinatesystems::coordinatesystemdimension(), messages::error(), and messages::fatal().

Referenced by compvecpot(), and magneticw1solver().

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

integer, parameter elementdescription::max_element_nodes = 256

The documentation for this module was generated from the following file: