Elmer FEM solver
Elmer is an open source finite element software for multiphysical problems
|
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) |
subroutine elementdescription::addelementdescription | ( | type(elementtype_t), target | element, |
integer, dimension(:) | BasisTerms | ||
) |
Add an element description to global list of element types.
basisterms | List of terms in the basis function that should be included for this element type. |
element | Structure holding element type description |
References compute1dpbasis(), and linearalgebra::invertmatrix().
Referenced by initializeelementdescriptions().
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().
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().
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().
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))).
element | Element structure |
nodes | Element nodal coordinates. |
u | 1st Local coordinate at which to calculate the basis function. |
v | 2nd local coordinate. |
w | 3rd local coordinate. |
detj | Square root of determinant of element coordinate system metric |
basis | Basis function values at (u,v,w) |
dbasisdx | Global first derivatives of basis functions at (u,v,w) |
ddbasisddx | Global second derivatives of basis functions at (u,v,w) if requested |
basisdegree | Degree of each basis function in Basis(:) vector. May be used with P element basis functions |
secondderivatives | Are the second derivatives needed? (still present for historical reasons) |
bubbles | Are 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)))
ndofs | Number of active nodes in element |
elm | Element structure |
nodes | element nodal coordinates |
metric | Contravariant metric tensor |
dlbasisdx | Derivatives of element basis function with respect to local coordinates |
detg | SQRT of determinant of element coordinate metric |
ltogmap | Mapping between local and global coordinates |
References coordinatesystems::coordinatesystemdimension(), messages::error(), g(), messages::info(), and invertmatrix3x3().
Referenced by elementinfo().
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().
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.
element | element structure |
u | Point at which to evaluate the partial derivative |
x | Nodal values of the quantity whose partial derivative we want to know |
Referenced by globaltolocal(), freesurface::meancurvature(), freesurface::moveboundary(), and normalvector().
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().
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().
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().
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().
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().
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().
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().
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().
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().
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().
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().
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.
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().
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().
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().
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.
element | element structure |
u | Point at which to evaluate the value |
x | Nodal values of the quantity whose value we want to know |
Referenced by interpolateinelement().
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.
element | element structure |
u | Point at which to evaluate the partial derivative |
v | Point at which to evaluate the partial derivative |
x | Nodal values of the quantity whose partial derivative we want to know |
Referenced by interpolateinelement(), and meshutils::splitmeshequal().
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().
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().
integer function, dimension(:,:), pointer elementdescription::lgetedgemap | ( | integer | ElementFamily) |
Referenced by elementdiameter(), and rtandbdmelementinfo().
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().
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().
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().
subroutine elementdescription::nodalbasisfunctions1d | ( | real(kind=dp), dimension(:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u | ||
) |
element | element structure |
u | Point at which to evaluate the value |
y | value of the quantity y = x(u) |
Referenced by nodalbasisfunctions().
subroutine elementdescription::nodalbasisfunctions2d | ( | real(kind=dp), dimension(:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u, | ||
real(kind=dp) | v | ||
) |
element | element structure |
u | Point at which to evaluate the value |
v | Point at which to evaluate the value |
y | value of the quantity y = x(u,v) |
Referenced by nodalbasisfunctions().
subroutine elementdescription::nodalbasisfunctions3d | ( | real(kind=dp), dimension(:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u, | ||
real(kind=dp) | v, | ||
real(kind=dp) | w | ||
) |
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().
subroutine elementdescription::nodalfirstderivatives1d | ( | real(kind=dp), dimension(:,:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u | ||
) |
u | Point at which to evaluate the partial derivative |
y | value of the quantity y = / |
element | element structure |
Referenced by nodalfirstderivatives().
subroutine elementdescription::nodalfirstderivatives2d | ( | real(kind=dp), dimension(:,:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u, | ||
real(kind=dp) | v | ||
) |
subroutine elementdescription::nodalfirstderivatives3d | ( | real(kind=dp), dimension(:,:) | y, |
type(element_t) | element, | ||
real(kind=dp) | u, | ||
real(kind=dp) | v, | ||
real(kind=dp) | w | ||
) |
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().
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.
nn | The number of classic nodes used in the element mapping |
element | Element structure |
nodes | Data corresponding to the classic element nodes |
f | The gradient of the element mapping |
detf | The determinant of the gradient matrix (or the Jacobian matrix) |
dlbasisdx | Derivatives of nodal basis functions with respect to local coordinates |
Referenced by rtandbdmelementinfo().
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().
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 | ||
) |
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.
element | Element structure |
nodes | Data corresponding to the classic element nodes |
u | 1st Local coordinate at which the basis functions are evaluated |
v | 2nd local coordinate |
w | 3rd local coordinate |
f | The gradient of the element map K=f(k), with k the reference element |
detf | The determinant of the gradient matrix (or the Jacobian matrix) |
basis | Standard nodal basis functions evaluated at (u,v,w) |
rtbasis | Basis functions for spanning RT/BDM approximation on reference element |
divrtbasis | The divergence of RT/BDM basis functions with respect to the local coordinates |
bdm | If .TRUE., a basis for BDM space is constructed |
basisdegree | This a dummy parameter at the moment |
References crossproduct(), pelementbase::dtetranodalpbasis(), pelementbase::dtrianglenodalpbasis(), messages::fatal(), lgetedgemap(), piolatransformationdata(), pelementbase::tetranodalpbasis(), and pelementbase::trianglenodalpbasis().
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.
element | element structure |
u | Point at which to evaluate the partial derivative |
x | Nodal values of the quantity whose partial derivative we want to know |
Referenced by globalsecondderivatives().
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().
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().
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().
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().
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....
Referenced by brickinside().
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 | ||
) |
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().
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().
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().
integer, parameter elementdescription::max_element_nodes = 256 |