v0.5.86
Classes | Typedefs | Functions | Variables
MoFEM Namespace Reference

implementation of Data Operators for Forces and Sources More...

Classes

struct  AccelerationCubitBcData
 Definition of the acceleration bc data structure. More...
 
struct  BaseFEDofEntity
 keeps basic information about indexed dofs for the finite element More...
 
struct  BaseFunction
 Base class if inherited used to calculate base functions. More...
 
struct  BaseFunctionCtx
 Base class used to exchange data between element data structures and class calculating base functions. More...
 
struct  BasicEntity
 this struct keeps basic methods for moab entity More...
 
struct  BasicEntityData
 Basic data. like access to moab interface and basic tag handlers. More...
 
struct  BasicMethod
 Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data between MoFEM and user functions. It stores information about multi-indices. More...
 
struct  BasicMethodPtr
 
struct  BitFEId_mi_tag
 
struct  BitFieldId_mi_tag
 
struct  BitFieldId_space_mi_tag
 
struct  BitLevelCoupler
 Interface set parent for vertices, edges, triangles and tetrahedrons. More...
 
struct  BitProblemId_mi_tag
 
struct  BitRefManager
 Managing BitRefLevels . More...
 
struct  Block_BodyForces
 Body force data structure. More...
 
struct  BlockData
 
struct  BlockSetAttributes
 Arbitrary block attributes data structure. More...
 
struct  CfgCubitBcData
 Definition of the cfd_bc data structure. More...
 
struct  ComposedProblemsData
 
struct  Composite_Cubit_msId_And_MeshSetType_mi_tag
 
struct  Composite_Ent_and_ShortId_mi_tag
 
struct  Composite_EntType_and_ParentEntType_mi_tag
 
struct  Composite_EntType_and_Space_mi_tag
 
struct  Composite_mi_tag
 
struct  Composite_Name_And_Ent_And_EntDofIdx_mi_tag
 
struct  Composite_Name_And_Ent_mi_tag
 
struct  Composite_Name_And_Part_mi_tag
 
struct  Composite_Name_And_Type_mi_tag
 
struct  Composite_Name_Ent_And_Part_mi_tag
 
struct  Composite_Name_Ent_Order_And_CoeffIdx_mi_tag
 
struct  Composite_Name_Type_And_Side_Number_mi_tag
 
struct  Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag
 
struct  Composite_ParentEnt_And_EntType_mi_tag
 
struct  Composite_Part_And_Order_mi_tag
 
struct  Composite_SeriesID_And_Step_mi_tag
 
struct  Composite_SeriesName_And_Step_mi_tag
 
struct  Composite_SeriesName_And_Time_mi_tag
 
struct  Composite_Unique_mi_tag
 
struct  CoordSys
 Structure for Coordinate system of two-point tensorScientific computing applications deal in physical quantities expressed as tensors: scalars such as temperature, vectors such as velocity, and second-order tensors such as stress. In practice, these are formally tensor fields: a tensor field assigns a tensor to each point in a mathematical space (typically a Euclidean space or manifold). More...
 
struct  CoordSysName_mi_tag
 
struct  CoordSystemsManager
 
struct  Core
 Core (interface) classThis is the implementation of abstract MoFEM::Interface class. Similarly to the convention used in MoAB, we use small letters to name function of purely abstract classes. This is an exception used only here. For more details about naming functions see Coding practice. More...
 
struct  CreateRowComressedADJMatrix
 Create compressed matrix. More...
 
struct  CubitMeshSets
 this struct keeps basic methods for moab meshset about material and boundary conditions More...
 
struct  CubitMeshSets_change_add_bit_to_cubit_bc_type
 change meshset type More...
 
struct  CubitMeshSets_change_attributes
 
struct  CubitMeshSets_change_attributes_data_structure
 
struct  CubitMeshSets_change_bc_data_structure
 
struct  CubitMeshSets_change_name
 change meshset name More...
 
struct  CubitMeshSets_mask_meshset_mi_tag
 
struct  CubitMeshSets_mi_tag
 MultiIndex Tag for field id. More...
 
struct  CubitMeshSets_name
 
struct  CutMeshInterface
 Interface to cut meshes. More...
 
struct  DataForcesAndSourcesCore
 data structure for finite element entityIt keeps that about indices of degrees of freedom, dofs data, base functions functions, entity side number, type of entities, approximation order, etc. More...
 
struct  DataOperator
 base operator to do operations at Gauss Pt. level More...
 
struct  DefaultElementAdjacency
 default adjacency map More...
 
struct  DerivedDataForcesAndSourcesCore
 this class derive data form other data structure More...
 
struct  DisplacementCubitBcData
 Definition of the displacement bc data structure. More...
 
struct  DMCtx
 PETSc Discrete Manager data structure. More...
 
struct  Dof_shared_ptr_change
 
struct  DofEntity
 keeps information about DOF on the entity More...
 
struct  DofEntity_active_change
 
struct  EdgeElementForcesAndSourcesCore
 Edge finite elementUser is implementing own operator at Gauss points level, by own object derived from EdgeElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing objects to rowOpPtrVector and rowColOpPtrVector. More...
 
struct  EdgePolynomialBase
 Calculate base functions on tetrahedral. More...
 
struct  Ent_Ent_mi_tag
 
struct  Ent_FiniteElement_mi_tag
 
struct  Ent_mi_tag
 
struct  Ent_Owner_mi_tag
 
struct  Ent_ParallelStatus
 
struct  EntDofIdx_mi_tag
 
struct  EntFiniteElement
 Finite element data for entitiy. More...
 
struct  Entity_update_part_proc
 
struct  Entity_update_pcomm_data
 
struct  EntMethod
 Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange data between MoFEM and user functions. It stores information about multi-indices. More...
 
struct  EntPolynomialBaseCtx
 Class used to pass element data to calculate base functions on tet,triangle,edge. More...
 
struct  EntType_mi_tag
 
struct  EqBit
 
struct  FaceElementForcesAndSourcesCore
 Face finite elementUser is implementing own operator at Gauss point level, by own object derived from FaceElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing objects to OpPtrVector. More...
 
struct  FatPrismElementForcesAndSourcesCore
 FatPrism finite elementUser is implementing own operator at Gauss points level, by own object derived from FatPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing objects to rowOpPtrVector and rowColOpPtrVector. More...
 
struct  FatPrismPolynomialBase
 Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (that not perfect) More...
 
struct  FatPrismPolynomialBaseCtx
 Class used to pass element data to calculate base functions on fat prism. More...
 
struct  FEDofEntity
 keeps information about indexed dofs for the finite element More...
 
struct  FEEnt_mi_tag
 
struct  FEMethod
 structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices, residuals, load vectors etc. It is low level class however in some class users looking for speed and efficiency, can use it directly. More...
 
struct  FENumeredDofEntity
 keeps information about indexed dofs for the finite element More...
 
struct  Field
 Provide data structure for (tensor) field approximation.The Field is intended to provide support for fields, with a strong bias towards supporting first and best the capabilities required for scientific computing applications. Since we work with discrete spaces, data structure has to carry information about type of approximation space, its regularity. More...
 
struct  FieldBlas
 Section manager is used to create sections . More...
 
struct  FieldChangeCoordinateSystem
 Set field coordinate system. More...
 
struct  FieldEntity
 Struct keeps handle to entity in the field. More...
 
struct  FieldEntity_change_order
 structure to chane FieldEntity order More...
 
struct  FieldEntityEntFiniteElementAdjacencyMap
 FieldEntityEntFiniteElementAdjacencyMap of mofem finite element and entities. More...
 
struct  FieldEntityEntFiniteElementAdjacencyMap_change_ByWhat
 
struct  FieldName_mi_tag
 MultiIndex Tag for field name. More...
 
struct  FieldSeries
 Structure for recording (time) series. More...
 
struct  FieldSeriesStep
 Structure for keeping time and step. More...
 
struct  FiniteElement
 user adjacency function table More...
 
struct  FiniteElement_change_bit_add
 Add field to data. More...
 
struct  FiniteElement_change_bit_off
 Unset field from data. More...
 
struct  FiniteElement_col_change_bit_add
 Add field to column. More...
 
struct  FiniteElement_col_change_bit_off
 Unset field from column. More...
 
struct  FiniteElement_Meshset_mi_tag
 
struct  FiniteElement_name_mi_tag
 
struct  FiniteElement_row_change_bit_add
 Add field to row. More...
 
struct  FiniteElement_row_change_bit_off
 Unset field from row. More...
 
struct  FlatPrismElementForcesAndSourcesCore
 FlatPrism finite elementUser is implementing own operator at Gauss points level, by own object derived from FlatPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing objects to rowOpPtrVector and rowColOpPtrVector. More...
 
struct  FlatPrismPolynomialBase
 Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (that not perfect) More...
 
struct  FlatPrismPolynomialBaseCtx
 Class used to pass element data to calculate base functions on flat prism. More...
 
struct  ForceCubitBcData
 Definition of the force bc data structure. More...
 
struct  ForcesAndSourcesCore
 structure to get information form mofem into DataForcesAndSourcesCore More...
 
struct  GenericAttributeData
 Generic attribute data structure. More...
 
struct  GenericCubitBcData
 Generic bc data structure. More...
 
struct  HashBit
 
struct  HeatFluxCubitBcData
 Definition of the heat flux bc data structure. More...
 
struct  HelmholtzElement
 struture grouping operators and data used for helmholtz problemsIn order to assemble matrices and right hand side vectors, the loops over elements, enetities over that elememnts and finally loop over intergration points are executed. More...
 
struct  Idx_mi_tag
 
struct  IdxDataType
 
struct  IdxDataTypePtr
 
struct  Interface
 InterfaceThis interface is used by user to:
More...
 
struct  interface_DofEntity
 Interface to DofEntity. More...
 
struct  interface_EntFiniteElement
 interface to EntFiniteElement More...
 
struct  interface_Field
 Pointer interface for MoFEM::Field. More...
 
struct  interface_FieldEntity
 Interface to FieldEntityinterface to FieldEntity. More...
 
struct  interface_FieldSeries
 
struct  interface_FiniteElement
 Inetface for FE. More...
 
struct  interface_NumeredDofEntity
 interface to NumeredDofEntity More...
 
struct  interface_NumeredEntFiniteElement
 interface for NumeredEntFiniteElement More...
 
struct  interface_RefElement
 intrface to RefElement More...
 
struct  interface_RefEntity
 interface to RefEntity More...
 
struct  ISManager
 Section manager is used to create sections . More...
 
struct  JacobiPolynomial
 Calculating Legendre base functions. More...
 
struct  JacobiPolynomialCtx
 Class used to give arguments to Legendre base functions. More...
 
struct  KernelLobattoPolynomial
 Calculating Lobatto base functions. More...
 
struct  KernelLobattoPolynomialCtx
 Class used to give arguments to Kernel Lobatto base functions. More...
 
struct  KeyFromKey
 
struct  KspCtx
 Interface for linear (KSP) solver. More...
 
struct  KspMethod
 data structure for ksp (linear solver) contextStruture stores context data which are set in functions run by PETSc SNES functions. More...
 
struct  LegendrePolynomial
 Calculating Legendre base functions. More...
 
struct  LegendrePolynomialCtx
 Class used to give arguments to Legendre base functions. More...
 
struct  LobattoPolynomial
 Calculating Lobatto base functions. More...
 
struct  LobattoPolynomialCtx
 Class used to give arguments to Lobatto base functions. More...
 
struct  LtBit
 
struct  Mat_Elastic
 Elastic material data structure. More...
 
struct  Mat_Elastic_EberleinHolzapfel1
 Mat_Elastic with Fibres. More...
 
struct  Mat_Elastic_TransIso
 Transverse Isotropic material data structure. More...
 
struct  Mat_Interf
 Linear interface data structure. More...
 
struct  Mat_Moisture
 moisture transport material data structure More...
 
struct  Mat_Thermal
 Thermal material data structure. More...
 
struct  MedInterface
 Interface for load MED files. More...
 
struct  MeshRefinement
 Mesh refinement interface. More...
 
struct  Meshset_mi_tag
 
struct  MeshsetsManager
 Interface for managing meshsets containing materials and boundary conditions. More...
 
struct  MoFEMException
 
struct  MOFEMuuid
 MoFEM interface unique ID. More...
 
struct  NodeMergerInterface
 merge node from two bit levels More...
 
struct  NormElement
 finite element to appeximate analytical solution on surface More...
 
struct  NumeredDofEntity
 keeps information about indexed dofs for the problemFIXME: Is too many iterator, this has to be manage more efficiently, some iterators could be moved to multi_indices views. More...
 
struct  NumeredDofEntity_local_idx_change
 
struct  NumeredDofEntity_mofem_index_change
 
struct  NumeredDofEntity_mofem_part_and_all_index_change
 
struct  NumeredDofEntity_part_change
 
struct  NumeredEntFiniteElement
 Partitioned (Indexed) Finite Element in Problem. More...
 
struct  NumeredEntFiniteElement_change_part
 Change finite element part. More...
 
struct  OpCalculateHdivVectorField
 Get vector field for H-div approximation. More...
 
struct  OpCalculateHdivVectorField_General
 Get vector field for H-div approximation. More...
 
struct  OpCalculateHdivVectorField_General< Tensor_Dim, double, ublas::row_major, DoubleAllacator >
 Get vector field for H-div approximation. More...
 
struct  OpCalculateInvJacForFace
 Calculate inverse of jacobian for face element. More...
 
struct  OpCalculateInvJacForFatPrism
 Calculate inverse of jacobian for face element. More...
 
struct  OpCalculateInvJacForFlatPrism
 Calculate inverse of jacobian for face element. More...
 
struct  OpCalculateScalarFieldGradient
 Get field gradients at integration pts for scalar filed rank 0, i.e. vector field. More...
 
struct  OpCalculateScalarFieldGradient_General
 Evaluate field gradient values for scalar field, i.e. gradient is tensor rank 1 (vector) More...
 
struct  OpCalculateScalarFieldGradient_General< Tensor_Dim, double, ublas::row_major, DoubleAllacator >
 Evaluate field gradient values for scalar field, i.e. gradient is tensor rank 1 (vector), specialization. More...
 
struct  OpCalculateScalarFieldValues
 Get value at integration points for scalar field. More...
 
struct  OpCalculateScalarFieldValues_General
 Calculate field values for tenor field rank 0, i.e. scalar field. More...
 
struct  OpCalculateTensor2FieldValues_General
 Calculate field values for tenor field rank 2. More...
 
struct  OpCalculateTensor2FieldValues_General< Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllacator >
 
struct  OpCalculateVectorFieldGradient
 Get field gradients at integration pts for scalar filed rank 0, i.e. vector field. More...
 
struct  OpCalculateVectorFieldGradient_General
 Evaluate field gradient values for vector field, i.e. gradient is tensor rank 2. More...
 
struct  OpCalculateVectorFieldGradient_General< Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllacator >
 
struct  OpCalculateVectorFieldValues
 Get values at integration pts for tensor filed rank 1, i.e. vector field. More...
 
struct  OpCalculateVectorFieldValues_General
 Calculate field values for tenor field rank 1, i.e. vector field. More...
 
struct  OpCalculateVectorFieldValues_General< Tensor_Dim, double, ublas::row_major, DoubleAllacator >
 Calculate field values (template specialization) for tensor field rank 1, i.e. vector field. More...
 
struct  OpGetCoordsAndNormalsOnFace
 Calculate normals at Gauss points of triangle element. More...
 
struct  OpGetCoordsAndNormalsOnPrism
 calculate normals at Gauss points of triangle element More...
 
struct  OpGetDataAndGradient
 Get field values and gradients at Gauss points. More...
 
struct  OpGetHoTangentOnEdge
 Calculate tangent vector on edge form HO geometry approximation. More...
 
struct  OpSetContravariantPiolaTransform
 apply contravariant (Piola) transfer to Hdiv space More...
 
struct  OpSetContravariantPiolaTransoformOnTriangle
 transform Hdiv base fluxes from reference element to physical triangle More...
 
struct  OpSetCovariantPiolaTransform
 apply covariant transfer to Hcurl space More...
 
struct  OpSetCovariantPiolaTransoformOnEdge
 transform Hcurl base fluxes from reference element to physical edge More...
 
struct  OpSetCovariantPiolaTransoformOnTriangle
 transform Hcurl base fluxes from reference element to physical triangle More...
 
struct  OpSetHoContravariantPiolaTransform
 Apply contravariant (Piola) transfer to Hdiv space for HO geometr. More...
 
struct  OpSetHoCovariantPiolaTransform
 Apply covariant (Piola) transfer to Hcurl space for HO geometry. More...
 
struct  OpSetHoInvJacH1
 transform local reference derivatives of shape function to global derivatives if higher order geometry is given More...
 
struct  OpSetHoInvJacHdivAndHcurl
 transform local reference derivatives of shape function to global derivatives if higher order geometry is given More...
 
struct  OpSetInvJacH1
 Transform local reference derivatives of shape function to global derivatives. More...
 
struct  OpSetInvJacH1ForFace
 Transform local reference derivatives of shape functions to global derivatives. More...
 
struct  OpSetInvJacH1ForFatPrism
 Transform local reference derivatives of shape functions to global derivatives. More...
 
struct  OpSetInvJacH1ForFlatPrism
 Transform local reference derivatives of shape functions to global derivatives. More...
 
struct  OpSetInvJacHcurlFace
 brief Transform local reference derivatives of shape function to global derivatives for face More...
 
struct  OpSetInvJacHdivAndHcurl
 brief Transform local reference derivatives of shape function to global derivatives More...
 
struct  Order_mi_tag
 MultiIndex Tag for field order. More...
 
struct  PairNameFEMethodPtr
 
struct  ParentEntType_mi_tag
 
struct  Part_mi_tag
 
struct  PetscGlobalIdx_mi_tag
 
struct  PetscLocalIdx_mi_tag
 
struct  PressureCubitBcData
 Definition of the pressure bc data structure. More...
 
struct  PrismInterface
 Make interface on given faces and create flat prism in that space. More...
 
struct  PrismsFromSurfaceInterface
 merge node from two bit levels More...
 
struct  Problem
 keeps basic data about problemThis is low level structure with information about problem, what elements compose problem and what DOFs are on rows and columns. More...
 
struct  Problem_mi_tag
 
struct  ProblemChangeRefLevelBitAdd
 add ref level to problem More...
 
struct  ProblemChangeRefLevelBitDofMaskSet
 set prof dof bit ref mask More...
 
struct  ProblemChangeRefLevelBitSet
 set ref level to problem More...
 
struct  ProblemClearNumeredFiniteElementsChange
 clear problem finite elements More...
 
struct  ProblemFiniteElementChangeBitAdd
 add finite element to problem More...
 
struct  ProblemFiniteElementChangeBitUnSet
 remove finite element from problem More...
 
struct  ProblemsManager
 Problem manager is used to build and partition problems . More...
 
struct  ProblemZeroNbColsChange
 zero nb. of DOFs in col More...
 
struct  ProblemZeroNbRowsChange
 zero nb. of DOFs in row More...
 
struct  Proc_mi_tag
 
struct  Projection10NodeCoordsOnField
 Projection of edge entities with one mid-node on hierarchical basis. More...
 
struct  ProjectionFieldOn10NodeTet
 
struct  RefElement
 keeps data about abstract refined finite element More...
 
struct  RefElement_change_parent
 change parentUsing this function with care. Some other multi-indices can deponent on this. More...
 
struct  RefElement_EDGE
 keeps data about abstract EDGE finite element More...
 
struct  RefElement_MESHSET
 keeps data about abstract MESHSET finite element More...
 
struct  RefElement_PRISM
 keeps data about abstract PRISM finite element More...
 
struct  RefElement_TET
 keeps data about abstract TET finite element More...
 
struct  RefElement_TRI
 keeps data about abstract TRI finite element More...
 
struct  RefElement_VERTEX
 keeps data about abstract VERTEX finite element More...
 
struct  RefEntity
 Struct keeps handle to refined handle. More...
 
struct  RefEntity_change_add_bit
 ref mofem entity, change bit More...
 
struct  RefEntity_change_and_bit
 ref mofem entity, change bit More...
 
struct  RefEntity_change_left_shift
 ref mofem entity, left shift More...
 
struct  RefEntity_change_parent
 change parentUse this function with care. Some other multi-indices can deponent on this. More...
 
struct  RefEntity_change_remove_parent
 ref mofem entity, remove parent More...
 
struct  RefEntity_change_right_shift
 ref mofem entity, right shift More...
 
struct  RefEntity_change_set_bit
 ref mofem entity, change bit More...
 
struct  RefEntity_change_set_nth_bit
 ref mofem entity, change bit More...
 
struct  RefEntity_change_xor_bit
 ref mofem entity, change bit More...
 
struct  SeriesID_mi_tag
 
struct  SeriesName_mi_tag
 
struct  SeriesRecorder
 
struct  SideNumber_mi_tag
 
struct  Simple
 Simple interface for fast problem set-up. More...
 
struct  SnesCtx
 Interface for nonlinear (SNES) solver. More...
 
struct  SnesMethod
 data structure for snes (nonlinear solver) contextStructure stores context data which are set in functions run by PETSc SNES functions. More...
 
struct  TemperatureCubitBcData
 Definition of the temperature bc data structure. More...
 
struct  TetGenInterface
 use TetGen to generate mesh More...
 
struct  TetPolynomialBase
 Calculate base functions on tetrahedral. More...
 
struct  TriPolynomialBase
 Calculate base functions on triangle. More...
 
struct  TsCtx
 Interface for Time Stepping (TS) solver. More...
 
struct  TSMethod
 data structure for TS (time stepping) contextStructure stores context data which are set in functions run by PETSc Time Stepping functions. More...
 
struct  Unique_Ent_mi_tag
 
struct  Unique_FiniteElement_mi_tag
 
struct  Unique_mi_tag
 
struct  UnknownInterface
 base class for all interface classes More...
 
struct  VecManager
 Section manager is used to create sections . More...
 
struct  VelocityCubitBcData
 Definition of the velocity bc data structure. More...
 
struct  Version
 
struct  VertexElementForcesAndSourcesCore
 Vertex finite elementUser is implementing own operator at Gauss points level, by own object derived from VertexElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing objects to rowOpPtrVector and rowColOpPtrVector. More...
 
struct  VolumeElementForcesAndSourcesCore
 Volume finite elementUser is implementing own operator at Gauss point level, by class derived from VolumeElementForcesAndSourcesCore::UserDataOperator. Arbitrary number of operator can be added by pushing objects to OpPtrVector. More...
 
struct  VolumeElementForcesAndSourcesCoreOnSide
 Volume element used to integrate on skeleton. More...
 

Typedefs

typedef ErrorCode MoABErrorCode
 
typedef int DofIdx
 Index of DOF. More...
 
typedef int FEIdx
 Index of the element. More...
 
typedef int EntIdx
 Index of DOF on the entity. More...
 
typedef int EntPart
 Partition owning entity. More...
 
typedef double FieldData
 Field data type. More...
 
typedef int ApproximationOrder
 Approximation on the entity. More...
 
typedef int FieldCoefficientsNumber
 Number of field coefficients. More...
 
typedef uint128_t UId
 Unique Id. More...
 
typedef int ShortId
 Unique Id in the field. More...
 
typedef std::bitset< BITREFEDGES_SIZEBitRefEdges
 
typedef std::bitset< BITREFLEVEL_SIZEBitRefLevel
 Bit structure attached to each entity identifying to what mesh entity is attached. More...
 
typedef std::bitset< BITFIELDID_SIZEBitFieldId
 
typedef std::bitset< BITFEID_SIZEBitFEId
 
typedef std::bitset< BITPROBLEMID_SIZEBitProblemId
 
typedef std::bitset< BITINTERFACEUID_SIZEBitIntefaceId
 
typedef std::bitset< 32 > CubitBCType
 
typedef std::vector< int, std::allocator< int > > IntAllacator
 
typedef std::vector< double, std::allocator< double > > DoubleAllacator
 
typedef ublas::vector< int, IntAllacatorVectorInt
 
typedef ublas::vector< double, DoubleAllacatorVectorDouble
 
typedef ublas::matrix< double, ublas::row_major, DoubleAllacatorMatrixDouble
 
typedef ublas::matrix< double, ublas::row_major, ublas::bounded_array< double, 9 > > MatrixDouble3by3
 
typedef ublas::vector< double, ublas::bounded_array< double, 3 > > VectorDouble3
 
typedef ublas::vector< double, ublas::bounded_array< double, 9 > > VectorDouble9
 
typedef ublas::vector< double, ublas::shallow_array_adaptor< double > > VectorAdaptor
 
typedef ublas::matrix< double, ublas::row_major, ublas::shallow_array_adaptor< double > > MatrixAdaptor
 
typedef ublas::vector< int, ublas::shallow_array_adaptor< int > > VectorIntAdaptor
 
typedef std::vector< boost::shared_ptr< MatrixDouble > > ShapeFunctionBasesVector
 
typedef ublas::unbounded_array< boost::shared_ptr< const FEDofEntity >, std::allocator< boost::shared_ptr< const FEDofEntity > >> DofsAllacator
 
typedef ublas::vector< boost::shared_ptr< const FEDofEntity >, DofsAllacatorVectorDofs
 
typedef CubitMeshSet_multiIndex::index< CubitMeshSets_mi_tag >::type CubitMeshsetByType
 
typedef CubitMeshSet_multiIndex::index< CubitMeshSets_mask_meshset_mi_tag >::type CubitMeshsetByMask
 
typedef CubitMeshSet_multiIndex::index< CubitMeshSets_name >::type CubitMeshsetByName
 
typedef CubitMeshSet_multiIndex::index< CubitMeshSets_mi_tag >::type CubitMeshsetById
 
typedef multi_index_container< CubitMeshSets, indexed_by< hashed_unique< tag< Meshset_mi_tag >, member< CubitMeshSets, EntityHandle,&CubitMeshSets::meshset > >, ordered_non_unique< tag< CubitMeshSets_mi_tag >, const_mem_fun< CubitMeshSets, unsigned long int,&CubitMeshSets::getBcTypeULong > >, ordered_non_unique< tag< CubitMeshSets_mask_meshset_mi_tag >, const_mem_fun< CubitMeshSets, unsigned long int,&CubitMeshSets::getMaksedBcTypeULong > >, ordered_non_unique< tag< CubitMeshSets_name >, const_mem_fun< CubitMeshSets, std::string,&CubitMeshSets::getName > >, hashed_unique< tag< Composite_Cubit_msId_And_MeshSetType_mi_tag >, composite_key< CubitMeshSets, const_mem_fun< CubitMeshSets, int,&CubitMeshSets::getMeshsetId >, const_mem_fun< CubitMeshSets, unsigned long int,&CubitMeshSets::getMaksedBcTypeULong > > > >> CubitMeshSet_multiIndex
 Stores data about meshsets (see CubitMeshSets) storing data about boundary conditions, interfaces, sidesets, nodests, blocksets. More...
 
typedef multi_index_container< boost::shared_ptr< CoordSys >, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< CoordSys, EntityHandle,&CoordSys::meshSet > >, ordered_unique< tag< CoordSysName_mi_tag >, const_mem_fun< CoordSys, boost::string_ref,&CoordSys::getNameRef > > > > CoordSys_multiIndex
 
typedef DofEntity_multiIndex::index< FieldName_mi_tag >::type DofEntityByFieldName
 Dof entity multi-index by field name. More...
 
typedef DofEntity_multiIndex::index< Composite_Name_And_Ent_mi_tag >::type DofEntityByNameAndEnt
 Dof entity multi-index by field name and entity. More...
 
typedef DofEntity_multiIndex::index< Composite_Name_And_Type_mi_tag >::type DofEntityByNameAndType
 Dof entity multi-index by field name and entity type. More...
 
typedef multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< const_mem_fun< DofEntity, UId,&DofEntity::getGlobalUniqueId > > > > DofEntity_multiIndex_uid_view
 multi-index view on DofEntity by uid More...
 
typedef multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< const_mem_fun< DofEntity, UId,&DofEntity::getGlobalUniqueId > >, ordered_non_unique< const_mem_fun< DofEntity, char,&DofEntity::getActive > > > > DofEntity_multiIndex_active_view
 multi-index view on DofEntity activity More...
 
typedef multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_non_unique< const_mem_fun< DofEntity, ApproximationOrder,&DofEntity::getDofOrder > > > > DofEntity_multiIndex_order_view
 multi-index view on DofEntity order More...
 
typedef multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_non_unique< const_mem_fun< DofEntity::interface_type_RefEntity, EntityType,&DofEntity::getEntType > > > > DofEntity_multiIndex_ent_type_view
 multi-index view on DofEntity type More...
 
typedef FEDofEntity_multiIndex::index< FieldName_mi_tag >::type FEDofEntityByFieldName
 Finite element DoF multi-index by field name. More...
 
typedef FEDofEntity_multiIndex::index< Composite_Name_And_Ent_mi_tag >::type FEDofEntityByNameAndEnt
 Dof entity multi-index by field name and entity. More...
 
typedef FEDofEntity_multiIndex::index< Composite_Name_And_Type_mi_tag >::type FEDofEntityByNameAndType
 Dof entity multi-index by field name and entity type. More...
 
typedef FENumeredDofEntity_multiIndex::index< FieldName_mi_tag >::type FENumeredDofEntityByFieldName
 Finite element numbered DoF multi-index by field name. More...
 
typedef FENumeredDofEntity_multiIndex::index< Composite_Name_And_Ent_mi_tag >::type FENumeredDofEntityByNameAndEnt
 Dof entity multi-index by field name and entity. More...
 
typedef FENumeredDofEntity_multiIndex::index< Composite_Name_And_Type_mi_tag >::type FENumeredDofEntityByNameAndType
 Dof entity multi-index by field name and entity type. More...
 
typedef FENumeredDofEntity_multiIndex::index< Unique_mi_tag >::type FENumeredDofEntityByUId
 Dof entity multi-index by UId. More...
 
typedef FENumeredDofEntity_multiIndex::index< Ent_mi_tag >::type FENumeredDofEntityByEnt
 Numbered DoF multi-index by entity. More...
 
typedef NumeredDofEntity_multiIndex::index< FieldName_mi_tag >::type NumeredDofEntityByFieldName
 Numbered DoF multi-index by field name. More...
 
typedef NumeredDofEntity_multiIndex::index< Unique_mi_tag >::type NumeredDofEntityByUId
 Numbered DoF multi-index by UId. More...
 
typedef NumeredDofEntity_multiIndex::index< PetscLocalIdx_mi_tag >::type NumeredDofEntityByLocalIdx
 Numbered DoF multi-index by local index. More...
 
typedef NumeredDofEntity_multiIndex::index< Ent_mi_tag >::type NumeredDofEntityByEnt
 Numbered DoF multi-index by entity. More...
 
typedef NumeredDofEntity_multiIndex::index< Composite_Name_Ent_And_Part_mi_tag >::type NumeredDofEntityByNameEntAndPart
 Numbered DoF multi-index by name entity and partition. More...
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< const_mem_fun< NumeredDofEntity::interface_type_DofEntity, const UId, &NumeredDofEntity::getGlobalUniqueId > > >> NumeredDofEntity_multiIndex_uid_view_ordered
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< hashed_unique< const_mem_fun< NumeredDofEntity, DofIdx,&NumeredDofEntity::getDofIdx > > >> NumeredDofEntity_multiIndex_idx_view_hashed
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_non_unique< const_mem_fun< NumeredDofEntity, DofIdx,&NumeredDofEntity::getPetscLocalDofIdx > > >> NumeredDofEntity_multiIndex_petsc_local_dof_view_ordered_non_unique
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_non_unique< const_mem_fun< NumeredDofEntity, DofIdx,&NumeredDofEntity::getPetscGlobalDofIdx > > >> NumeredDofEntity_multiIndex_petsc_global_dof_view_ordered_non_unique
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_non_unique< const_mem_fun< NumeredDofEntity::interface_type_DofEntity, FieldCoefficientsNumber,&NumeredDofEntity::getDofCoeffIdx > > >> NumeredDofEntity_multiIndex_coeff_idx_ordered_non_unique
 
typedef multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< member< NumeredDofEntity, const DofIdx,&NumeredDofEntity::petscGloablDofIdx > > > > NumeredDofEntity_multiIndex_global_index_view
 
typedef multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, member< RefEntity::BasicEntity, EntityHandle,&RefEntity::ent > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle,&RefEntity::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< RefEntity::BasicEntity, EntityType,&RefEntity::getEntType > >, ordered_non_unique< tag< ParentEntType_mi_tag >, const_mem_fun< RefEntity, EntityType,&RefEntity::getParentEntType > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity::BasicEntity, EntityType,&RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType,&RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityHandle,&RefEntity::getParentEnt >, const_mem_fun< RefEntity::BasicEntity, EntityType,&RefEntity::getEntType > > > >> RefEntity_multiIndex
 
typedef multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< hashed_unique< const_mem_fun< RefEntity, EntityHandle,&RefEntity::getParentEnt > >, hashed_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< boost::shared_ptr< RefEntity >, const_mem_fun< RefEntity, EntityHandle,&RefEntity::getRefEnt >, const_mem_fun< RefEntity, EntityHandle,&RefEntity::getParentEnt > > > >> RefEntity_multiIndex_view_by_parent_entity
 multi-index view of RefEntity by parent entity More...
 
typedef FieldEntity_multiIndex::index< FieldName_mi_tag >::type FieldEntityByFieldName
 Entity nulti index by field name. More...
 
typedef multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< sequenced<>, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity, EntityHandle,&FieldEntity::getEnt > > >> FieldEntity_multiIndex_ent_view
 
typedef interface_RefElement< RefElementptrWrapperRefElement
 
typedef multi_index_container< ptrWrapperRefElement, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityHandle, &ptrWrapperRefElement::getRefEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityHandle, &ptrWrapperRefElement::getParentEnt > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityType, &ptrWrapperRefElement::getEntType > >, ordered_non_unique< tag< Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag >, composite_key< ptrWrapperRefElement, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityHandle, &ptrWrapperRefElement::getParentEnt >, const_mem_fun< ptrWrapperRefElement::interface_type_RefElement, int, &ptrWrapperRefElement::getBitRefEdgesUlong > > >, hashed_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< ptrWrapperRefElement, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityHandle, &ptrWrapperRefElement::getRefEnt >, const_mem_fun< ptrWrapperRefElement::interface_type_RefEntity, EntityHandle, &ptrWrapperRefElement::getParentEnt > > > >> RefElement_multiIndex
 
typedef boost::function< PetscErrorCode(Interface &moab, const Field &field_ptr, const EntFiniteElement &fe_ptr, Range &adjacency)> ElementAdjacencyFunct
 user adjacency function More...
 
typedef EntFiniteElement_multiIndex::index< FiniteElement_name_mi_tag >::type EntFiniteElementbyName
 Entity finite element multi-index by finite element name. More...
 
typedef NumeredEntFiniteElement_multiIndex::index< FiniteElement_name_mi_tag >::type NumeredEntFiniteElementbyName
 Entity finite element multi-index by finite element name. More...
 
typedef NumeredEntFiniteElement_multiIndex::index< Composite_Name_And_Part_mi_tag >::type NumeredEntFiniteElementbyNameAndPart
 Entity finite element multi-index by finite element name and partition. More...
 
typedef int(* FieldOrderTable[MBMAXTYPE]) (const int order)
 user adjacency function table More...
 
typedef int(* FieldOrderFunct) (const int order)
 user adjacency function More...
 
typedef multi_index_container< boost::shared_ptr< Field >, indexed_by< ordered_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &,&Field::getId >, LtBit< BitFieldId > >> > Field_multiIndex_view
 
typedef multi_index_container< FieldSeries, indexed_by< ordered_unique< tag< SeriesID_mi_tag >, const_mem_fun< FieldSeries, EntityID,&FieldSeries::get_meshset_id > >, ordered_unique< tag< SeriesName_mi_tag >, const_mem_fun< FieldSeries, boost::string_ref,&FieldSeries::getNameRef > > >> Series_multiIndex
 Series multi index. More...
 
typedef multi_index_container< FieldSeriesStep, indexed_by< ordered_unique< tag< Composite_SeriesID_And_Step_mi_tag >, composite_key< FieldSeriesStep, const_mem_fun< FieldSeriesStep::interface_type_FieldSeries, EntityID, &FieldSeriesStep::get_meshset_id >, member< FieldSeriesStep, int,&FieldSeriesStep::step_number > > >, ordered_unique< tag< Composite_SeriesName_And_Step_mi_tag >, composite_key< FieldSeriesStep, const_mem_fun< FieldSeriesStep::interface_type_FieldSeries, boost::string_ref, &FieldSeriesStep::getNameRef >, member< FieldSeriesStep, int,&FieldSeriesStep::step_number > > >, ordered_non_unique< tag< SeriesName_mi_tag >, const_mem_fun< FieldSeriesStep::interface_type_FieldSeries, boost::string_ref, &FieldSeriesStep::getNameRef > >, ordered_non_unique< tag< Composite_SeriesName_And_Time_mi_tag >, composite_key< FieldSeriesStep, const_mem_fun< FieldSeriesStep::interface_type_FieldSeries, boost::string_ref, &FieldSeriesStep::getNameRef >, const_mem_fun< FieldSeriesStep, double,&FieldSeriesStep::get_time > > > >> SeriesStep_multiIndex
 Step multi index. More...
 
typedef std::vector< PairNameFEMethodPtrFEMethodsSequence
 
typedef std::vector< BasicMethodPtrBasicMethodsSequence
 

Functions

template<class X >
std::string toString (X x)
 
PetscErrorCode PetscOptionsGetInt (PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
 
PetscErrorCode PetscOptionsGetReal (PetscOptions *, const char pre[], const char name[], PetscReal *dval, PetscBool *set)
 
PetscErrorCode PetscOptionsGetScalar (PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
 
PetscErrorCode PetscOptionsGetString (PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
 
PetscErrorCode PetscOptionsGetBool (PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
 
PetscErrorCode PetscOptionsGetRealArray (PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
 
PetscErrorCode PetscOptionsGetEList (PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
 
PetscErrorCode PetscOptionsGetIntArray (PetscOptions options, const char pre[], const char name[], PetscInt dvalue[], PetscInt *nmax, PetscBool *set)
 
PetscErrorCode PetscOptionsGetScalarArray (PetscOptions options, const char pre[], const char name[], PetscScalar dvalue[], PetscInt *nmax, PetscBool *set)
 
void tet_type_6 (moab::Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
int tet_type_5 (moab::Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
int tet_type_4 (const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
int tet_type_3 (const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
int tet_type_2 (const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
void tet_type_1 (const EntityHandle *conn, const int split_edge, const EntityHandle edge_new_node, EntityHandle *new_tets_conn)
 
PetscErrorCode tri_type_3 (const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tris_conn)
 
PetscErrorCode prism_type_1 (const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
 
PetscErrorCode prism_type_2 (const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
 
PetscErrorCode prism_type_3 (const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
 
PetscErrorCode Hcurl_EdgeBaseFunctions_MBTET (int *sense, int *p, double *N, double *diffN, double *edgeN[], double *diff_edgeN[], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Edge based H-curl base functions on tetrahedral. More...
 
PetscErrorCode Hcurl_EdgeBaseFunctions_MBTET_ON_EDGE (int sense, int p, double *N, double *diffN, double *edgeN, double *diff_edgeN, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Edge based H-curl base functions on edge. More...
 
PetscErrorCode Hcurl_EdgeBaseFunctions_MBTET_ON_FACE (int *sense, int *p, double *N, double *diffN, double *edgeN[], double *diff_edgeN[], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Edge based H-curl base functions on face. More...
 
PetscErrorCode Hcurl_EdgeBasedFaceFunctions_MBTET (int *faces_nodes, int *p, double *N, double *diffN, double *phi_f_e[4][3], double *diff_phi_f_e[4][3], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face edge base functions of Hcurl space on tetrahedral. More...
 
PetscErrorCode Hcurl_EdgeBasedFaceFunctions_MBTET_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face edge base functions of Hcurl space. More...
 
PetscErrorCode Hcurl_BubbleFaceFunctions_MBTET (int *faces_nodes, int *p, double *N, double *diffN, double *phi_f[4], double *diff_phi_f[4], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face edge base functions of Hcurl space on face on tetrahedral. More...
 
PetscErrorCode Hcurl_BubbleFaceFunctions_MBTET_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face edge base functions of Hcurl space on face. More...
 
PetscErrorCode Hcurl_FaceInteriorFunctions_MBTET (int *faces_nodes, int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face base interior function. More...
 
PetscErrorCode Hcurl_VolumeInteriorFunctions_MBTET (int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Volume interior function. More...
 
PetscErrorCode Hcurl_FaceFunctions_MBTET (int *face_nodes, int *p, double *N, double *diffN, double *phi_f[4], double *diff_phi_f[4], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face H-curl functions. More...
 
PetscErrorCode Hcurl_FaceFunctions_MBTET_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face H-curl functions. More...
 
PetscErrorCode Hcurl_VolumeFunctions_MBTET (int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 H-curl volume base functions. More...
 
PetscErrorCode Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET (int *faces_nodes, int *p, double *N, double *diffN, double *phi_f_e[4][3], double *diff_phi_f_e[4][3], int gdim, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Hdiv base functions, Edge-based face functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Hdiv base functions, Edge-based face functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Ainsworth_FaceBubbleShapeFunctions (int *faces_nodes, int *p, double *N, double *diffN, double *phi_f[], double *diff_phi_f[], int gdim, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face bubble functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Face bubble functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Ainsworth_EdgeBasedVolumeShapeFunctions_MBTET (int p, double *N, double *diffN, double *phi_v_e[6], double *diff_phi_v_e[6], int gdim, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Hdiv base function, Edge-based interior (volume) functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Ainsworth_FaceBasedVolumeShapeFunctions_MBTET (int p, double *N, double *diffN, double *phi_v_f[], double *diff_phi_v_f[], int gdim, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 
PetscErrorCode Hdiv_Ainsworth_VolumeBubbleShapeFunctions_MBTET (int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int gdim, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Interior bubble functions by Ainsworth [1]. More...
 
PetscErrorCode Hdiv_Demkowicz_Face_MBTET_ON_FACE (int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb)
 
PetscErrorCode Hdiv_Demkowicz_Interior_MBTET (int p, double *N, double *diffN, int p_face[], double *phi_f[4], double *diff_phi_f[4], double *phi_v, double *diff_phi_v, int gdim)
 
void tet_type_6 (Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
int tet_type_5 (Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
 
template<int Tensor_Dim, class T , class L , class A >
PetscErrorCode invertTensor3by3 (ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
 Calculate inverse of tensor rank 2 at integration points. More...
 
template<>
PetscErrorCode invertTensor3by3< 3, double, ublas::row_major, DoubleAllacator > (MatrixDouble &jac_data, VectorDouble &det_data, MatrixDouble &inv_jac_data)
 
template<class T1 , class T2 >
PetscErrorCode determinantTensor3by3 (T1 &t, T2 &det)
 Calculate determinant. More...
 
template<class T1 , class T2 , class T3 >
PetscErrorCode invertTensor3by3 (T1 &t, T2 &det, T3 &inv_t)
 Calculate matrix inverse. More...
 
template<class T , class A >
FTensor::Tensor0< T * > getTensor0FormData (ublas::vector< T, A > &data)
 Get tensor rank 0 (scalar) form data vectorExample how to use it. More...
 
template<>
FTensor::Tensor0< double * > getTensor0FormData< double, DoubleAllacator > (ublas::vector< double, DoubleAllacator > &data)
 
template<int Tensor_Dim, class T , class L , class A >
FTensor::Tensor1< T *, Tensor_Dim > getTensor1FormData (ublas::matrix< T, L, A > &data)
 Get tensor rank 1 (vector) form data matrix. More...
 
template<int Tensor_Dim>
FTensor::Tensor1< double *, Tensor_Dim > getTensor1FormData (MatrixDouble &data)
 Get tensor rank 1 (vector) form data matrix (specialization) More...
 
template<>
FTensor::Tensor1< double *, 3 > getTensor1FormData< 3, double, ublas::row_major, DoubleAllacator > (MatrixDouble &data)
 
template<>
FTensor::Tensor1< double *, 2 > getTensor1FormData< 2, double, ublas::row_major, DoubleAllacator > (MatrixDouble &data)
 
template<int Tensor_Dim0, int Tensor_Dim1, class T , class L , class A >
FTensor::Tensor2< T *, Tensor_Dim0, Tensor_Dim1 > getTensor2FormData (ublas::matrix< T, L, A > &data)
 Get tensor rank 2 (matrix) form data matrix. More...
 
template<>
FTensor::Tensor2< double *, 3, 3 > getTensor2FormData (MatrixDouble &data)
 Get tensor rank 2 (matrix) form data matrix (specialization) More...
 
template<>
FTensor::Tensor2< double *, 3, 2 > getTensor2FormData (MatrixDouble &data)
 Get tensor rank 2 (matrix) form data matrix (specialization) More...
 
template<int Tensor_Dim0, int Tensor_Dim1>
FTensor::Tensor2< double *, Tensor_Dim0, Tensor_Dim1 > getTensor2FormData (MatrixDouble &data)
 Get tensor rank 2 (matrix) form data matrix (specialization) More...
 
template<>
FTensor::Tensor2< double *, 3, 3 > getTensor2FormData< 3, 3, double, ublas::row_major, DoubleAllacator > (MatrixDouble &data)
 
template<>
FTensor::Tensor2< double *, 3, 2 > getTensor2FormData< 3, 2, double, ublas::row_major, DoubleAllacator > (MatrixDouble &data)
 
template<class T >
void cOnstructor (DataForcesAndSourcesCore *data, EntityType type, T)
 
std::ostream & operator<< (std::ostream &os, const DataForcesAndSourcesCore::EntData &e)
 
std::ostream & operator<< (std::ostream &os, const DataForcesAndSourcesCore &e)
 
template<typename DOFMULTIINDEX >
static int getMaxOrder (const DOFMULTIINDEX &dof_multi_index)
 
static PetscErrorCode get_porblem_row_indices (const ForcesAndSourcesCore *fe_ptr, const EntityType type, const int side, const std::string field_name, VectorInt &indices)
 
static PetscErrorCode get_porblem_col_indices (const ForcesAndSourcesCore *fe_ptr, const EntityType type, const int side, const std::string field_name, VectorInt &indices)
 
static void error_printf_hilight (void)
 
static void error_printf_normal (void)
 
static PetscErrorCode mofem_error_handler (MPI_Comm comm, int line, const char *fun, const char *file, PetscErrorCode n, PetscErrorType p, const char *mess, void *ctx)
 
PetscErrorCode print_verison ()
 
static med_geometrie_element moab2med_element_type (const EntityType type)
 
std::ostream & operator<< (std::ostream &os, const MedInterface::FieldData &field_data)
 
template<class T >
static double determinant (T &t)
 
template<class T >
static double volume_length_quality (T *coords)
 
int fNBENTITY_GENERIC (int P)
 
int fNBENTITYSET_NOFIELD (int P)
 
int fNBVERTEX_L2 (int P)
 
int fNBVOLUMETET_L2 (int P)
 
int fNBFACETRI_L2 (int P)
 
int fNBEDGE_L2 (int P)
 
int fNBVERTEX_H1 (int P)
 number of approx. functions for H1 space on vertex More...
 
int fNBEDGE_H1 (int P)
 number of approx. functions for H1 space on edge More...
 
int fNBFACETRI_H1 (int P)
 number of approx. functions for H1 space on face More...
 
int fNBFACEQUAD_H1 (int P)
 
int fNBVOLUMETET_H1 (int P)
 number of approx. functions for H1 space on volume More...
 
int fNBVOLUMEPRISM_H1 (int P)
 
int fNBVERTEX_HCURL (int P)
 number of approx. functions for HCURL space on vertex More...
 
int fNBEDGE_HCURL (int P)
 
int fNBFACETRI_HCURL (int P)
 
int fNBVOLUMETET_HCURL (int P)
 
int fNBVERTEX_HDIV (int P)
 number of approx. functions for HDIV space on vertex More...
 
int fNBEDGE_HDIV (int P)
 number of approx. functions for HDIV space on edge More...
 
int fNBFACETRI_AINSWORTH_HDIV (int P)
 number of approx. functions for HDIV space on face More...
 
int fNBVOLUMETET_AINSWORTH_HDIV (int P)
 number of approx. functions for HDIV space on volume More...
 
int fNBFACETRI_DEMKOWICZ_HDIV (int P)
 
int fNBVOLUMETET_DEMKOWICZ_HDIV (int P)
 number of approx. functions for HDIV space on volume More...
 
PetscErrorCode test_moab (Interface &moab, const EntityHandle ent)
 Test MoAB entity handle if has structure as is assumed by MoFEM. More...
 
struct __attribute__ ((__packed__)) SideNumber
 keeps information about side number for the finite element More...
 
std::ostream & operator<< (std::ostream &os, const CubitMeshSets &e)
 
std::ostream & operator<< (std::ostream &os, const DisplacementCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const ForceCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const VelocityCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const AccelerationCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const TemperatureCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const PressureCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const HeatFluxCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const CfgCubitBcData &e)
 
std::ostream & operator<< (std::ostream &os, const BlockSetAttributes &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Elastic &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Elastic_EberleinHolzapfel1 &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Thermal &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Moisture &e)
 
std::ostream & operator<< (std::ostream &os, const Block_BodyForces &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Elastic_TransIso &e)
 
std::ostream & operator<< (std::ostream &os, const Mat_Interf &e)
 
std::ostream & operator<< (std::ostream &os, const Field &e)
 
std::ostream & operator<< (std::ostream &os, const FieldEntityEntFiniteElementAdjacencyMap &e)
 
std::ostream & operator<< (std::ostream &os, const DofEntity &e)
 
std::ostream & operator<< (std::ostream &os, const NumeredDofEntity &e)
 
std::ostream & operator<< (std::ostream &os, const FEDofEntity &e)
 
std::ostream & operator<< (std::ostream &os, const FENumeredDofEntity &e)
 
void * get_tag_ptr (SequenceManager *sequence_manager, Tag th, EntityHandle ent, int *tag_size)
 
PetscErrorCode getPatentEnt (Interface &moab, Range ents, std::vector< EntityHandle > vec_patent_ent)
 
std::ostream & operator<< (std::ostream &os, const RefEntity &e)
 
std::ostream & operator<< (std::ostream &os, const FieldEntity &e)
 
std::ostream & operator<< (std::ostream &os, const RefElement &e)
 
std::ostream & operator<< (std::ostream &os, const RefElement_TET &e)
 
std::ostream & operator<< (std::ostream &os, const RefElement_TRI &e)
 
std::ostream & operator<< (std::ostream &os, const RefElement_EDGE &e)
 
std::ostream & operator<< (std::ostream &os, const RefElement_VERTEX &e)
 
std::ostream & operator<< (std::ostream &os, const FiniteElement &e)
 
std::ostream & operator<< (std::ostream &os, const EntFiniteElement &e)
 
template<typename FE_DOFS , typename MOFEM_DOFS , typename MOFEM_DOFS_VIEW >
static PetscErrorCode get_fe_dof_view (const FE_DOFS &fe_dofs_view, const MOFEM_DOFS &mofem_dofs, MOFEM_DOFS_VIEW &mofem_dofs_view, const int operation_type)
 
std::ostream & operator<< (std::ostream &os, const Problem &e)
 
std::ostream & operator<< (std::ostream &os, const FieldSeries &e)
 
std::ostream & operator<< (std::ostream &os, const FieldSeriesStep &e)
 
PetscErrorCode KspRhs (KSP ksp, Vec f, void *ctx)
 Run over elements in the lists. More...
 
PetscErrorCode KspMat (KSP ksp, Mat A, Mat B, void *ctx)
 Run over elenents in the list. More...
 
PetscErrorCode SnesRhs (SNES snes, Vec x, Vec f, void *ctx)
 This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver. More...
 
PetscErrorCode SnesMat (SNES snes, Vec x, Mat A, Mat B, void *ctx)
 This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver. More...
 
PetscErrorCode SNESMoFEMSetAssmblyType (SNES snes, MatAssemblyType type)
 Set assembly type at the end of SnesMat. More...
 
PetscErrorCode SNESMoFEMSetBehavior (SNES snes, MoFEMTypes bh)
 Set behavior if finite element in sequence does not exist. More...
 
PetscErrorCode f_TSSetIFunction (TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
 
PetscErrorCode f_TSSetIJacobian (TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
 
PetscErrorCode f_TSMonitorSet (TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
 

Variables

DEPRECATED typedef DataForcesAndSourcesCore DataForcesAndSurcesCore
 
DEPRECATED typedef DerivedDataForcesAndSourcesCore DerivedDataForcesAndSurcesCore
 
DEPRECATED typedef EdgeElementForcesAndSourcesCore EdgeElementForcesAndSurcesCore
 Use EdgeElementForcesAndSourcesCore. More...
 
DEPRECATED typedef FaceElementForcesAndSourcesCore FaceElementForcesAndSurcesCore
 
DEPRECATED typedef FlatPrismElementForcesAndSourcesCore FlatPrismElementForcesAndSurcesCore
 USe FlatPrismElementForcesAndSourcesCore. More...
 
DEPRECATED typedef ForcesAndSourcesCore ForcesAndSurcesCore
 
DEPRECATED typedef VertexElementForcesAndSourcesCore VertexElementForcesAndSurcesCore
 
DEPRECATED typedef VolumeElementForcesAndSourcesCore VolumeElementForcesAndSurcesCore
 
DEPRECATED typedef Problem MoFEMProblem
 
static MoABErrorCode rval
 
static PetscErrorCode ierr
 
const EntityHandle no_handle = 0
 
static const MOFEMuuid IDD_UNKNOWN_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(UNKNOWN_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_TET_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(TET_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_TRI_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(TRI_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(EDGE_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_FATPRISM_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(FATPRISM_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_FLATPRISM_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(FLATPRISM_BASE_FUNCTION_INTERFACE))
 
static const int edges_conn [] = { 0,1, 1,2, 2,0, 0,3, 1,3, 2,3 }
 
static const int oposite_edge [] = { 5, 3, 4, 1, 2, 0 }
 
static const int edge_permutations [6][6] = { { 0,1,2, 3,4,5 }, { 1,2,0, 4,5,3 }, { 2,0,1, 5,3,4 }, { 3,4,0, 2,5,1 }, { 4,5,1, 0,3,2 }, { 5,3,2, 1,4,0 } }
 
static const int edge_mirror_cross [6] = { 0,3,4,1,2,5 }
 
static const int edge_mirror_vertical [6] = { 0,4,3,2,1,5 }
 
static const int cyclic_node_rotate_face_3 [3][4] = { { 3,1,0,2 }, { 0,1,2,3 }, { 2,1,3,0 } }
 
static const int cyclic_edge_rotate_face_3 [3][6] = { { 4,0,3,5,1,2 }, { 0,1,2,3,4,5 }, { 1,4,5,2,0,3 } }
 
static const char edge_bits_mark [] = { 1, 2, 4, 8, 16, 32 }
 
static const MOFEMuuid IDD_JACOBI_BASE_FUNCTION
 
static const MOFEMuuid IDD_LEGENDRE_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(LEGENDRE_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_LOBATTO_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(LOBATTO_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_KERNEL_BASE_FUNCTION = MOFEMuuid(BitIntefaceId(KERNEL_BASE_FUNCTION_INTERFACE))
 
static const MOFEMuuid IDD_MOFEMBitLevelCoupler = MOFEMuuid( BitIntefaceId(BITLEVELCOUPLER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMBitRefManager = MOFEMuuid( BitIntefaceId(BITREFMANAGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMCoordsSystemsManager = MOFEMuuid( BitIntefaceId(COORDSSYSTEMMANAGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMCutMesh = MOFEMuuid( BitIntefaceId(CUTMESH_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMFieldBlas = MOFEMuuid( BitIntefaceId(FIELDBLAS_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMInterface = MOFEMuuid( BitIntefaceId(CORE_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMISManager = MOFEMuuid( BitIntefaceId(ISMANAGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMKspMethod = MOFEMuuid( BitIntefaceId(KSP_METHOD) )
 
static const MOFEMuuid IDD_MOFEMSnesMethod = MOFEMuuid( BitIntefaceId(SNES_METHOD) )
 
static const MOFEMuuid IDD_MOFEMTsMethod = MOFEMuuid( BitIntefaceId(TS_METHOD) )
 
static const MOFEMuuid IDD_MOFEMBasicMethod = MOFEMuuid( BitIntefaceId(BASIC_METHOD) )
 
static const MOFEMuuid IDD_MOFEMFEMethod = MOFEMuuid( BitIntefaceId(FE_METHOD) )
 
static const MOFEMuuid IDD_MOFEMEntMethod = MOFEMuuid( BitIntefaceId(ENT_METHOD) )
 
static const MOFEMuuid IDD_MOFEMMedInterface = MOFEMuuid( BitIntefaceId(MED_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMMeshRefine = MOFEMuuid( BitIntefaceId(MESH_REFINE) )
 
static const MOFEMuuid IDD_MOFEMMeshsetsManager = MOFEMuuid( BitIntefaceId(MESHSETSMANAGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMNodeMerger = MOFEMuuid( BitIntefaceId(NODEMERGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMPrismInterface = MOFEMuuid( BitIntefaceId(PRISM_INTEFACE) )
 
static const MOFEMuuid IDD_MOFEMPrismsFromSurface = MOFEMuuid(BitIntefaceId(PRISMSFROMSURFACE_INTERFACE))
 
static const MOFEMuuid IDD_MOFEMProblemsManager = MOFEMuuid( BitIntefaceId(PROBLEMSMANAGER_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMSeriesRecorder = MOFEMuuid( BitIntefaceId(SERIES_RECORDER) )
 
static const MOFEMuuid IDD_MOFEMSimple = MOFEMuuid( BitIntefaceId(SIMPLE_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMTetGegInterface = MOFEMuuid( BitIntefaceId(TETGEN_INTERFACE) )
 
static const MOFEMuuid IDD_MOFEMUnknown = MOFEMuuid( BitIntefaceId(UNKNOWNINTERFACE) )
 
static const MOFEMuuid IDD_MOFEMVEC = MOFEMuuid( BitIntefaceId(VECMANAGER_INTERFACE) )
 
const int prism_adj_edges [] = { 6,7,8, -1,-1,-1, 0,1,2 }
 
const int prism_edges_conn [6][2] = { {0,1},{1,2},{2,0}, {3,4}, {4,5}, {5,3} }
 
static moab::Error error
 
static const MOFEMuuid IDD_DMCTX = MOFEMuuid(BitIntefaceId(DMCTX_INTERFACE))
 

Detailed Description

implementation of Data Operators for Forces and Sources

name space of MoFEM library functions and classes

Create adjacent matrices using different indices.

file DataOperators.cpp

The MoFEM package is copyrighted by Lukasz Kaczmarczyk. It can be freely used for educational and research purposes by other institutions. If you use this softwre pleas cite my work.

MoFEM is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

MoFEM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with MoFEM. If not, see http://www.gnu.org/licenses/

MoFEM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with MoFEM. If not, see http://www.gnu.org/licenses/

Typedef Documentation

◆ ApproximationOrder

Approximation on the entity.

Definition at line 33 of file Common.hpp.

◆ BasicMethodsSequence

Definition at line 54 of file AuxPETSc.hpp.

◆ BitFEId

typedef std::bitset<BITFEID_SIZE> MoFEM::BitFEId

Definition at line 51 of file Common.hpp.

◆ BitFieldId

typedef std::bitset<BITFIELDID_SIZE> MoFEM::BitFieldId

Definition at line 50 of file Common.hpp.

◆ BitIntefaceId

Definition at line 53 of file Common.hpp.

◆ BitProblemId

typedef std::bitset<BITPROBLEMID_SIZE> MoFEM::BitProblemId

Definition at line 52 of file Common.hpp.

◆ BitRefEdges

typedef std::bitset<BITREFEDGES_SIZE> MoFEM::BitRefEdges

Definition at line 43 of file Common.hpp.

◆ BitRefLevel

typedef std::bitset<BITREFLEVEL_SIZE> MoFEM::BitRefLevel

Bit structure attached to each entity identifying to what mesh entity is attached.

Examples:
MagneticElement.hpp, Remodeling.hpp, unsaturated_transport.cpp, and UnsaturatedFlow.hpp.

Definition at line 48 of file Common.hpp.

◆ CoordSys_multiIndex

typedef multi_index_container< boost::shared_ptr<CoordSys>, indexed_by< ordered_unique< tag<Meshset_mi_tag>, member<CoordSys,EntityHandle,&CoordSys::meshSet> >, ordered_unique< tag<CoordSysName_mi_tag>, const_mem_fun<CoordSys,boost::string_ref,&CoordSys::getNameRef> > > > MoFEM::CoordSys_multiIndex

Definition at line 140 of file CoordSysMultiIndices.hpp.

◆ CubitBCType

bc & material meshsets

Definition at line 125 of file Common.hpp.

◆ CubitMeshSet_multiIndex

Stores data about meshsets (see CubitMeshSets) storing data about boundary conditions, interfaces, sidesets, nodests, blocksets.

Parameters
Meshset_mi_tagindex by meshset handle
CubitMeshSets_mi_tagindex by bc type, see CubitBC
CubitMeshSets_mask_meshset_mi_tagindex by NODESET, SIDESET, BLOCKSET only
CubitMeshSets_nameindex by meshset name
Composite_Cubit_msId_And_MeshSetType_mi_tagindex by meshset id and type NODESET, SIDESET or BLOCKSET

Example:

MeshsetsManager *m_mng;
ierr = m_field.query_interface(m_mng); CHKERRQ(ierr);
CubitMeshSet_multiIndex &meshsets_index = m_mng->etMeshsetsMultindex();
CubitMeshSet_multiIndex::index<CubitMeshSets_mask_meshset_mi_tag>::type::iterator mit,hi_mit;
mit = meshsets_index.get<CubitMeshSets_mask_meshset_mi_tag>().lower_bound(BLOCKSET);
hi_mit = meshsets_index.get<CubitMeshSets_mask_meshset_mi_tag>().upper_bound(BLOCKSET);
// Make a loop over all BLOCKSET
for(;mit!=hi_mit;mit++) {
int id = mit->getMeshsetId(); // get blockset id
EntityHandle handle = mit->getMeshset(); // get block meshset
std::vector< double > attributes;
// get block attributes
ierr = mit->getAttributes(attributes); CHKERRQ(ierr);
// do something
}

Definition at line 406 of file BCMultiIndices.hpp.

◆ CubitMeshsetById

typedef CubitMeshSet_multiIndex::index<CubitMeshSets_mi_tag>::type MoFEM::CubitMeshsetById

Definition at line 30 of file MeshsetsManager.hpp.

◆ CubitMeshsetByMask

typedef CubitMeshSet_multiIndex::index<CubitMeshSets_mask_meshset_mi_tag>::type MoFEM::CubitMeshsetByMask

Definition at line 26 of file MeshsetsManager.hpp.

◆ CubitMeshsetByName

typedef CubitMeshSet_multiIndex::index<CubitMeshSets_name>::type MoFEM::CubitMeshsetByName

Definition at line 28 of file MeshsetsManager.hpp.

◆ CubitMeshsetByType

typedef CubitMeshSet_multiIndex::index<CubitMeshSets_mi_tag>::type MoFEM::CubitMeshsetByType

Definition at line 24 of file MeshsetsManager.hpp.

◆ DofIdx

typedef int MoFEM::DofIdx

Index of DOF.

Definition at line 28 of file Common.hpp.

◆ DofsAllacator

typedef ublas::unbounded_array< boost::shared_ptr<const FEDofEntity>, std::allocator<boost::shared_ptr<const FEDofEntity> >> MoFEM::DofsAllacator

Definition at line 32 of file DataStructures.hpp.

◆ DoubleAllacator

typedef std::vector<double,std::allocator<double> > MoFEM::DoubleAllacator

Definition at line 131 of file Common.hpp.

◆ EntIdx

typedef int MoFEM::EntIdx

Index of DOF on the entity.

Definition at line 30 of file Common.hpp.

◆ EntPart

Partition owning entity.

Definition at line 31 of file Common.hpp.

◆ FEIdx

typedef int MoFEM::FEIdx

Index of the element.

Definition at line 29 of file Common.hpp.

◆ FEMethodsSequence

Definition at line 53 of file AuxPETSc.hpp.

◆ Field_multiIndex_view

typedef multi_index_container< boost::shared_ptr<Field>, indexed_by< ordered_unique< tag<BitFieldId_mi_tag>, const_mem_fun<Field,const BitFieldId&,&Field::getId>, LtBit<BitFieldId> >> > MoFEM::Field_multiIndex_view

Definition at line 451 of file FieldMultiIndices.hpp.

◆ FieldCoefficientsNumber

Number of field coefficients.

Definition at line 34 of file Common.hpp.

◆ FieldData

Field data type.

Definition at line 32 of file Common.hpp.

◆ FieldEntity_multiIndex_ent_view

typedef multi_index_container< boost::shared_ptr<FieldEntity>, indexed_by< sequenced<>, ordered_non_unique< tag<Ent_mi_tag>, const_mem_fun<FieldEntity,EntityHandle,&FieldEntity::getEnt> > >> MoFEM::FieldEntity_multiIndex_ent_view

Definition at line 919 of file EntsMultiIndices.hpp.

◆ IntAllacator

typedef std::vector<int,std::allocator<int> > MoFEM::IntAllacator

Definition at line 130 of file Common.hpp.

◆ MatrixAdaptor

typedef ublas::matrix<double,ublas::row_major,ublas::shallow_array_adaptor<double> > MoFEM::MatrixAdaptor
Examples:
SmallStrainPlasticity.hpp.

Definition at line 143 of file Common.hpp.

◆ MatrixDouble

◆ MatrixDouble3by3

typedef ublas::matrix<double,ublas::row_major,ublas::bounded_array<double,9> > MoFEM::MatrixDouble3by3

Definition at line 137 of file Common.hpp.

◆ MoABErrorCode

typedef ErrorCode MoFEM::MoABErrorCode

Definition at line 23 of file Common.hpp.

◆ NumeredDofEntity_multiIndex_coeff_idx_ordered_non_unique

Definition at line 772 of file DofsMultiIndices.hpp.

◆ NumeredDofEntity_multiIndex_global_index_view

typedef multi_index_container< boost::shared_ptr<NumeredDofEntity>, indexed_by< ordered_unique< member<NumeredDofEntity,const DofIdx,&NumeredDofEntity::petscGloablDofIdx> > > > MoFEM::NumeredDofEntity_multiIndex_global_index_view

Definition at line 876 of file DofsMultiIndices.hpp.

◆ NumeredDofEntity_multiIndex_idx_view_hashed

typedef multi_index_container< boost::shared_ptr<NumeredDofEntity>, indexed_by< hashed_unique< const_mem_fun<NumeredDofEntity,DofIdx,&NumeredDofEntity::getDofIdx> > >> MoFEM::NumeredDofEntity_multiIndex_idx_view_hashed

Definition at line 746 of file DofsMultiIndices.hpp.

◆ NumeredDofEntity_multiIndex_petsc_global_dof_view_ordered_non_unique

typedef multi_index_container< boost::shared_ptr<NumeredDofEntity>, indexed_by< ordered_non_unique< const_mem_fun<NumeredDofEntity,DofIdx,&NumeredDofEntity::getPetscGlobalDofIdx> > >> MoFEM::NumeredDofEntity_multiIndex_petsc_global_dof_view_ordered_non_unique

Definition at line 763 of file DofsMultiIndices.hpp.

◆ NumeredDofEntity_multiIndex_petsc_local_dof_view_ordered_non_unique

typedef multi_index_container< boost::shared_ptr<NumeredDofEntity>, indexed_by< ordered_non_unique< const_mem_fun<NumeredDofEntity,DofIdx,&NumeredDofEntity::getPetscLocalDofIdx> > >> MoFEM::NumeredDofEntity_multiIndex_petsc_local_dof_view_ordered_non_unique

Definition at line 754 of file DofsMultiIndices.hpp.

◆ NumeredDofEntity_multiIndex_uid_view_ordered

typedef multi_index_container< boost::shared_ptr<NumeredDofEntity>, indexed_by< ordered_unique< const_mem_fun< NumeredDofEntity::interface_type_DofEntity, const UId, &NumeredDofEntity::getGlobalUniqueId > > >> MoFEM::NumeredDofEntity_multiIndex_uid_view_ordered

Definition at line 737 of file DofsMultiIndices.hpp.

◆ ptrWrapperRefElement

Definition at line 187 of file FEMultiIndices.hpp.

◆ ShapeFunctionBasesVector

typedef std::vector<boost::shared_ptr<MatrixDouble> > MoFEM::ShapeFunctionBasesVector

Definition at line 146 of file Common.hpp.

◆ ShortId

Unique Id in the field.

Definition at line 39 of file Common.hpp.

◆ UId

typedef uint128_t MoFEM::UId

Unique Id.

Definition at line 38 of file Common.hpp.

◆ VectorAdaptor

typedef ublas::vector<double,ublas::shallow_array_adaptor<double> > MoFEM::VectorAdaptor
Examples:
SmallStrainPlasticity.hpp.

Definition at line 142 of file Common.hpp.

◆ VectorDofs

typedef ublas::vector< boost::shared_ptr<const FEDofEntity>,DofsAllacator> MoFEM::VectorDofs

Definition at line 35 of file DataStructures.hpp.

◆ VectorDouble

typedef ublas::vector<double,DoubleAllacator > MoFEM::VectorDouble

◆ VectorDouble3

typedef ublas::vector<double,ublas::bounded_array<double,3> > MoFEM::VectorDouble3

Definition at line 138 of file Common.hpp.

◆ VectorDouble9

typedef ublas::vector<double,ublas::bounded_array<double,9> > MoFEM::VectorDouble9

Definition at line 139 of file Common.hpp.

◆ VectorInt

typedef ublas::vector<int,IntAllacator > MoFEM::VectorInt

Definition at line 134 of file Common.hpp.

◆ VectorIntAdaptor

typedef ublas::vector<int,ublas::shallow_array_adaptor<int> > MoFEM::VectorIntAdaptor

Definition at line 144 of file Common.hpp.

Function Documentation

◆ cOnstructor()

template<class T >
void MoFEM::cOnstructor ( DataForcesAndSourcesCore data,
EntityType  type,
 
)

Definition at line 142 of file DataStructures.cpp.

142  {
143 
144  data->dataOnEntities[MBENTITYSET].push_back(new T());
145 
146  switch (type) {
147  case MBTET:
148  data->dataOnEntities[MBVERTEX].push_back(new T());
149  for(int ee = 0;ee<6;ee++) {
150  data->dataOnEntities[MBEDGE].push_back(new T());
151  }
152  for(int ff = 0;ff<4;ff++) {
153  data->dataOnEntities[MBTRI].push_back(new T());
154  }
155  data->dataOnEntities[MBTET].push_back(new T());
156  break;
157  case MBTRI:
158  data->dataOnEntities[MBVERTEX].push_back(new T());
159  for(int ee = 0;ee<3;ee++) {
160  data->dataOnEntities[MBEDGE].push_back(new T());
161  }
162  data->dataOnEntities[MBTRI].push_back(new T());
163  break;
164  case MBEDGE:
165  data->dataOnEntities[MBVERTEX].push_back(new T());
166  data->dataOnEntities[MBEDGE].push_back(new T());
167  break;
168  case MBVERTEX:
169  data->dataOnEntities[MBVERTEX].push_back(new T());
170  break;
171  case MBPRISM:
172  data->dataOnEntities[MBVERTEX].push_back(new T());
173  for(int ee = 0;ee<9;ee++) {
174  data->dataOnEntities[MBEDGE].push_back(new T());
175  }
176  for(int ff = 0;ff<5;ff++) {
177  data->dataOnEntities[MBQUAD].push_back(new T());
178  }
179  for(int ff = 0;ff<5;ff++) {
180  data->dataOnEntities[MBTRI].push_back(new T());
181  }
182  data->dataOnEntities[MBPRISM].push_back(new T());
183  break;
184  default:
186  }
187 
188 }
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ determinant()

template<class T >
static double MoFEM::determinant ( T &  t)
static

Definition at line 56 of file NodeMerger.cpp.

56  {
57  return
58  +t(0,0)*t(1,1)*t(2,2) + t(1,0)*t(2,1)*t(0,2)
59  +t(2,0)*t(0,1)*t(1,2) - t(0,0)*t(2,1)*t(1,2)
60  -t(2,0)*t(1,1)*t(0,2) - t(1,0)*t(0,1)*t(2,2);
61 }

◆ error_printf_hilight()

static void MoFEM::error_printf_hilight ( void  )
static

Definition at line 267 of file Core.cpp.

267  {
268 #if defined(PETSC_HAVE_UNISTD_H) && defined(PETSC_USE_ISATTY)
269  if (PetscErrorPrintf == PetscErrorPrintfDefault) {
270  if (isatty(fileno(PETSC_STDERR))) fprintf(PETSC_STDERR,"\033[1;32m");
271  }
272 #endif
273 }

◆ error_printf_normal()

static void MoFEM::error_printf_normal ( void  )
static

Definition at line 275 of file Core.cpp.

275  {
276 #if defined(PETSC_HAVE_UNISTD_H) && defined(PETSC_USE_ISATTY)
277  if (PetscErrorPrintf == PetscErrorPrintfDefault) {
278  if (isatty(fileno(PETSC_STDERR))) fprintf(PETSC_STDERR,"\033[0;39m\033[0;49m");
279  }
280 #endif
281 }

◆ f_TSMonitorSet()

PetscErrorCode MoFEM::f_TSMonitorSet ( TS  ts,
PetscInt  step,
PetscReal  t,
Vec  u,
void *  ctx 
)
Examples:
UnsaturatedFlow.hpp.

Definition at line 187 of file TsCtx.cpp.

187  {
188  // PetscValidHeaderSpecific(ts,TS_CLASSID,1);
189  PetscFunctionBegin;
190  TsCtx* ts_ctx = (TsCtx*)ctx;
191  PetscLogEventBegin(ts_ctx->USER_EVENT_TsCtxRHSFunction,0,0,0,0);
192  ierr = VecGhostUpdateBegin(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
193  ierr = VecGhostUpdateEnd(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
194  ierr = ts_ctx->mField.query_interface<VecManager>()->setLocalGhostVector(
195  ts_ctx->problemName,COL,u,INSERT_VALUES,SCATTER_REVERSE
196  ); CHKERRQ(ierr);
197  //preproces
198  TsCtx::BasicMethodsSequence::iterator bit = ts_ctx->preProcess_Monitor.begin();
199  for(;bit!=ts_ctx->preProcess_Monitor.end();bit++) {
200  (*bit)->ts_u = u;
201  (*bit)->ts_t = t;
202  (*bit)->ts_step = step;
203  (*bit)->ts_F = PETSC_NULL;
204  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
205  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
206  ierr = ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
207  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE); CHKERRQ(ierr); CHKERRQ(ierr);
208  }
209  TsCtx::FEMethodsSequence::iterator lit = ts_ctx->loops_to_do_Monitor.begin();
210  for(;lit!=ts_ctx->loops_to_do_Monitor.end();lit++) {
211  lit->second->ts_u = u;
212  lit->second->ts_t = t;
213  lit->second->ts_step = step;
214  lit->second->ts_F = PETSC_NULL;
215  ierr = lit->second->setTsCtx(TSMethod::CTX_TSTSMONITORSET);
216  ierr = lit->second->setTs(ts); CHKERRQ(ierr);
217  ierr = ts_ctx->mField.loop_finite_elements(ts_ctx->problemName,lit->first,*(lit->second),ts_ctx->bH); CHKERRQ(ierr);
218  ierr = lit->second->setTsCtx(TSMethod::CTX_TSNONE);
219  }
220  //post process
221  bit = ts_ctx->postProcess_Monitor.begin();
222  for(;bit!=ts_ctx->postProcess_Monitor.end();bit++) {
223  (*bit)->ts_u = u;
224  (*bit)->ts_t = t;
225  (*bit)->ts_step = step;
226  (*bit)->ts_F = PETSC_NULL;
227  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
228  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
229  ierr = ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
230  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE); CHKERRQ(ierr);
231  }
232  PetscLogEventEnd(ts_ctx->USER_EVENT_TsCtxRHSFunction,0,0,0,0);
233  PetscFunctionReturn(0);
234 }
CHKERRQ(ierr)

◆ f_TSSetIFunction()

PetscErrorCode MoFEM::f_TSSetIFunction ( TS  ts,
PetscReal  t,
Vec  u,
Vec  u_t,
Vec  F,
void *  ctx 
)

Definition at line 56 of file TsCtx.cpp.

56  {
57  // PetscValidHeaderSpecific(ts,TS_CLASSID,1);
58  PetscFunctionBegin;
59  TsCtx* ts_ctx = (TsCtx*)ctx;
60  PetscLogEventBegin(ts_ctx->USER_EVENT_TsCtxIFunction,0,0,0,0);
61  ierr = VecGhostUpdateBegin(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
62  ierr = VecGhostUpdateEnd(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
63  ierr = VecGhostUpdateBegin(u_t,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
64  ierr = VecGhostUpdateEnd(u_t,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
65  ierr = ts_ctx->mField.query_interface<VecManager>()->setLocalGhostVector(
66  ts_ctx->problemName,COL,u,INSERT_VALUES,SCATTER_REVERSE
67  ); CHKERRQ(ierr);
68  ierr = VecZeroEntries(F); CHKERRQ(ierr);
69  ierr = VecGhostUpdateBegin(F,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
70  ierr = VecGhostUpdateEnd(F,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
71  int step;
72  ierr = TSGetTimeStepNumber(ts,&step); CHKERRQ(ierr);
73  //preprocess
74  TsCtx::BasicMethodsSequence::iterator bit = ts_ctx->preProcess_IFunction.begin();
75  for(;bit!=ts_ctx->preProcess_IFunction.end();bit++) {
76  (*bit)->ts_u = u;
77  (*bit)->ts_u_t = u_t;
78  (*bit)->ts_F = F;
79  (*bit)->ts_t = t;
80  (*bit)->ts_step = step;
81  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
82  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
83  ierr = ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
84  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE);
85  }
86  //fe loops
87  TsCtx::FEMethodsSequence::iterator lit = ts_ctx->loops_to_do_IFunction.begin();
88  for(;lit!=ts_ctx->loops_to_do_IFunction.end();lit++) {
89  lit->second->ts_u = u;
90  lit->second->ts_u_t = u_t;
91  lit->second->ts_F = F;
92  lit->second->ts_t = t;
93  lit->second->ts_step = step;
94  ierr = lit->second->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
95  ierr = lit->second->setTs(ts); CHKERRQ(ierr);
96  ierr = ts_ctx->mField.loop_finite_elements(ts_ctx->problemName,lit->first,*(lit->second),ts_ctx->bH); CHKERRQ(ierr);
97  ierr = lit->second->setTsCtx(TSMethod::CTX_TSNONE);
98  }
99  //post process
100  bit = ts_ctx->postProcess_IFunction.begin();
101  for(;bit!=ts_ctx->postProcess_IFunction.end();bit++) {
102  (*bit)->ts_u = u;
103  (*bit)->ts_u_t = u_t;
104  (*bit)->ts_F = F;
105  (*bit)->ts_t = t;
106  (*bit)->ts_step = step;
107  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIFUNCTION);
108  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
109  ierr = ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
110  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE);
111  }
112  ierr = VecGhostUpdateBegin(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
113  ierr = VecGhostUpdateEnd(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
114  ierr = VecAssemblyBegin(F); CHKERRQ(ierr);
115  ierr = VecAssemblyEnd(F); CHKERRQ(ierr);
116  PetscLogEventEnd(ts_ctx->USER_EVENT_TsCtxIFunction,0,0,0,0);
117  PetscFunctionReturn(0);
118 }
CHKERRQ(ierr)

◆ f_TSSetIJacobian()

PetscErrorCode MoFEM::f_TSSetIJacobian ( TS  ts,
PetscReal  t,
Vec  u,
Vec  u_t,
PetscReal  a,
Mat  A,
Mat  B,
void *  ctx 
)

Definition at line 119 of file TsCtx.cpp.

119  {
120  // PetscValidHeaderSpecific(ts,TS_CLASSID,1);
121  PetscFunctionBegin;
122  TsCtx* ts_ctx = (TsCtx*)ctx;
123  PetscLogEventBegin(ts_ctx->USER_EVENT_TsCtxIFunction,0,0,0,0);
124  ierr = VecGhostUpdateBegin(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
125  ierr = VecGhostUpdateEnd(u,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
126  ierr = VecGhostUpdateBegin(u_t,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
127  ierr = VecGhostUpdateEnd(u_t,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
128  ierr = ts_ctx->mField.query_interface<VecManager>()->setLocalGhostVector(
129  ts_ctx->problemName,COL,u,INSERT_VALUES,SCATTER_REVERSE
130  ); CHKERRQ(ierr);
131  if(ts_ctx->zeroMatrix) {
132  ierr = MatZeroEntries(B); CHKERRQ(ierr);
133  }
134  int step;
135  ierr = TSGetTimeStepNumber(ts,&step); CHKERRQ(ierr);
136  //preproces
137  TsCtx::BasicMethodsSequence::iterator bit = ts_ctx->preProcess_IJacobian.begin();
138  for(;bit!=ts_ctx->preProcess_IJacobian.end();bit++) {
139  (*bit)->ts_u = u;
140  (*bit)->ts_u_t = u_t;
141  (*bit)->ts_A = A;
142  (*bit)->ts_B = B;
143  (*bit)->ts_t = t;
144  (*bit)->ts_a = a;
145  (*bit)->ts_step = step;
146  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
147  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
148  ierr = ts_ctx->mField.problem_basic_method_preProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
149  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE); CHKERRQ(ierr); CHKERRQ(ierr);
150  }
151  TsCtx::FEMethodsSequence::iterator lit = ts_ctx->loops_to_do_IJacobian.begin();
152  for(;lit!=ts_ctx->loops_to_do_IJacobian.end();lit++) {
153  lit->second->ts_u = u;
154  lit->second->ts_u_t = u_t;
155  lit->second->ts_A = A;
156  lit->second->ts_B = B;
157  lit->second->ts_t = t;
158  lit->second->ts_a = a;
159  lit->second->ts_step = step;
160  ierr = lit->second->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
161  ierr = lit->second->setTs(ts); CHKERRQ(ierr);
162  ierr = ts_ctx->mField.loop_finite_elements(ts_ctx->problemName,lit->first,*(lit->second),ts_ctx->bH); CHKERRQ(ierr);
163  ierr = lit->second->setTsCtx(TSMethod::CTX_TSNONE); CHKERRQ(ierr);
164  }
165  //post process
166  bit = ts_ctx->postProcess_IJacobian.begin();
167  for(;bit!=ts_ctx->postProcess_IJacobian.end();bit++) {
168  (*bit)->ts_u = u;
169  (*bit)->ts_u_t = u_t;
170  (*bit)->ts_A = A;
171  (*bit)->ts_B = B;
172  (*bit)->ts_t = t;
173  (*bit)->ts_a = a;
174  (*bit)->ts_step = step;
175  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSSETIJACOBIAN);
176  ierr = (*bit)->setTs(ts); CHKERRQ(ierr);
177  ierr = ts_ctx->mField.problem_basic_method_postProcess(ts_ctx->problemName,*(*(bit))); CHKERRQ(ierr);
178  ierr = (*bit)->setTsCtx(TSMethod::CTX_TSNONE); CHKERRQ(ierr);
179  }
180  if(ts_ctx->zeroMatrix) {
181  ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
182  ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
183  }
184  PetscLogEventEnd(ts_ctx->USER_EVENT_TsCtxIFunction,0,0,0,0);
185  PetscFunctionReturn(0);
186 }
CHKERRQ(ierr)

◆ fNBEDGE_H1()

int MoFEM::fNBEDGE_H1 ( int  P)

number of approx. functions for H1 space on edge

Definition at line 39 of file CoreDataStructures.hpp.

39 { return NBEDGE_H1(P); }
#define NBEDGE_H1(P)
Number of dofs on edge for H1 space.

◆ fNBEDGE_HCURL()

int MoFEM::fNBEDGE_HCURL ( int  P)

Definition at line 50 of file CoreDataStructures.hpp.

50 { return NBEDGE_HCURL(P); }
#define NBEDGE_HCURL(P)
Number of base functions H curl on faces.

◆ fNBEDGE_HDIV()

int MoFEM::fNBEDGE_HDIV ( int  P)

number of approx. functions for HDIV space on edge

Definition at line 59 of file CoreDataStructures.hpp.

59 { (void)P; return NBEDGE_HDIV(P); }
#define NBEDGE_HDIV(P)

◆ fNBEDGE_L2()

int MoFEM::fNBEDGE_L2 ( int  P)

Definition at line 33 of file CoreDataStructures.hpp.

33 { return NBEDGE_L2(P); }
#define NBEDGE_L2(P)

◆ fNBENTITY_GENERIC()

int MoFEM::fNBENTITY_GENERIC ( int  P)

Definition at line 27 of file CoreDataStructures.hpp.

27 { (void)P; return 0; }

◆ fNBENTITYSET_NOFIELD()

int MoFEM::fNBENTITYSET_NOFIELD ( int  P)

Definition at line 28 of file CoreDataStructures.hpp.

28 { (void)P; return 1; }

◆ fNBFACEQUAD_H1()

int MoFEM::fNBFACEQUAD_H1 ( int  P)

Definition at line 42 of file CoreDataStructures.hpp.

42 { return NBFACEQUAD_H1(P); }
#define NBFACEQUAD_H1(P)

◆ fNBFACETRI_AINSWORTH_HDIV()

int MoFEM::fNBFACETRI_AINSWORTH_HDIV ( int  P)

number of approx. functions for HDIV space on face

Definition at line 61 of file CoreDataStructures.hpp.

61 { return NBFACETRI_AINSWORTH_HDIV(P); }
#define NBFACETRI_AINSWORTH_HDIV(P)

◆ fNBFACETRI_DEMKOWICZ_HDIV()

int MoFEM::fNBFACETRI_DEMKOWICZ_HDIV ( int  P)

Definition at line 65 of file CoreDataStructures.hpp.

65 { return NBFACETRI_DEMKOWICZ_HDIV(P); }
#define NBFACETRI_DEMKOWICZ_HDIV(P)

◆ fNBFACETRI_H1()

int MoFEM::fNBFACETRI_H1 ( int  P)

number of approx. functions for H1 space on face

Definition at line 41 of file CoreDataStructures.hpp.

41 { return NBFACETRI_H1(P); }
#define NBFACETRI_H1(P)
Number of dofs on face for H1 space.

◆ fNBFACETRI_HCURL()

int MoFEM::fNBFACETRI_HCURL ( int  P)

Definition at line 51 of file CoreDataStructures.hpp.

51 { return NBFACETRI_HCURL(P); }
#define NBFACETRI_HCURL(P)

◆ fNBFACETRI_L2()

int MoFEM::fNBFACETRI_L2 ( int  P)

Definition at line 32 of file CoreDataStructures.hpp.

32 { return NBFACETRI_L2(P); }
#define NBFACETRI_L2(P)

◆ fNBVERTEX_H1()

int MoFEM::fNBVERTEX_H1 ( int  P)

number of approx. functions for H1 space on vertex

Definition at line 37 of file CoreDataStructures.hpp.

37 { return (P==1) ? 1 : 0; }

◆ fNBVERTEX_HCURL()

int MoFEM::fNBVERTEX_HCURL ( int  P)

number of approx. functions for HCURL space on vertex

Definition at line 49 of file CoreDataStructures.hpp.

49 { (void)P; return 0; }

◆ fNBVERTEX_HDIV()

int MoFEM::fNBVERTEX_HDIV ( int  P)

number of approx. functions for HDIV space on vertex

zero number of digrees of freedom on vertex for that space

Definition at line 57 of file CoreDataStructures.hpp.

57 { (void)P; return 0; }

◆ fNBVERTEX_L2()

int MoFEM::fNBVERTEX_L2 ( int  P)

Definition at line 30 of file CoreDataStructures.hpp.

30 { (void)P; return 0; }

◆ fNBVOLUMEPRISM_H1()

int MoFEM::fNBVOLUMEPRISM_H1 ( int  P)

Definition at line 45 of file CoreDataStructures.hpp.

45 { return NBVOLUMEPRISM_H1(P); }
#define NBVOLUMEPRISM_H1(P)

◆ fNBVOLUMETET_AINSWORTH_HDIV()

int MoFEM::fNBVOLUMETET_AINSWORTH_HDIV ( int  P)

number of approx. functions for HDIV space on volume

Definition at line 63 of file CoreDataStructures.hpp.

63 { return NBVOLUMETET_AINSWORTH_HDIV(P); }
#define NBVOLUMETET_AINSWORTH_HDIV(P)

◆ fNBVOLUMETET_DEMKOWICZ_HDIV()

int MoFEM::fNBVOLUMETET_DEMKOWICZ_HDIV ( int  P)

number of approx. functions for HDIV space on volume

Definition at line 67 of file CoreDataStructures.hpp.

67 { return NBVOLUMETET_DEMKOWICZ_HDIV(P); }
#define NBVOLUMETET_DEMKOWICZ_HDIV(P)

◆ fNBVOLUMETET_H1()

int MoFEM::fNBVOLUMETET_H1 ( int  P)

number of approx. functions for H1 space on volume

Definition at line 44 of file CoreDataStructures.hpp.

44 { return NBVOLUMETET_H1(P); }
#define NBVOLUMETET_H1(P)
Number of dofs on volume for H1 space.

◆ fNBVOLUMETET_HCURL()

int MoFEM::fNBVOLUMETET_HCURL ( int  P)

Definition at line 52 of file CoreDataStructures.hpp.

52 { return NBVOLUMETET_HCURL(P); }
#define NBVOLUMETET_HCURL(P)

◆ fNBVOLUMETET_L2()

int MoFEM::fNBVOLUMETET_L2 ( int  P)

Definition at line 31 of file CoreDataStructures.hpp.

31 { return NBVOLUMETET_L2(P); }
#define NBVOLUMETET_L2(P)
Number of dofs for L2 space.

◆ get_fe_dof_view()

template<typename FE_DOFS , typename MOFEM_DOFS , typename MOFEM_DOFS_VIEW >
static PetscErrorCode MoFEM::get_fe_dof_view ( const FE_DOFS &  fe_dofs_view,
const MOFEM_DOFS &  mofem_dofs,
MOFEM_DOFS_VIEW &  mofem_dofs_view,
const int  operation_type 
)
static

Definition at line 940 of file FEMultiIndices.cpp.

945  {
946  PetscFunctionBegin;
947  typename boost::multi_index::index<MOFEM_DOFS,Unique_mi_tag>::type::iterator mofem_it,mofem_it_end;
948  typename FE_DOFS::iterator it,it_end;
949  if(operation_type==moab::Interface::UNION) {
950  mofem_it = mofem_dofs.template get<Unique_mi_tag>().begin();
951  mofem_it_end = mofem_dofs.template get<Unique_mi_tag>().end();
952  it = fe_dofs_view.begin();
953  it_end = fe_dofs_view.end();
954  for(;it!=it_end;it++) {
955  const UId &globalUid = (*it)->getGlobalUniqueId();
956  if(mofem_it != mofem_it_end) {
957  if((*mofem_it)->getGlobalUniqueId() != globalUid) {
958  mofem_it = mofem_dofs.template get<Unique_mi_tag>().find(globalUid);
959  } // else lucky guess
960  } else {
961  mofem_it = mofem_dofs.template get<Unique_mi_tag>().find(globalUid);
962  }
963  if(mofem_it != mofem_it_end) {
964  mofem_dofs_view.insert(mofem_dofs_view.end(),*mofem_it);
965  mofem_it++;
966  }
967  }
968  } else {
969  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"not implemented");
970  }
971  PetscFunctionReturn(0);
972 }
uint128_t UId
Unique Id.
Definition: Common.hpp:38

◆ get_porblem_col_indices()

static PetscErrorCode MoFEM::get_porblem_col_indices ( const ForcesAndSourcesCore fe_ptr,
const EntityType  type,
const int  side,
const std::string  field_name,
VectorInt indices 
)
static

Definition at line 1174 of file ForcesAndSourcesCore.cpp.

1175  {
1176  PetscFunctionBegin;;
1177 
1178 
1179 
1180  switch(type) {
1181  case MBVERTEX:
1182  ierr = fe_ptr->getProblemNodesColIndices(field_name,indices); CHKERRQ(ierr);
1183  break;
1184  default:
1185  ierr = fe_ptr->getProblemTypeColIndices(field_name,type,side,indices); CHKERRQ(ierr);
1186  }
1187 
1188  PetscFunctionReturn(0);
1189 }
CHKERRQ(ierr)

◆ get_porblem_row_indices()

static PetscErrorCode MoFEM::get_porblem_row_indices ( const ForcesAndSourcesCore fe_ptr,
const EntityType  type,
const int  side,
const std::string  field_name,
VectorInt indices 
)
static

Definition at line 1157 of file ForcesAndSourcesCore.cpp.

1158  {
1159  PetscFunctionBegin;
1160 
1161 
1162 
1163  switch(type) {
1164  case MBVERTEX:
1165  ierr = fe_ptr->getProblemNodesRowIndices(field_name,indices); CHKERRQ(ierr);
1166  break;
1167  default:
1168  ierr = fe_ptr->getProblemTypeRowIndices(field_name,type,side,indices); CHKERRQ(ierr);
1169  }
1170 
1171  PetscFunctionReturn(0);
1172 }
CHKERRQ(ierr)

◆ get_tag_ptr()

void* MoFEM::get_tag_ptr ( SequenceManager *  sequence_manager,
Tag  th,
EntityHandle  ent,
int tag_size 
)

Definition at line 44 of file EntsMultiIndices.cpp.

44  {
45  ApproximationOrder *ret_val;
46  if(th->get_storage_type()==MB_TAG_SPARSE) {
47  rval = static_cast<SparseTag*>(th)->get_data(
48  sequence_manager,&error,&ent,1,(const void**)&ret_val,tag_size
49  ); MOAB_THROW(rval);
50  return ret_val;
51  } else {
52  rval = static_cast<DenseTag*>(th)->get_data(
53  sequence_manager,&error,&ent,1,(const void**)&ret_val,tag_size
54  ); MOAB_THROW(rval);
55  return ret_val;
56  }
57 }
static MoABErrorCode rval
Definition: Common.hpp:25
#define MOAB_THROW(a)
Definition: definitions.h:376
static moab::Error error
int ApproximationOrder
Approximation on the entity.
Definition: Common.hpp:33

◆ getMaxOrder()

template<typename DOFMULTIINDEX >
static int MoFEM::getMaxOrder ( const DOFMULTIINDEX &  dof_multi_index)
static

Definition at line 150 of file ForcesAndSourcesCore.cpp.

152  {
153  int max_order = 0;
154  typename DOFMULTIINDEX::iterator dit,hi_dit;
155  dit = dof_multi_index.begin();
156  hi_dit = dof_multi_index.end();
157  for(;dit!=hi_dit;dit++) {
158  if((*dit)->getEntDofIdx()!=0) continue;
159  int dit_max_order = (*dit)->getMaxOrder();
160  max_order = (max_order>dit_max_order) ? max_order : dit_max_order;
161  }
162  return max_order;
163 }

◆ getPatentEnt()

PetscErrorCode MoFEM::getPatentEnt ( Interface moab,
Range  ents,
std::vector< EntityHandle vec_patent_ent 
)

Definition at line 123 of file EntsMultiIndices.cpp.

123  {
124 
125  PetscFunctionBegin;
126  Tag th_ref_parent_handle;
127  rval = moab.tag_get_handle("_RefParentHandle",th_ref_parent_handle); CHKERRQ_MOAB(rval);
128  vec_patent_ent.resize(ents.size());
129  rval = moab.tag_get_data(th_ref_parent_handle,ents,&*vec_patent_ent.begin()); CHKERRQ_MOAB(rval);
130  PetscFunctionReturn(0);
131 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:362
static MoABErrorCode rval
Definition: Common.hpp:25

◆ getTensor0FormData< double, DoubleAllacator >()

Definition at line 66 of file DataStructures.cpp.

68  {
69  return FTensor::Tensor0<double*>(&*data.data().begin());
70 }
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ getTensor1FormData< 2, double, ublas::row_major, DoubleAllacator >()

Definition at line 83 of file DataStructures.cpp.

85  {
86  if(data.size1()!=2) {
87  THROW_MESSAGE("Wrong size of data matrix");
88  }
89  return FTensor::Tensor1<double*,2>(&data(0,0),&data(1,0));
90 }
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:389
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ getTensor1FormData< 3, double, ublas::row_major, DoubleAllacator >()

Definition at line 73 of file DataStructures.cpp.

75  {
76  if(data.size1()!=3) {
77  THROW_MESSAGE("Wrong size of data matrix");
78  }
79  return FTensor::Tensor1<double*,3>(&data(0,0),&data(1,0),&data(2,0));
80 }
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:389
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ getTensor2FormData< 3, 2, double, ublas::row_major, DoubleAllacator >()

Definition at line 105 of file DataStructures.cpp.

107  {
108  if(data.size1()!=6) {
109  THROW_MESSAGE("Wrong size of data matrix");
110  }
112  &data(0,0),&data(1,0),&data(2,0),&data(3,0),&data(4,0),&data(5,0)
113  );
114 }
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:389
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ getTensor2FormData< 3, 3, double, ublas::row_major, DoubleAllacator >()

Definition at line 93 of file DataStructures.cpp.

95  {
96  if(data.size1()!=9) {
97  THROW_MESSAGE("Wrong size of data matrix");
98  }
100  &data(0,0),&data(1,0),&data(2,0),&data(3,0),&data(4,0),&data(5,0),&data(6,0),&data(7,0),&data(8,0)
101  );
102 }
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:389
T data[(layout==column_major) ? Tensor_Dim0 :Tensor_Dim1][(layout==column_major) ? Tensor_Dim1 :Tensor_Dim0]

◆ Hcurl_BubbleFaceFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_BubbleFaceFunctions_MBTET ( int faces_nodes,
int p,
double N,
double diffN,
double phi_f[4],
double diff_phi_f[4],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face edge base functions of Hcurl space on face on tetrahedral.

On each face we have P*(P-1) base functions and are 4 faces.

See NBFACETRI_EDGE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Narray shape functions evaluated at each integration point
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 499 of file Hcurl.cpp.

502  {
503 
504  PetscFunctionBegin;
505 
508 
509  // double coords[] = { 0,0,0, 1,0,0, 0,1,0, 0,0,1 };
510  // FTensor::Tensor1<double*,3> t_coords[4] = {
511  // FTensor::Tensor1<double*,3>(&coords[0],&coords[ 1],&coords[ 2]),
512  // FTensor::Tensor1<double*,3>(&coords[3],&coords[ 4],&coords[ 5]),
513  // FTensor::Tensor1<double*,3>(&coords[6],&coords[ 7],&coords[ 8]),
514  // FTensor::Tensor1<double*,3>(&coords[9],&coords[10],&coords[11])
515  // };
516  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
517  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
518  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
519  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
520  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
521  };
522  FTensor::Tensor1<double,3> t_diff_ksi0i,t_diff_ksi0j;
523  FTensor::Tensor1<double,3> diff_beta_0ij;
524 
527 
528  for(int ff = 0;ff!=4;ff++) {
529 
530  if(NBFACETRI_FACE_HCURL(p[ff])==0) continue;
531 
532  int n0 = faces_nodes[3*ff+0];
533  int n1 = faces_nodes[3*ff+1];
534  int n2 = faces_nodes[3*ff+2];
535 
536  // tou_0i(i) = t_coords[n1](i)-t_coords[n0](i);
537  // tou_0j(i) = t_coords[n2](i)-t_coords[n0](i);
538  tou_0i(i) = t_node_diff_ksi[n1](i)-t_node_diff_ksi[n0](i);
539  tou_0j(i) = t_node_diff_ksi[n2](i)-t_node_diff_ksi[n0](i);
540 
541  t_diff_ksi0i(i) = t_node_diff_ksi[n1](i)-t_node_diff_ksi[n0](i);
542  t_diff_ksi0j(i) = t_node_diff_ksi[n2](i)-t_node_diff_ksi[n0](i);
543 
544  double psi_l_0i[p[ff]+1],diff_psi_l_0i[3*p[ff]+3];
545  double psi_l_0j[p[ff]+1],diff_psi_l_0j[3*p[ff]+3];
546 
548  &phi_f[ff][0],&phi_f[ff][1],&phi_f[ff][2],3
549  );
550  FTensor::Tensor2<double*,3,3> t_diff_phi_f(
551  &diff_phi_f[ff][0],&diff_phi_f[ff][3],&diff_phi_f[ff][6],
552  &diff_phi_f[ff][1],&diff_phi_f[ff][4],&diff_phi_f[ff][7],
553  &diff_phi_f[ff][2],&diff_phi_f[ff][5],&diff_phi_f[ff][8],9
554  );
556 
557  for(int ii = 0;ii!=nb_integration_pts;ii++) {
558 
559  const int node_shift = ii*4;
560  const double beta_0ij = N[node_shift+n0]*N[node_shift+n1]*N[node_shift+n2];
561  diff_beta_0ij(i) =
562  t_node_diff_ksi[n0](i)*N[node_shift+n1]*N[node_shift+n2]+
563  N[node_shift+n0]*t_node_diff_ksi[n1](i)*N[node_shift+n2]+
564  N[node_shift+n0]*N[node_shift+n1]*t_node_diff_ksi[n2](i);
565 
566  const double ksi_0i = N[node_shift+n1]-N[node_shift+n0];
567  ierr = base_polynomials(p[ff],ksi_0i,&t_diff_ksi0i(0),psi_l_0i,diff_psi_l_0i,3); CHKERRQ(ierr);
568  const double ksi_0j = N[node_shift+n2]-N[node_shift+n0];
569  ierr = base_polynomials(p[ff],ksi_0j,&t_diff_ksi0j(0),psi_l_0j,diff_psi_l_0j,3); CHKERRQ(ierr);
570 
571  int cc = 0;
572  for(int oo = 0;oo<=(p[ff]-3);oo++) {
573  FTensor::Tensor1<double*,3> t_diff_psi_l_0i(
574  &diff_psi_l_0i[0],&diff_psi_l_0i[p[ff]+1],&diff_psi_l_0i[2*p[ff]+2],1
575  );
576  for(int pp0 = 0;pp0<=oo;pp0++) {
577  const int pp1 = oo-pp0;
578  if(pp1>=0) {
579  FTensor::Tensor1<double*,3> t_diff_psi_l_0j(
580  &diff_psi_l_0j[pp1],&diff_psi_l_0j[p[ff]+1+pp1],&diff_psi_l_0j[2*p[ff]+2+pp1],1
581  );
582  // base functions
583  const double a = beta_0ij*psi_l_0i[pp0]*psi_l_0j[pp1];
584  t_phi_f(i) = a*tou_0i(i);
585  ++t_phi_f;
586  ++cc;
587  t_phi_f(i) = a*tou_0j(i);
588  ++t_phi_f;
589  ++cc;
590  // diff base functions
591  t_b(j) =
592  diff_beta_0ij(j)*psi_l_0i[pp0]*psi_l_0j[pp1]+
593  beta_0ij*t_diff_psi_l_0i(j)*psi_l_0j[pp1]+
594  beta_0ij*psi_l_0i[pp0]*t_diff_psi_l_0j(j);
595  t_diff_phi_f(i,j) = t_b(j)*tou_0i(i);
596  ++t_diff_phi_f;
597  t_diff_phi_f(i,j) = t_b(j)*tou_0j(i);
598  ++t_diff_phi_f;
599  ++t_diff_psi_l_0i;
600  }
601  }
602  }
603  const int nb_base_fun_on_face = NBFACETRI_FACE_HCURL(p[ff]);
604  if(cc!=nb_base_fun_on_face) {
605  SETERRQ2(
606  PETSC_COMM_SELF,
608  "Wrong number of base functions %d != %d",
609  cc,nb_base_fun_on_face
610  );
611  }
612 
613  }
614 
615  }
616  PetscFunctionReturn(0);
617 }
CHKERRQ(ierr)
#define NBFACETRI_FACE_HCURL(P)

◆ Hcurl_BubbleFaceFunctions_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hcurl_BubbleFaceFunctions_MBTET_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f,
double diff_phi_f,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face edge base functions of Hcurl space on face.

On each face we have P*(P-1) base functions and are 4 faces.

See NBFACETRI_EDGE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Narray shape functions evaluated at each integration point
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 619 of file Hcurl.cpp.

622  {
623 
624  PetscFunctionBegin;
625 
626  double zero = 0;
627  FTensor::Tensor2<double*,3,3> t_node_diff_ksi(
628  &diffN[0],&diffN[ 1],&zero,
629  &diffN[2],&diffN[ 3],&zero,
630  &diffN[4],&diffN[ 5],&zero
631  );
632 
635 
636  if(NBFACETRI_FACE_HCURL(p)==0) PetscFunctionReturn(0);
637 
641 
642  const int node0 = faces_nodes[0];
643  const int node1 = faces_nodes[1];
644  const int node2 = faces_nodes[2];
645 
648 
649  tou_0i(i) = t_node_diff_ksi(N1,i)-t_node_diff_ksi(N0,i);;
650  tou_0j(i) = t_node_diff_ksi(N2,i)-t_node_diff_ksi(N0,i);;
651 
652  double psi_l_0i[p+1],psi_l_0j[p+1];
653  double diff_psi_l_0i[2*p+2],diff_psi_l_0j[2*p+2];
654  FTensor::Tensor1<double*,3> t_phi_f(&phi_f[0],&phi_f[1],&phi_f[2],3);
655  FTensor::Tensor2<double*,3,2> t_diff_phi_f(
656  &diff_phi_f[HCURL0_0],&diff_phi_f[HCURL0_1],
657  &diff_phi_f[HCURL1_0],&diff_phi_f[HCURL1_1],
658  &diff_phi_f[HCURL2_0],&diff_phi_f[HCURL2_1],6
659  );
660 
661  double diff_ksi_0i[] = {
662  t_node_diff_ksi(node1,0)-t_node_diff_ksi(node0,0),
663  t_node_diff_ksi(node1,1)-t_node_diff_ksi(node0,1)
664  };
665  double diff_ksi_0j[] = {
666  t_node_diff_ksi(node2,0)-t_node_diff_ksi(node0,0),
667  t_node_diff_ksi(node2,1)-t_node_diff_ksi(node0,1)
668  };
669 
670  for(int ii = 0;ii!=nb_integration_pts;ii++) {
671 
672  const int node_shift = ii*3;
673  const double n0 = N[node_shift+node0];
674  const double n1 = N[node_shift+node1];
675  const double n2 = N[node_shift+node2];
676 
677  const double beta_0ij = n0*n1*n2;
678  FTensor::Tensor1<double,2> diff_beta_0ij(
679  t_node_diff_ksi(node0,0)*n1*n2+n0*t_node_diff_ksi(node1,0)*n2+n0*n1*t_node_diff_ksi(node2,0),
680  t_node_diff_ksi(node0,1)*n1*n2+n0*t_node_diff_ksi(node1,1)*n2+n0*n1*t_node_diff_ksi(node2,1)
681  );
682 
683  const double ksi_0i = N[node_shift+node1]-N[node_shift+node0];
684  ierr = base_polynomials(p,ksi_0i,diff_ksi_0i,psi_l_0i,diff_psi_l_0i,2); CHKERRQ(ierr);
685  const double ksi_0j = N[node_shift+node2]-N[node_shift+node0];
686  ierr = base_polynomials(p,ksi_0j,diff_ksi_0j,psi_l_0j,diff_psi_l_0j,2); CHKERRQ(ierr);
687 
688  int cc = 0;
690  for(int oo = 0;oo<=(p-3);oo++) {
691  for(int pp0 = 0;pp0<=oo;pp0++) {
692  const int pp1 = oo-pp0;
693  if(pp1>=0) {
694  FTensor::Tensor1<double,2> t_diff_psi_l_0i(
695  diff_psi_l_0i[0+pp0],diff_psi_l_0i[p+1+pp0]
696  );
697  FTensor::Tensor1<double,2> t_diff_psi_l_0j(
698  diff_psi_l_0j[0+pp1],diff_psi_l_0j[p+1+pp1]
699  );
700  const double a = beta_0ij*psi_l_0i[pp0]*psi_l_0j[pp1];
701  t_diff_a(j) =
702  diff_beta_0ij(j)*psi_l_0i[pp0]*psi_l_0j[pp1]+
703  beta_0ij*psi_l_0i[pp0]*t_diff_psi_l_0j(j)+
704  beta_0ij*psi_l_0j[pp1]*t_diff_psi_l_0i(j);
705 
706  t_phi_f(i) = a*tou_0i(i);
707  t_diff_phi_f(i,j) = tou_0i(i)*t_diff_a(j);
708  ++t_phi_f;
709  ++t_diff_phi_f;
710  ++cc;
711  t_phi_f(i) = a*tou_0j(i);
712  t_diff_phi_f(i,j) = tou_0j(i)*t_diff_a(j);
713  ++t_phi_f;
714  ++t_diff_phi_f;
715  ++cc;
716  }
717  }
718  }
719 
720  const int nb_base_fun_on_face = NBFACETRI_FACE_HCURL(p);
721  if(cc!=nb_base_fun_on_face) {
722  SETERRQ2(
723  PETSC_COMM_SELF,
725  "Wrong number of base functions %d != %d",
726  cc,nb_base_fun_on_face
727  );
728  }
729 
730  }
731 
732  PetscFunctionReturn(0);
733 }
CHKERRQ(ierr)
#define NBFACETRI_FACE_HCURL(P)

◆ Hcurl_EdgeBasedFaceFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_EdgeBasedFaceFunctions_MBTET ( int faces_nodes,
int p,
double N,
double diffN,
double phi_f_e[4][3],
double diff_phi_f_e[4][3],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face edge base functions of Hcurl space on tetrahedral.

On each edge we have (P-1) base functions, and each face has 3 edges and are 4 faces on tetrahedral.

See NBFACETRI_EDGE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Narray shape functions evaluated at each integration point
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 328 of file Hcurl.cpp.

331  {
332 
333  PetscFunctionBegin;
334  const int edges[3][2] = { {0,1}, {1,2}, {2,0} };
335 
338 
339  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
340  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
341  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
342  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
343  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
344  };
345  FTensor::Tensor1<double,3> t_edge_diff_ksi;
346  FTensor::Tensor1<double,3> t_diff_beta_e;
347 
348  for(int ff = 0;ff!=4;ff++) {
349 
350  const int o_nodes[3] = {
351  faces_nodes[3*ff+2],faces_nodes[3*ff+0],faces_nodes[3*ff+1]
352  };
353  FTensor::Tensor1<double*,3> t_opposite_node_diff[3] = {
354  FTensor::Tensor1<double*,3>(&diffN[3*o_nodes[0]+0],&diffN[3*o_nodes[0]+1],&diffN[3*o_nodes[0]+2]),
355  FTensor::Tensor1<double*,3>(&diffN[3*o_nodes[1]+0],&diffN[3*o_nodes[1]+1],&diffN[3*o_nodes[1]+2]),
356  FTensor::Tensor1<double*,3>(&diffN[3*o_nodes[2]+0],&diffN[3*o_nodes[2]+1],&diffN[3*o_nodes[2]+2])
357  };
358  double psi_l[p[ff]+1],diff_psi_l[3*p[ff]+3];
359 
360  const int nb_base_fun_on_face = NBFACETRI_EDGE_HCURL(p[ff]);
361  // cerr << nb_base_fun_on_face << " " << p[ff] << endl;
362  if(nb_base_fun_on_face==0) continue;
363 
364  for(int ee = 0;ee!=3;ee++) {
365 
366  FTensor::Tensor1<double*,3> t_face_edge_base(
367  &phi_f_e[ff][ee][0],&phi_f_e[ff][ee][1],&phi_f_e[ff][ee][2],3
368  );
369  FTensor::Tensor2<double*,3,3> t_diff_face_edge_base(
370  &diff_phi_f_e[ff][ee][0],&diff_phi_f_e[ff][ee][3],&diff_phi_f_e[ff][ee][6],
371  &diff_phi_f_e[ff][ee][1],&diff_phi_f_e[ff][ee][4],&diff_phi_f_e[ff][ee][7],
372  &diff_phi_f_e[ff][ee][2],&diff_phi_f_e[ff][ee][5],&diff_phi_f_e[ff][ee][8],9
373  );
374  const int en[] = { faces_nodes[3*ff+edges[ee][0]],faces_nodes[3*ff+edges[ee][1]] };
375  t_edge_diff_ksi(i) = t_node_diff_ksi[en[1]](i)-t_node_diff_ksi[en[0]](i);
376 
377  for(int ii = 0;ii!=nb_integration_pts;ii++) {
378 
379  const int node_shift = ii*4;
380  const double n[] = {
381  N[node_shift+faces_nodes[3*ff+edges[ee][0]]],
382  N[node_shift+faces_nodes[3*ff+edges[ee][1]]]
383  };
384  const double ksi_0i = n[1]-n[0];
385  ierr = base_polynomials(p[ff],ksi_0i,&t_edge_diff_ksi(0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
386  FTensor::Tensor1<double*,3> t_diff_psi_l(&diff_psi_l[0],&diff_psi_l[p[ff]+1],&diff_psi_l[2*p[ff]+2],1);
387 
388  const double beta_e = n[0]*n[1];
389  t_diff_beta_e(j) = t_node_diff_ksi[en[0]](j)*n[1]+n[0]*t_node_diff_ksi[en[1]](j);
390 
391  for(int ll = 0;ll!=nb_base_fun_on_face;ll++) {
392  // if(ll == nb_base_fun_on_face-1) cerr << psi_l[ll] << endl;
393 
394  t_face_edge_base(i) = beta_e*psi_l[ll]*t_opposite_node_diff[ee](i);
395  ++t_face_edge_base;
396 
397  t_diff_face_edge_base(i,j) =
398  (beta_e*t_diff_psi_l(j)+t_diff_beta_e(j)*psi_l[ll])*t_opposite_node_diff[ee](i);
399 
400  ++t_diff_face_edge_base;
401  ++t_diff_psi_l;
402  }
403 
404  }
405  }
406  }
407 
408  PetscFunctionReturn(0);
409 }
#define NBFACETRI_EDGE_HCURL(P)
CHKERRQ(ierr)

◆ Hcurl_EdgeBasedFaceFunctions_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hcurl_EdgeBasedFaceFunctions_MBTET_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f_e[3],
double diff_phi_f_e[3],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face edge base functions of Hcurl space.

On each edge we have (P-1) base functions, and each face has 3 edges and are 4 faces on tetrahedral.

See NBFACETRI_EDGE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Narray shape functions evaluated at each integration point
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 411 of file Hcurl.cpp.

414  {
415 
416  PetscFunctionBegin;
417 
418  const int nb_base_fun_on_face = NBFACETRI_EDGE_HCURL(p);
419  if(nb_base_fun_on_face==0) PetscFunctionReturn(0);
420 
421  const int edges[3][2] = { {0,1}, {1,2}, {2,0} };
422 
425 
426  const int o_nodes[3] = { 2,0,1 };
427  FTensor::Tensor2<double,3,3> t_opposite_node_diff(
428  diffN[2*o_nodes[0]+0],diffN[2*o_nodes[0]+1],0,
429  diffN[2*o_nodes[1]+0],diffN[2*o_nodes[1]+1],0,
430  diffN[2*o_nodes[2]+0],diffN[2*o_nodes[2]+1],0
431  );
432  double psi_l[p+1];
433  double diff_psi_l[2*p+2];
434 
435  FTensor::Tensor1<double*,3> t_face_edge_base[3] = {
436  FTensor::Tensor1<double*,3>(&phi_f_e[0][HCURL0],&phi_f_e[0][HCURL1],&phi_f_e[0][HCURL2],3),
437  FTensor::Tensor1<double*,3>(&phi_f_e[1][HCURL0],&phi_f_e[1][HCURL1],&phi_f_e[1][HCURL2],3),
438  FTensor::Tensor1<double*,3>(&phi_f_e[2][HCURL0],&phi_f_e[2][HCURL1],&phi_f_e[2][HCURL2],3),
439  };
440  FTensor::Tensor2<double*,3,2> t_diff_face_edge_base[3] = {
442  &diff_phi_f_e[0][HCURL0_0],&diff_phi_f_e[0][HCURL0_1],
443  &diff_phi_f_e[0][HCURL1_0],&diff_phi_f_e[0][HCURL1_1],
444  &diff_phi_f_e[0][HCURL2_0],&diff_phi_f_e[0][HCURL2_1],6
445  ),
447  &diff_phi_f_e[1][HCURL0_0],&diff_phi_f_e[1][HCURL0_1],
448  &diff_phi_f_e[1][HCURL1_0],&diff_phi_f_e[1][HCURL1_1],
449  &diff_phi_f_e[1][HCURL2_0],&diff_phi_f_e[1][HCURL2_1],6
450  ),
452  &diff_phi_f_e[2][HCURL0_0],&diff_phi_f_e[2][HCURL0_1],
453  &diff_phi_f_e[2][HCURL1_0],&diff_phi_f_e[2][HCURL1_1],
454  &diff_phi_f_e[2][HCURL2_0],&diff_phi_f_e[2][HCURL2_1],6
455  )
456  };
457 
458  for(int ee = 0;ee!=3;ee++) {
459 
460  const int node0 = faces_nodes[edges[ee][0]];
461  const int node1 = faces_nodes[edges[ee][1]];
462  double diff_ksi0i[] = {
463  diffN[2*node1+0]-diffN[2*node0+0],
464  diffN[2*node1+1]-diffN[2*node0+1]
465  };
466 
467  for(int ii = 0;ii!=nb_integration_pts;ii++) {
468 
469  const int node_shift = ii*3;
470  const double n0 = N[node_shift+node0];
471  const double n1 = N[node_shift+node1];
472  const double ksi_0i = n1-n0;
473  ierr = base_polynomials(p,ksi_0i,diff_ksi0i,psi_l,diff_psi_l,2); CHKERRQ(ierr);
474  const double beta_e = n0*n1;
475  FTensor::Tensor1<double,2> t_diff_beta_e(
476  diffN[2*node0+0]*n1+n0*diffN[2*node1+0],
477  diffN[2*node0+1]*n1+n0*diffN[2*node1+1]
478  );
479  FTensor::Tensor1<double*,2> t_diff_psi_l(
480  &diff_psi_l[0],&diff_psi_l[p+1],1
481  );
482 
483  for(int ll = 0;ll!=nb_base_fun_on_face;ll++) {
484  t_face_edge_base[ee](i) = beta_e*psi_l[ll]*t_opposite_node_diff(ee,i);
485  t_diff_face_edge_base[ee](i,j) =
486  beta_e*t_opposite_node_diff(ee,i)*t_diff_psi_l(j)+
487  psi_l[ll]*t_opposite_node_diff(ee,i)*t_diff_beta_e(j);
488  ++t_face_edge_base[ee];
489  ++t_diff_face_edge_base[ee];
490  ++t_diff_psi_l;
491  }
492 
493  }
494  }
495 
496  PetscFunctionReturn(0);
497 }
#define NBFACETRI_EDGE_HCURL(P)
CHKERRQ(ierr)

◆ Hcurl_EdgeBaseFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_EdgeBaseFunctions_MBTET ( int sense,
int p,
double N,
double diffN,
double edgeN[],
double diff_edgeN[],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Edge based H-curl base functions on tetrahedral.

Function generates hierarchical base of h-curl comforting functions on tetrahedral edge. For more details see [2].

On each tetrahedral's edge we have P+1 functions. See NBEDGE_HCURL

Parameters
sensesense fo edge (i.e. unique orientation)
parray of oder for each edge
Narray shape functions evaluated at each integration point
diffNderivatives of shape functions
edgeNbase functions on edges
diff_edgeNderivatives of edge shape functions
nb_integration_ptsnumber of integration points
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 21 of file Hcurl.cpp.

24  {
25  PetscFunctionBegin;
26 
27  const int edges_nodes[6][2] = {
28  {0,1}, {1,2}, {2,0}, {0,3}, {1,3}, {2,3}
29  };
30  int P[6];
31  for(int ee = 0;ee!=6; ee++) P[ee] = NBEDGE_HCURL(p[ee]);
32 
35 
36  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
37  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
38  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
39  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
40  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
41  };
42  double edge_diff_ksi[6][3];
43  FTensor::Tensor1<double*,3> t_edge_diff_ksi[6] = {
44  FTensor::Tensor1<double*,3>(&edge_diff_ksi[0][0],&edge_diff_ksi[0][1],&edge_diff_ksi[0][2]),
45  FTensor::Tensor1<double*,3>(&edge_diff_ksi[1][0],&edge_diff_ksi[1][1],&edge_diff_ksi[1][2]),
46  FTensor::Tensor1<double*,3>(&edge_diff_ksi[2][0],&edge_diff_ksi[2][1],&edge_diff_ksi[2][2]),
47  FTensor::Tensor1<double*,3>(&edge_diff_ksi[3][0],&edge_diff_ksi[3][1],&edge_diff_ksi[3][2]),
48  FTensor::Tensor1<double*,3>(&edge_diff_ksi[4][0],&edge_diff_ksi[4][1],&edge_diff_ksi[4][2]),
49  FTensor::Tensor1<double*,3>(&edge_diff_ksi[5][0],&edge_diff_ksi[5][1],&edge_diff_ksi[5][2])
50  };
51  for(int ee = 0;ee!=6;ee++) {
52  t_edge_diff_ksi[ee](i) = (
53  t_node_diff_ksi[edges_nodes[ee][1]](i)-t_node_diff_ksi[edges_nodes[ee][0]](i)
54  )*sense[ee];
55  }
56 
57  FTensor::Tensor1<double*,3> t_edge_n[6] = {
58  FTensor::Tensor1<double*,3>(&edge_n[0][0],&edge_n[0][1],&edge_n[0][2],3),
59  FTensor::Tensor1<double*,3>(&edge_n[1][0],&edge_n[1][1],&edge_n[1][2],3),
60  FTensor::Tensor1<double*,3>(&edge_n[2][0],&edge_n[2][1],&edge_n[2][2],3),
61  FTensor::Tensor1<double*,3>(&edge_n[3][0],&edge_n[3][1],&edge_n[3][2],3),
62  FTensor::Tensor1<double*,3>(&edge_n[4][0],&edge_n[4][1],&edge_n[4][2],3),
63  FTensor::Tensor1<double*,3>(&edge_n[5][0],&edge_n[5][1],&edge_n[5][2],3)
64  };
65  FTensor::Tensor2<double*,3,3> t_diff_edge_n[6] = {
67  &diff_edge_n[0][0],&diff_edge_n[0][3],&diff_edge_n[0][6],
68  &diff_edge_n[0][1],&diff_edge_n[0][4],&diff_edge_n[0][7],
69  &diff_edge_n[0][2],&diff_edge_n[0][5],&diff_edge_n[0][8],9
70  ),
72  &diff_edge_n[1][0],&diff_edge_n[1][3],&diff_edge_n[1][6],
73  &diff_edge_n[1][1],&diff_edge_n[1][4],&diff_edge_n[1][7],
74  &diff_edge_n[1][2],&diff_edge_n[1][5],&diff_edge_n[1][8],9
75  ),
77  &diff_edge_n[2][0],&diff_edge_n[2][3],&diff_edge_n[2][6],
78  &diff_edge_n[2][1],&diff_edge_n[2][4],&diff_edge_n[2][7],
79  &diff_edge_n[2][2],&diff_edge_n[2][5],&diff_edge_n[2][8],9
80  ),
82  &diff_edge_n[3][0],&diff_edge_n[3][3],&diff_edge_n[3][6],
83  &diff_edge_n[3][1],&diff_edge_n[3][4],&diff_edge_n[3][7],
84  &diff_edge_n[3][2],&diff_edge_n[3][5],&diff_edge_n[3][8],9
85  ),
87  &diff_edge_n[4][0],&diff_edge_n[4][3],&diff_edge_n[4][6],
88  &diff_edge_n[4][1],&diff_edge_n[4][4],&diff_edge_n[4][7],
89  &diff_edge_n[4][2],&diff_edge_n[4][5],&diff_edge_n[4][8],9
90  ),
92  &diff_edge_n[5][0],&diff_edge_n[5][3],&diff_edge_n[5][6],
93  &diff_edge_n[5][1],&diff_edge_n[5][4],&diff_edge_n[5][7],
94  &diff_edge_n[5][2],&diff_edge_n[5][5],&diff_edge_n[5][8],9
95  )
96  };
97  FTensor::Tensor1<double,3> t_psi_e_0,t_psi_e_1;
98  FTensor::Tensor2<double,3,3> t_diff_psi_e_0,t_diff_psi_e_1;
99 
100  for(int ii = 0;ii!=nb_integration_pts;ii++) {
101 
102  const int node_shift = ii*4;
103  for(int ee = 0;ee!=6;ee++) {
104 
105  if(P[ee]==0) continue;
106 
107  t_psi_e_0(i) =
108  (N[node_shift+edges_nodes[ee][1]]*t_node_diff_ksi[edges_nodes[ee][0]](i)-
109  N[node_shift+edges_nodes[ee][0]]*t_node_diff_ksi[edges_nodes[ee][1]](i))*sense[ee];
110  t_diff_psi_e_0(i,j) =
111  (t_node_diff_ksi[edges_nodes[ee][1]](j)*t_node_diff_ksi[edges_nodes[ee][0]](i)-
112  t_node_diff_ksi[edges_nodes[ee][0]](j)*t_node_diff_ksi[edges_nodes[ee][1]](i))*sense[ee];
113 
114  t_psi_e_1(i) =
115  N[node_shift+edges_nodes[ee][1]]*t_node_diff_ksi[edges_nodes[ee][0]](i)+
116  N[node_shift+edges_nodes[ee][0]]*t_node_diff_ksi[edges_nodes[ee][1]](i);
117  t_diff_psi_e_1(i,j) =
118  t_node_diff_ksi[edges_nodes[ee][1]](j)*t_node_diff_ksi[edges_nodes[ee][0]](i)+
119  t_node_diff_ksi[edges_nodes[ee][0]](j)*t_node_diff_ksi[edges_nodes[ee][1]](i);
120 
121  (t_edge_n[ee])(i) = t_psi_e_0(i);
122  ++(t_edge_n[ee]);
123  (t_edge_n[ee])(i) = t_psi_e_1(i);
124  ++(t_edge_n[ee]);
125 
126  (t_diff_edge_n[ee])(i,j) = t_diff_psi_e_0(i,j);
127  ++(t_diff_edge_n[ee]);
128  (t_diff_edge_n[ee])(i,j) = t_diff_psi_e_1(i,j);
129  ++(t_diff_edge_n[ee]);
130 
131  if(p[ee]>1) {
132 
133  const double ksi_0i = (N[node_shift+edges_nodes[ee][1]]-N[node_shift+edges_nodes[ee][0]])*sense[ee];
134  double psi_l[p[ee]+1],diff_psi_l[3*p[ee]+3];
135  ierr = base_polynomials(p[ee],ksi_0i,&edge_diff_ksi[ee][0],psi_l,diff_psi_l,3); CHKERRQ(ierr);
136  FTensor::Tensor1<double*,3> t_diff_psi_l(&diff_psi_l[0],&diff_psi_l[p[ee]+1],&diff_psi_l[2*p[ee]+2],1);
137 
138  for(int ll = 2;ll!=P[ee];ll++) {
139 
140  const double a = (double)(2*ll+1)/(double)(ll+1);
141  const double b = (double)(ll)/(double)(ll+1);
142 
143  (t_edge_n[ee])(i) = a*psi_l[ll-1]*t_psi_e_1(i)-b*psi_l[ll-2]*t_psi_e_0(i);
144  ++(t_edge_n[ee]);
145 
146  (t_diff_edge_n[ee])(i,j) = -b*(t_diff_psi_l(j)*t_psi_e_0(i)+psi_l[ll-2]*t_diff_psi_e_0(i,j));
147  ++t_diff_psi_l;
148  (t_diff_edge_n[ee])(i,j) += a*(t_diff_psi_l(j)*t_psi_e_1(i)+psi_l[ll-1]*t_diff_psi_e_1(i,j));
149  ++(t_diff_edge_n[ee]);
150 
151  }
152  }
153 
154  }
155  }
156 
157  PetscFunctionReturn(0);
158 }
#define NBEDGE_HCURL(P)
Number of base functions H curl on faces.
CHKERRQ(ierr)

◆ Hcurl_EdgeBaseFunctions_MBTET_ON_EDGE()

PetscErrorCode MoFEM::Hcurl_EdgeBaseFunctions_MBTET_ON_EDGE ( int  sense,
int  p,
double N,
double diffN,
double edgeN,
double diff_edgeN,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Edge based H-curl base functions on edge.

Function generates hierarchical base of h-curl comforting functions on tetrahedral edge. For more details see [2].

On each edge we have P+1 functions. See NBEDGE_HCURL

Parameters
sensesense fo edge (i.e. unique orientation)
parray of oder for each edge
Narray shape functions evaluated at each integration point
diffNderivatives of shape functions
edgeNbase functions on edges
diff_edgeNderivatives of edge shape functions
nb_integration_ptsnumber of integration points
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 160 of file Hcurl.cpp.

163  {
164  PetscFunctionBegin;
165 
166 
167  if(NBEDGE_HCURL(p)==0) PetscFunctionReturn(0);
168  if(diff_edge_n!=NULL) {
169  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Calculation of derivatives not implemented");
170  }
171 
173  FTensor::Tensor1<double,3> t_node_diff_ksi[2];
174  t_node_diff_ksi[0](0) = diffN[0];
175  t_node_diff_ksi[0](1) = 0;
176  t_node_diff_ksi[0](2) = 0;
177  t_node_diff_ksi[1](0) = diffN[1];
178  t_node_diff_ksi[1](1) = 0;
179  t_node_diff_ksi[1](2) = 0;
180 
181  FTensor::Tensor1<double*,3> t_edge_n(&edge_n[0],&edge_n[1],&edge_n[2],3);
182  FTensor::Tensor1<double,3> t_psi_e_0,t_psi_e_1;
183 
184  for(int ii = 0;ii!=nb_integration_pts;ii++) {
185 
186  const int node_shift = ii*2;
187 
188  t_psi_e_0(i) = (N[node_shift+1]*t_node_diff_ksi[0](i)-N[node_shift+0]*t_node_diff_ksi[1](i))*sense;
189  t_psi_e_1(i) = N[node_shift+1]*t_node_diff_ksi[0](i)+N[node_shift+0]*t_node_diff_ksi[1](i);
190  t_edge_n(i) = t_psi_e_0(i);
191  ++t_edge_n;
192  t_edge_n(i) = t_psi_e_1(i);
193  ++t_edge_n;
194 
195  if(p>1) {
196  const double ksi_0i = (N[node_shift+1]-N[node_shift+0])*sense;
197  double psi_l[p+1];
198  ierr = base_polynomials(p,ksi_0i,NULL,psi_l,NULL,3); CHKERRQ(ierr);
199  for(int ll = 2;ll!=NBEDGE_HCURL(p);ll++) {
200  const double a = (double)(2*ll+1)/(double)(ll+1);
201  const double b = (double)(ll)/(double)(ll+1);
202  t_edge_n(i) = a*psi_l[ll-1]*t_psi_e_1(i)-b*psi_l[ll-2]*t_psi_e_0(i);
203  ++t_edge_n;
204  }
205  }
206 
207  }
208 
209  PetscFunctionReturn(0);
210 }
#define NBEDGE_HCURL(P)
Number of base functions H curl on faces.
CHKERRQ(ierr)

◆ Hcurl_EdgeBaseFunctions_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hcurl_EdgeBaseFunctions_MBTET_ON_FACE ( int sense,
int p,
double N,
double diffN,
double edgeN[],
double diff_edgeN[],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Edge based H-curl base functions on face.

Function generates hierarchical base of h-curl comforting functions on tetrahedral edge. For more details see [2].

On each face's edge we have P+1 functions. See NBEDGE_HCURL

Parameters
sensesense fo edge (i.e. unique orientation)
parray of oder for each edge
Narray shape functions evaluated at each integration point
diffNderivatives of shape functions
edgeNbase functions on edges
diff_edgeNderivatives of edge shape functions
nb_integration_ptsnumber of integration points
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 212 of file Hcurl.cpp.

215  {
216 
217  PetscFunctionBegin;
218 
219  // TODO This is not by atom tests properly
220 
221  const int edges_nodes[3][2] = { {0,1}, {1,2}, {2,0} };
222  int P[3];
223  for(int ee = 0;ee<3; ee++) P[ee] = NBEDGE_HCURL(p[ee]);
224 
227 
228  FTensor::Tensor1<double,3> t_node_diff_ksi[3] = {
229  FTensor::Tensor1<double,3>(diffN[0],diffN[1],0),
230  FTensor::Tensor1<double,3>(diffN[2],diffN[3],0),
231  FTensor::Tensor1<double,3>(diffN[4],diffN[5],0),
232  };
233  FTensor::Tensor1<double,2> t_2d_diff_ksi[3] = {
234  FTensor::Tensor1<double,2>(diffN[0],diffN[1]),
235  FTensor::Tensor1<double,2>(diffN[2],diffN[3]),
236  FTensor::Tensor1<double,2>(diffN[4],diffN[5])
237  };
238 
239  FTensor::Tensor1<double*,3> t_edge_n[3] = {
240  FTensor::Tensor1<double*,3>(&edge_n[0][0],&edge_n[0][1],&edge_n[0][2],3),
241  FTensor::Tensor1<double*,3>(&edge_n[1][0],&edge_n[1][1],&edge_n[1][2],3),
242  FTensor::Tensor1<double*,3>(&edge_n[2][0],&edge_n[2][1],&edge_n[2][2],3)
243  };
244  FTensor::Tensor2<double*,3,2> t_diff_edge_n[3] = {
246  &diff_edge_n[0][HCURL0_0],&diff_edge_n[0][HCURL0_1],
247  &diff_edge_n[0][HCURL1_0],&diff_edge_n[0][HCURL1_1],
248  &diff_edge_n[0][HCURL2_0],&diff_edge_n[0][HCURL2_1],6
249  ),
251  &diff_edge_n[1][HCURL0_0],&diff_edge_n[1][HCURL0_1],
252  &diff_edge_n[1][HCURL1_0],&diff_edge_n[1][HCURL1_1],
253  &diff_edge_n[1][HCURL2_0],&diff_edge_n[1][HCURL2_1],6
254  ),
256  &diff_edge_n[2][HCURL0_0],&diff_edge_n[2][HCURL0_1],
257  &diff_edge_n[2][HCURL1_0],&diff_edge_n[2][HCURL1_1],
258  &diff_edge_n[2][HCURL2_0],&diff_edge_n[2][HCURL2_1],6
259  )
260  };
261 
262  FTensor::Tensor1<double,3> t_psi_e_0,t_psi_e_1;
263  FTensor::Tensor2<double,3,2> t_diff_psi_e_0,t_diff_psi_e_1;
264 
265  for(int ee = 0;ee!=3;ee++) {
266 
267  if(P[ee]==0) continue;
268  const int node0 = edges_nodes[ee][0];
269  const int node1 = edges_nodes[ee][1];
270  const int sense_edge = sense[ee];
271 
272  t_diff_psi_e_0(i,j) =
273  (t_node_diff_ksi[node0](i)*t_2d_diff_ksi[node1](j)-
274  t_node_diff_ksi[node1](i)*t_2d_diff_ksi[node0](j))*sense_edge;
275  t_diff_psi_e_1(i,j) =
276  t_node_diff_ksi[node0](i)*t_2d_diff_ksi[node1](j)+
277  t_node_diff_ksi[node1](i)*t_2d_diff_ksi[node0](j);
278 
279  for(int ii = 0;ii!=nb_integration_pts;ii++) {
280 
281  const int node_shift = ii*3;
282  const double n0 = N[node_shift+node0];
283  const double n1 = N[node_shift+node1];
284 
285  t_psi_e_0(i) = (n1*t_node_diff_ksi[node0](i)-n0*t_node_diff_ksi[node1](i))*sense_edge;
286  t_psi_e_1(i) = n1*t_node_diff_ksi[node0](i)+n0*t_node_diff_ksi[node1](i);
287 
288  (t_edge_n[ee])(i) = t_psi_e_0(i);
289  (t_diff_edge_n[ee])(i,j) = t_diff_psi_e_0(i,j);
290  ++(t_edge_n[ee]);
291  ++(t_diff_edge_n[ee]);
292  (t_edge_n[ee])(i) = t_psi_e_1(i);
293  (t_diff_edge_n[ee])(i,j) = t_diff_psi_e_1(i,j);
294  ++(t_edge_n[ee]);
295  ++(t_diff_edge_n[ee]);
296 
297  if(p[ee]>1) {
298  const double ksi_0i = (n1-n0)*sense_edge;
299  double diff_ksi_0i[] = {
300  ((t_2d_diff_ksi[node1])(0)-(t_2d_diff_ksi[node0])(0))*sense_edge,
301  ((t_2d_diff_ksi[node1])(1)-(t_2d_diff_ksi[node0])(1))*sense_edge
302  };
303 
304  double psi_l[p[ee]+1],diff_psi_l[2*p[ee]+2];
305  ierr = base_polynomials(p[ee],ksi_0i,diff_ksi_0i,psi_l,diff_psi_l,2); CHKERRQ(ierr);
306  FTensor::Tensor1<double*,2> t_diff_psi_ll_m1(&diff_psi_l[0+2-1],&diff_psi_l[p[ee]+1+2-1],1);
307  FTensor::Tensor1<double*,2> t_diff_psi_ll_m2(&diff_psi_l[0+2-2],&diff_psi_l[p[ee]+1+2-2],1);
308  for(int ll = 2;ll!=P[ee];ll++) {
309  const double a = (double)(2*ll+1)/(double)(ll+1);
310  const double b = (double)(ll)/(double)(ll+1);
311  (t_edge_n[ee])(i) = a*psi_l[ll-1]*t_psi_e_1(i)-b*psi_l[ll-2]*t_psi_e_0(i);
312  (t_diff_edge_n[ee])(i,j) =
313  a*t_psi_e_1(i)*t_diff_psi_ll_m1(j)+a*psi_l[ll-1]*t_diff_psi_e_1(i,j)-
314  b*t_psi_e_0(i)*t_diff_psi_ll_m2(j)-b*psi_l[ll-2]*t_diff_psi_e_0(i,j);
315  ++(t_edge_n[ee]);
316  ++(t_diff_edge_n[ee]);
317  ++t_diff_psi_ll_m1;
318  ++t_diff_psi_ll_m2;
319  }
320  }
321 
322  }
323  }
324 
325  PetscFunctionReturn(0);
326 }
#define NBEDGE_HCURL(P)
Number of base functions H curl on faces.
CHKERRQ(ierr)

◆ Hcurl_FaceFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_FaceFunctions_MBTET ( int face_nodes,
int p,
double N,
double diffN,
double phi_f[4],
double diff_phi_f[4],
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face H-curl functions.

Face H-curl functions are set of Eddge-Based Face functions (Hcurl_EdgeBasedFaceFunctions_MBTET) and Bubble-Face functions (Hcurl_BubbleFaceFunctions_MBTET).

See NBVOLUMETET_FACE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Nnodal shape functions
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 1018 of file Hcurl.cpp.

1021  {
1022 
1023  PetscFunctionBegin;
1024 
1025  try {
1026 
1027  MatrixDouble base_face_edge_functions[4];
1028  MatrixDouble diff_base_face_edge_functions[4];
1029  double *phi_f_e[4][3];
1030  double *diff_phi_f_e[4][3];
1031  for(int ff = 0;ff!=4;ff++) {
1032  if(NBFACETRI_EDGE_HCURL(p[ff])==0) {
1033  for(int ee = 0;ee!=3;ee++) {
1034  phi_f_e[ff][ee] = NULL;
1035  diff_phi_f_e[ff][ee] = NULL;
1036  }
1037  } else {
1038  base_face_edge_functions[ff].resize(3,3*NBFACETRI_EDGE_HCURL(p[ff])*nb_integration_pts);
1039  diff_base_face_edge_functions[ff].resize(3,9*NBFACETRI_EDGE_HCURL(p[ff])*nb_integration_pts);
1040  // base_face_edge_functions[ff].clear();
1041  // diff_base_face_edge_functions[ff].clear();
1042  for(int ee = 0;ee!=3;ee++) {
1043  phi_f_e[ff][ee] = &base_face_edge_functions[ff](ee,0);
1044  diff_phi_f_e[ff][ee] = &diff_base_face_edge_functions[ff](ee,0);
1045  }
1046  }
1047  }
1049  face_nodes,p,N,diffN,phi_f_e,diff_phi_f_e,nb_integration_pts,base_polynomials
1050  ); CHKERRQ(ierr);
1051 
1052  VectorDouble base_face_bubble_functions[4];
1053  VectorDouble diff_base_face_bubble_functions[4];
1054  double *phi_f_f[4];
1055  double *diff_phi_f_f[4];
1056  for(int ff=0;ff!=4;ff++) {
1057  int nb_dofs = NBFACETRI_FACE_HCURL(p[ff]);
1058  if(nb_dofs==0) {
1059  phi_f_f[ff] = NULL;
1060  diff_phi_f_f[ff] = NULL;
1061  } else {
1062  base_face_bubble_functions[ff].resize(3*nb_dofs*nb_integration_pts,false);
1063  diff_base_face_bubble_functions[ff].resize(9*nb_dofs*nb_integration_pts,false);
1064  phi_f_f[ff] = &*base_face_bubble_functions[ff].data().begin();
1065  diff_phi_f_f[ff] = &*diff_base_face_bubble_functions[ff].data().begin();
1066  }
1067  }
1069  face_nodes,p,N,diffN,phi_f_f,diff_phi_f_f,nb_integration_pts,base_polynomials
1070  ); CHKERRQ(ierr);
1071 
1074 
1075  for(int ff = 0;ff!=4;ff++) {
1076 
1077  if(NBFACETRI_EDGE_HCURL(p[ff])==0) continue;
1078 
1079  FTensor::Tensor1<double*,3> t_face_edge_base[] = {
1080  FTensor::Tensor1<double*,3>(&phi_f_e[ff][0][0],&phi_f_e[ff][0][1],&phi_f_e[ff][0][2],3),
1081  FTensor::Tensor1<double*,3>(&phi_f_e[ff][1][0],&phi_f_e[ff][1][1],&phi_f_e[ff][1][2],3),
1082  FTensor::Tensor1<double*,3>(&phi_f_e[ff][2][0],&phi_f_e[ff][2][1],&phi_f_e[ff][2][2],3)
1083  };
1084  FTensor::Tensor2<double*,3,3> t_diff_face_edge_base[]= {
1086  &diff_phi_f_e[ff][0][0],&diff_phi_f_e[ff][0][3],&diff_phi_f_e[ff][0][6],
1087  &diff_phi_f_e[ff][0][1],&diff_phi_f_e[ff][0][4],&diff_phi_f_e[ff][0][7],
1088  &diff_phi_f_e[ff][0][2],&diff_phi_f_e[ff][0][5],&diff_phi_f_e[ff][0][8],9
1089  ),
1091  &diff_phi_f_e[ff][1][0],&diff_phi_f_e[ff][1][3],&diff_phi_f_e[ff][1][6],
1092  &diff_phi_f_e[ff][1][1],&diff_phi_f_e[ff][1][4],&diff_phi_f_e[ff][1][7],
1093  &diff_phi_f_e[ff][1][2],&diff_phi_f_e[ff][1][5],&diff_phi_f_e[ff][1][8],9
1094  ),
1096  &diff_phi_f_e[ff][2][0],&diff_phi_f_e[ff][2][3],&diff_phi_f_e[ff][2][6],
1097  &diff_phi_f_e[ff][2][1],&diff_phi_f_e[ff][2][4],&diff_phi_f_e[ff][2][7],
1098  &diff_phi_f_e[ff][2][2],&diff_phi_f_e[ff][2][5],&diff_phi_f_e[ff][2][8],9
1099  )
1100  };
1101 
1102  FTensor::Tensor1<double*,3> t_face_base(&phi_f[ff][0],&phi_f[ff][1],&phi_f[ff][2],3);
1103  FTensor::Tensor2<double*,3,3> t_diff_face_base(
1104  &diff_phi_f[ff][0],&diff_phi_f[ff][3],&diff_phi_f[ff][6],
1105  &diff_phi_f[ff][1],&diff_phi_f[ff][4],&diff_phi_f[ff][7],
1106  &diff_phi_f[ff][2],&diff_phi_f[ff][5],&diff_phi_f[ff][8],9
1107  );
1108 
1109  if(NBFACETRI_FACE_HCURL(p[ff])>0) {
1110  FTensor::Tensor1<double*,3> t_face_face_base(&phi_f_f[ff][0],&phi_f_f[ff][1],&phi_f_f[ff][2],3);
1111  FTensor::Tensor2<double*,3,3> t_diff_face_face_base(
1112  &diff_phi_f_f[ff][0],&diff_phi_f_f[ff][3],&diff_phi_f_f[ff][6],
1113  &diff_phi_f_f[ff][1],&diff_phi_f_f[ff][4],&diff_phi_f_f[ff][7],
1114  &diff_phi_f_f[ff][2],&diff_phi_f_f[ff][5],&diff_phi_f_f[ff][8],9
1115  );
1116  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1117  int cc = 0;
1118  for(int oo = 0;oo<=p[ff];oo++) {
1119  // Face-edge base
1120  if(oo>1) {
1121  for(int ee = 0;ee!=3;ee++) {
1122  for(int ll = NBFACETRI_EDGE_HCURL(oo-1);ll!=NBFACETRI_EDGE_HCURL(oo);ll++) {
1123  t_face_base(i) = t_face_edge_base[ee](i);
1124  ++cc;
1125  ++t_face_base;
1126  ++t_face_edge_base[ee];
1127  t_diff_face_base(i,j) = t_diff_face_edge_base[ee](i,j);
1128  ++t_diff_face_base;
1129  ++t_diff_face_edge_base[ee];
1130  // cerr << oo << " " << ll << " " << cc << " " << NBFACETRI_EDGE_HCURL(oo) << endl;
1131  }
1132  }
1133  }
1134  // Face-face base
1135  if(oo>2) {
1136  for(int ll = NBFACETRI_FACE_HCURL(oo-1);ll!=NBFACETRI_FACE_HCURL(oo);ll++) {
1137  t_face_base(i) = t_face_face_base(i);
1138  ++cc;
1139  ++t_face_base;
1140  ++t_face_face_base;
1141  t_diff_face_base(i,j) = t_diff_face_face_base(i,j);
1142  ++t_diff_face_base;
1143  ++t_diff_face_face_base;
1144  }
1145  }
1146  }
1147  // check consistency
1148  const int nb_base_fun_on_face = NBFACETRI_HCURL(p[ff]);
1149  if(cc!=nb_base_fun_on_face) {
1150  SETERRQ2(
1151  PETSC_COMM_SELF,
1153  "Wrong number of base functions %d != %d",
1154  cc,nb_base_fun_on_face
1155  );
1156  }
1157  }
1158  } else {
1159  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1160  int cc = 0;
1161  for(int oo = 0;oo<=p[ff];oo++) {
1162  // Face-edge base
1163  if(oo>1) {
1164  for(int ee = 0;ee!=3;ee++) {
1165  for(int ll = NBFACETRI_EDGE_HCURL(oo-1);ll!=NBFACETRI_EDGE_HCURL(oo);ll++) {
1166  t_face_base(i) = t_face_edge_base[ee](i);
1167  ++cc;
1168  ++t_face_base;
1169  ++t_face_edge_base[ee];
1170  t_diff_face_base(i,j) = t_diff_face_edge_base[ee](i,j);
1171  ++t_diff_face_base;
1172  ++t_diff_face_edge_base[ee];
1173  // cerr << oo << " " << ll << " " << cc << " " << NBFACETRI_EDGE_HCURL(oo) << endl;
1174  }
1175  }
1176  }
1177  }
1178  // check consistency
1179  const int nb_base_fun_on_face = NBFACETRI_HCURL(p[ff]);
1180  if(cc!=nb_base_fun_on_face) {
1181  SETERRQ2(
1182  PETSC_COMM_SELF,
1184  "Wrong number of base functions %d != %d",
1185  cc,nb_base_fun_on_face
1186  );
1187  }
1188  }
1189  }
1190  }
1191 
1192  } catch (MoFEMException const &e) {
1193  SETERRQ(PETSC_COMM_SELF,e.errorCode,e.errorMessage);
1194  } catch (std::exception& ex) {
1195  std::ostringstream ss;
1196  ss << "thorw in method: " << ex.what() << " at line " << __LINE__ << " in file " << __FILE__;
1197  SETERRQ(PETSC_COMM_SELF,MOFEM_STD_EXCEPTION_THROW,ss.str().c_str());
1198  }
1199 
1200  PetscFunctionReturn(0);
1201 }
MoFEMErrorCodes errorCode
Definition: Common.hpp:105
PetscErrorCode Hcurl_BubbleFaceFunctions_MBTET(int *faces_nodes, int *p, double *N, double *diffN, double *phi_f[4], double *diff_phi_f[4], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face edge base functions of Hcurl space on face on tetrahedral.
Definition: Hcurl.cpp:499
ublas::matrix< double, ublas::row_major, DoubleAllacator > MatrixDouble
Definition: Common.hpp:136
#define NBFACETRI_EDGE_HCURL(P)
char errorMessage[255]
Definition: Common.hpp:106
#define NBFACETRI_HCURL(P)
PetscErrorCode Hcurl_EdgeBasedFaceFunctions_MBTET(int *faces_nodes, int *p, double *N, double *diffN, double *phi_f_e[4][3], double *diff_phi_f_e[4][3], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face edge base functions of Hcurl space on tetrahedral.
Definition: Hcurl.cpp:328
CHKERRQ(ierr)
#define NBFACETRI_FACE_HCURL(P)
ublas::vector< double, DoubleAllacator > VectorDouble
Definition: Common.hpp:135

◆ Hcurl_FaceFunctions_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hcurl_FaceFunctions_MBTET_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f,
double diff_phi_f,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face H-curl functions.

Face H-curl functions are set of Eddge-Based Face functions (Hcurl_EdgeBasedFaceFunctions_MBTET) and Bubble-Face functions (Hcurl_BubbleFaceFunctions_MBTET).

See NBVOLUMETET_FACE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Nnodal shape functions
diffNderivatives of nodal shape functions
phi_f[4]calculated shape functions for each face
diff_phi_v[4]derivatives of shape functions for each face
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 1203 of file Hcurl.cpp.

1206  {
1207 
1208  PetscFunctionBegin;
1209 
1210  if(NBFACETRI_EDGE_HCURL(p)==0) PetscFunctionReturn(0);
1211 
1212  MatrixDouble base_face_edge_functions,diff_base_face_edge_functions;
1213  double *phi_f_e[3];
1214  double *diff_phi_f_e[3];
1215  base_face_edge_functions.resize(
1216  3,3*NBFACETRI_EDGE_HCURL(p)*nb_integration_pts
1217  );
1218  diff_base_face_edge_functions.resize(
1219  3,2*3*NBFACETRI_EDGE_HCURL(p)*nb_integration_pts
1220  );
1221  // base_face_edge_functions.clear();
1222  for(int ee = 0;ee!=3;ee++) {
1223  phi_f_e[ee] = &base_face_edge_functions(ee,0);
1224  diff_phi_f_e[ee] = &diff_base_face_edge_functions(ee,0);
1225  }
1227  faces_nodes,p,N,diffN,phi_f_e,diff_phi_f_e,nb_integration_pts,base_polynomials
1228  ); CHKERRQ(ierr);
1229 
1230  VectorDouble base_face_bubble_functions;
1231  VectorDouble diff_base_face_bubble_functions;
1232  double *phi_f_f,*diff_phi_f_f;
1233  base_face_bubble_functions.resize(3*NBFACETRI_FACE_HCURL(p)*nb_integration_pts);
1234  diff_base_face_bubble_functions.resize(2*3*NBFACETRI_FACE_HCURL(p)*nb_integration_pts);
1235  phi_f_f = &*base_face_bubble_functions.data().begin();
1236  diff_phi_f_f = &*diff_base_face_bubble_functions.data().begin();
1238  faces_nodes,p,N,diffN,phi_f_f,diff_phi_f_f,nb_integration_pts,base_polynomials
1239  ); CHKERRQ(ierr);
1240  // cerr << diff_base_face_bubble_functions << endl;
1241 
1244 
1245  FTensor::Tensor1<double*,3> t_face_edge_base[]= {
1246  FTensor::Tensor1<double*,3>(&phi_f_e[0][HCURL0],&phi_f_e[0][HCURL1],&phi_f_e[0][HCURL2],3),
1247  FTensor::Tensor1<double*,3>(&phi_f_e[1][HCURL0],&phi_f_e[1][HCURL1],&phi_f_e[1][HCURL2],3),
1248  FTensor::Tensor1<double*,3>(&phi_f_e[2][HCURL0],&phi_f_e[2][HCURL1],&phi_f_e[2][HCURL2],3)
1249  };
1250  FTensor::Tensor2<double*,3,2> t_diff_face_edge_base[]= {
1252  &diff_phi_f_e[0][HCURL0_0],&diff_phi_f_e[0][HCURL0_1],
1253  &diff_phi_f_e[0][HCURL1_0],&diff_phi_f_e[0][HCURL1_1],
1254  &diff_phi_f_e[0][HCURL2_0],&diff_phi_f_e[0][HCURL2_1],6
1255  ),
1257  &diff_phi_f_e[1][HCURL0_0],&diff_phi_f_e[1][HCURL0_1],
1258  &diff_phi_f_e[1][HCURL1_0],&diff_phi_f_e[1][HCURL1_1],
1259  &diff_phi_f_e[1][HCURL2_0],&diff_phi_f_e[1][HCURL2_1],6
1260  ),
1262  &diff_phi_f_e[2][HCURL0_0],&diff_phi_f_e[2][HCURL0_1],
1263  &diff_phi_f_e[2][HCURL1_0],&diff_phi_f_e[2][HCURL1_1],
1264  &diff_phi_f_e[2][HCURL2_0],&diff_phi_f_e[2][HCURL2_1],6
1265  )
1266  };
1267 
1268  FTensor::Tensor1<double*,3> t_face_base(&phi_f[0],&phi_f[1],&phi_f[2],3);
1269  FTensor::Tensor2<double*,3,2> t_diff_face_base(
1270  &diff_phi_f[HCURL0_0],&diff_phi_f[HCURL0_1],
1271  &diff_phi_f[HCURL1_0],&diff_phi_f[HCURL1_1],
1272  &diff_phi_f[HCURL2_0],&diff_phi_f[HCURL2_1],6
1273  );
1274 
1275  if(NBFACETRI_FACE_HCURL(p)>0) {
1276  FTensor::Tensor1<double*,3> t_face_face_base(
1277  &phi_f_f[HCURL0],&phi_f_f[HCURL1],&phi_f_f[HCURL2],3
1278  );
1279  FTensor::Tensor2<double*,3,2> t_diff_face_face_base(
1280  &diff_phi_f_f[HCURL0_0],&diff_phi_f_f[HCURL0_1],
1281  &diff_phi_f_f[HCURL1_0],&diff_phi_f_f[HCURL1_1],
1282  &diff_phi_f_f[HCURL2_0],&diff_phi_f_f[HCURL2_1],6
1283  );
1284  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1285  int cc = 0;
1286  for(int oo = 0;oo<=p;oo++) {
1287  // Face-edge base
1288  if(oo>1) {
1289  for(int ee = 0;ee!=3;ee++) {
1290  for(int ll = NBFACETRI_EDGE_HCURL(oo-1);ll!=NBFACETRI_EDGE_HCURL(oo);ll++) {
1291  t_face_base(i) = t_face_edge_base[ee](i);
1292  t_diff_face_base(i,j) = t_diff_face_edge_base[ee](i,j);
1293  ++cc;
1294  ++t_face_base;
1295  ++t_face_edge_base[ee];
1296  ++t_diff_face_base;
1297  ++t_diff_face_edge_base[ee];
1298  }
1299  }
1300  }
1301  // Face-face base
1302  if(oo>2) {
1303  for(int ll = NBFACETRI_FACE_HCURL(oo-1);ll!=NBFACETRI_FACE_HCURL(oo);ll++) {
1304  t_face_base(i) = t_face_face_base(i);
1305  t_diff_face_base(i,j) = t_diff_face_face_base(i,j);
1306  ++cc;
1307  ++t_face_base;
1308  ++t_face_face_base;
1309  ++t_diff_face_base;
1310  ++t_diff_face_face_base;
1311  }
1312  }
1313  }
1314  // check consistency
1315  const int nb_base_fun_on_face = NBFACETRI_HCURL(p);
1316  if(cc!=nb_base_fun_on_face) {
1317  SETERRQ2(
1318  PETSC_COMM_SELF,
1320  "Wrong number of base functions %d != %d",
1321  cc,nb_base_fun_on_face
1322  );
1323  }
1324  }
1325  } else {
1326  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1327  int cc = 0;
1328  for(int oo = 0;oo<=p;oo++) {
1329  // Face-edge base
1330  if(oo>1) {
1331  for(int ee = 0;ee!=3;ee++) {
1332  for(int ll = NBFACETRI_EDGE_HCURL(oo-1);ll!=NBFACETRI_EDGE_HCURL(oo);ll++) {
1333  t_face_base(i) = t_face_edge_base[ee](i);
1334  t_diff_face_base(i,j) = t_diff_face_edge_base[ee](i,j);
1335  ++cc;
1336  ++t_face_base;
1337  ++t_face_edge_base[ee];
1338  ++t_diff_face_base;
1339  ++t_diff_face_edge_base[ee];
1340  // cerr << oo << " " << ll << " " << cc << " " << NBFACETRI_EDGE_HCURL(oo) << endl;
1341  }
1342  }
1343  }
1344  }
1345  // check consistency
1346  const int nb_base_fun_on_face = NBFACETRI_HCURL(p);
1347  if(cc!=nb_base_fun_on_face) {
1348  SETERRQ2(
1349  PETSC_COMM_SELF,
1351  "Wrong number of base functions %d != %d",
1352  cc,nb_base_fun_on_face
1353  );
1354  }
1355  }
1356  }
1357 
1358  PetscFunctionReturn(0);
1359 }
ublas::matrix< double, ublas::row_major, DoubleAllacator > MatrixDouble
Definition: Common.hpp:136
#define NBFACETRI_EDGE_HCURL(P)
#define NBFACETRI_HCURL(P)
PetscErrorCode Hcurl_BubbleFaceFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face edge base functions of Hcurl space on face.
Definition: Hcurl.cpp:619
CHKERRQ(ierr)
#define NBFACETRI_FACE_HCURL(P)
ublas::vector< double, DoubleAllacator > VectorDouble
Definition: Common.hpp:135
PetscErrorCode Hcurl_EdgeBasedFaceFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face edge base functions of Hcurl space.
Definition: Hcurl.cpp:411

◆ Hcurl_FaceInteriorFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_FaceInteriorFunctions_MBTET ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_v,
double diff_phi_v,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face base interior function.

On each face we have P*(P-1)/2 and are 4 faces on Tetrahedral.

See NBFACETRI_FACE_HCURL

Parameters
face_nodesarray [4*3] of local indices of face nodes
papproximation order
Nnodal shape functions
diffNderivatives of nodal shape functions
phi_vcalculated shape functions
diff_phi_vderivatives of shape functions
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 735 of file Hcurl.cpp.

738  {
739 
740  PetscFunctionBegin;
741 
742  if(NBVOLUMETET_FACE_HCURL(p)==0) PetscFunctionReturn(0);
743 
744  const int face_opposite_nodes[] = { 2,0,1,3 };
745 
748 
749  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
750  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
751  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
752  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
753  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
754  };
755  FTensor::Tensor1<double,3> t_diff_ksi0i,t_diff_ksi0j;
756 
757  MatrixDouble m_psi_l_0i(4,p+1);
758  MatrixDouble m_psi_l_0j(4,p+1);
759  MatrixDouble m_diff_psi_l_0i(4,3*p+3);
760  MatrixDouble m_diff_psi_l_0j(4,3*p+3);
761 
762  double* psi_l_0i[] = {
763  &m_psi_l_0i(0,0),&m_psi_l_0i(1,0),&m_psi_l_0i(2,0),&m_psi_l_0i(3,0)
764  };
765  double* psi_l_0j[] = {
766  &m_psi_l_0j(0,0),&m_psi_l_0j(1,0),&m_psi_l_0j(2,0),&m_psi_l_0j(3,0)
767  };
768  double* diff_psi_l_0i[] = {
769  &m_diff_psi_l_0i(0,0),&m_diff_psi_l_0i(1,0),&m_diff_psi_l_0i(2,0),&m_diff_psi_l_0i(3,0)
770  };
771  double* diff_psi_l_0j[] = {
772  &m_diff_psi_l_0j(0,0),&m_diff_psi_l_0j(1,0),&m_diff_psi_l_0j(2,0),&m_diff_psi_l_0j(3,0)
773  };
774  double beta_f[4];
775 
776  FTensor::Tensor1<double,3> t_diff_beta_f[4];
777 
778  FTensor::Tensor1<double*,3> t_phi_v(&phi_v[0],&phi_v[1],&phi_v[2],3);
779  FTensor::Tensor2<double*,3,3> t_diff_phi_v(
780  &diff_phi_v[0],&diff_phi_v[3],&diff_phi_v[6],
781  &diff_phi_v[1],&diff_phi_v[4],&diff_phi_v[7],
782  &diff_phi_v[2],&diff_phi_v[5],&diff_phi_v[8],9
783  );
784 
785  for(int ii = 0;ii!=nb_integration_pts;ii++) {
786 
787  for(int ff = 0;ff!=4;ff++) {
788 
789  t_diff_ksi0i(i) =
790  t_node_diff_ksi[faces_nodes[3*ff+1]](i)-t_node_diff_ksi[faces_nodes[3*ff+0]](i);
791  t_diff_ksi0j(i) =
792  t_node_diff_ksi[faces_nodes[3*ff+2]](i)-t_node_diff_ksi[faces_nodes[3*ff+0]](i);
793 
794  const int node_shift = ii*4;
795 
796  beta_f[ff] =
797  N[node_shift+faces_nodes[3*ff+0]]*N[node_shift+faces_nodes[3*ff+1]]*
798  N[node_shift+faces_nodes[3*ff+2]];
799 
800  t_diff_beta_f[ff](j) =
801  t_node_diff_ksi[faces_nodes[3*ff+0]](j)*N[node_shift+faces_nodes[3*ff+1]]*N[node_shift+faces_nodes[3*ff+2]]+
802  N[node_shift+faces_nodes[3*ff+0]]*t_node_diff_ksi[faces_nodes[3*ff+1]](j)*N[node_shift+faces_nodes[3*ff+2]]+
803  N[node_shift+faces_nodes[3*ff+0]]*N[node_shift+faces_nodes[3*ff+1]]*t_node_diff_ksi[faces_nodes[3*ff+2]](j);
804 
805  const double ksi_0i =
806  N[node_shift+faces_nodes[3*ff+1]]-N[node_shift+faces_nodes[3*ff+0]];
807  ierr = base_polynomials(
808  p,ksi_0i,&t_diff_ksi0i(0),psi_l_0i[ff],diff_psi_l_0i[ff],3
809  ); CHKERRQ(ierr);
810 
811  const double ksi_0j =
812  N[node_shift+faces_nodes[3*ff+2]]-N[node_shift+faces_nodes[3*ff+0]];
813  ierr = base_polynomials(
814  p,ksi_0j,&t_diff_ksi0j(0),psi_l_0j[ff],diff_psi_l_0j[ff],3
815  ); CHKERRQ(ierr);
816 
817  }
818 
819  int cc = 0;
820  for(int oo = 0;oo<=(p-3);oo++) {
821  FTensor::Tensor1<double*,3> t_diff_psi_l_0i[] = {
822  FTensor::Tensor1<double*,3>(&diff_psi_l_0i[0][0],&diff_psi_l_0i[0][p+1],&diff_psi_l_0i[0][2*p+2],1),
823  FTensor::Tensor1<double*,3>(&diff_psi_l_0i[1][0],&diff_psi_l_0i[1][p+1],&diff_psi_l_0i[1][2*p+2],1),
824  FTensor::Tensor1<double*,3>(&diff_psi_l_0i[2][0],&diff_psi_l_0i[2][p+1],&diff_psi_l_0i[2][2*p+2],1),
825  FTensor::Tensor1<double*,3>(&diff_psi_l_0i[3][0],&diff_psi_l_0i[3][p+1],&diff_psi_l_0i[3][2*p+2],1),
826  };
827  for(int pp0 = 0;pp0<=oo;pp0++) {
828  const int pp1 = oo-pp0;
829  if(pp1>=0) {
830  for(int ff = 0;ff!=4;ff++) {
831  FTensor::Tensor1<double*,3> t_diff_psi_l_0j(
832  &m_diff_psi_l_0j(ff,pp1),
833  &m_diff_psi_l_0j(ff,p+1+pp1),
834  &m_diff_psi_l_0j(ff,2*p+2+pp1),1
835  );
836  const double t = psi_l_0i[ff][pp0]*psi_l_0j[ff][pp1];
837  const double a = beta_f[ff]*t;
838  t_phi_v(i) = a*t_node_diff_ksi[face_opposite_nodes[ff]](i);
839  ++t_phi_v;
840  ++cc;
841  t_diff_phi_v(i,j) = (
842  t_diff_beta_f[ff](j)*t+
843  beta_f[ff]*t_diff_psi_l_0i[ff](j)*psi_l_0j[ff][pp1]+
844  beta_f[ff]*psi_l_0i[ff][pp0]*t_diff_psi_l_0j(j)
845  )*t_node_diff_ksi[face_opposite_nodes[ff]](i);
846  ++t_diff_phi_v;
847  ++t_diff_psi_l_0i[ff];
848  }
849  }
850  }
851  }
852 
853  const int nb_base_fun_on_face = NBVOLUMETET_FACE_HCURL(p);
854  if(cc!=nb_base_fun_on_face) {
855  SETERRQ2(
856  PETSC_COMM_SELF,
858  "Wrong number of base functions %d != %d",
859  cc,nb_base_fun_on_face
860  );
861  }
862 
863  }
864 
865  PetscFunctionReturn(0);
866 }
ublas::matrix< double, ublas::row_major, DoubleAllacator > MatrixDouble
Definition: Common.hpp:136
CHKERRQ(ierr)
#define NBVOLUMETET_FACE_HCURL(P)

◆ Hcurl_VolumeFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_VolumeFunctions_MBTET ( int  p,
double N,
double diffN,
double phi_v,
double diff_phi_v,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

H-curl volume base functions.

Volume base functions are collection of iace interior base functions and volume interior base functions.

Parameters
papproximation order
Nnodal shape functions
diffNderivatives of nodal shape functions
phi_vcalculated shape functions
diff_phi_vderivatives of shape functions
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 1361 of file Hcurl.cpp.

1364  {
1365 
1366  PetscFunctionBegin;
1367 
1368  VectorDouble base_face_inetrior_functions(3*NBVOLUMETET_FACE_HCURL(p)*nb_integration_pts);
1369  VectorDouble diff_base_face_inetrior_functions(9*NBVOLUMETET_FACE_HCURL(p)*nb_integration_pts);
1370  // base_face_inetrior_functions.clear();
1371  // diff_base_face_inetrior_functions.clear();
1372  double *phi_v_f = &*base_face_inetrior_functions.data().begin();
1373  double *diff_phi_v_f = &*diff_base_face_inetrior_functions.data().begin();
1374  int faces_nodes[] = { 0,1,3, 1,2,3, 0,2,3, 0,1,2 };
1376  faces_nodes,p,N,diffN,phi_v_f,diff_phi_v_f,nb_integration_pts,base_polynomials
1377  ); CHKERRQ(ierr);
1378 
1379  VectorDouble base_interior_functions(3*NBVOLUMETET_TET_HCURL(p)*nb_integration_pts);
1380  VectorDouble diff_base_interior_functions(9*NBVOLUMETET_TET_HCURL(p)*nb_integration_pts);
1381  // base_interior_functions.clear();
1382  // diff_base_interior_functions.clear();
1383  double *phi_v_v = &*base_interior_functions.data().begin();
1384  double *diff_phi_v_v = &*diff_base_interior_functions.data().begin();
1386  p,N,diffN,phi_v_v,diff_phi_v_v,nb_integration_pts,base_polynomials
1387  ); CHKERRQ(ierr);
1388 
1391 
1392  FTensor::Tensor1<double*,3> t_face_interior(&phi_v_f[0],&phi_v_f[1],&phi_v_f[2],3);
1393  FTensor::Tensor2<double*,3,3> t_diff_face_interior(
1394  &diff_phi_v_f[0],&diff_phi_v_f[3],&diff_phi_v_f[6],
1395  &diff_phi_v_f[1],&diff_phi_v_f[4],&diff_phi_v_f[7],
1396  &diff_phi_v_f[2],&diff_phi_v_f[5],&diff_phi_v_f[8],9
1397  );
1398 
1399  FTensor::Tensor1<double*,3> t_phi_v(&phi_v[0],&phi_v[1],&phi_v[2],3);
1400  FTensor::Tensor2<double*,3,3> t_diff_phi_v(
1401  &diff_phi_v[0],&diff_phi_v[3],&diff_phi_v[6],
1402  &diff_phi_v[1],&diff_phi_v[4],&diff_phi_v[7],
1403  &diff_phi_v[2],&diff_phi_v[5],&diff_phi_v[8],9
1404  );
1405 
1406 
1407  if(NBVOLUMETET_TET_HCURL(p)>0) {
1408  FTensor::Tensor1<double*,3> t_volume_interior(&phi_v_v[0],&phi_v_v[1],&phi_v_v[2],3);
1409  FTensor::Tensor2<double*,3,3> t_diff_volume_interior(
1410  &diff_phi_v_v[0],&diff_phi_v_v[3],&diff_phi_v_v[6],
1411  &diff_phi_v_v[1],&diff_phi_v_v[4],&diff_phi_v_v[7],
1412  &diff_phi_v_v[2],&diff_phi_v_v[5],&diff_phi_v_v[8],9
1413  );
1414  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1415  int cc = 0;
1416  for(int oo = 0;oo<=p;oo++) {
1417  for(int ll = NBVOLUMETET_FACE_HCURL(oo-1);ll!=NBVOLUMETET_FACE_HCURL(oo);ll++) {
1418  t_phi_v(i) = t_face_interior(i);
1419  ++t_phi_v;
1420  ++t_face_interior;
1421  ++cc;
1422  t_diff_phi_v(i,j) = t_diff_face_interior(i,j);
1423  ++t_diff_phi_v;
1424  ++t_diff_face_interior;
1425  }
1426  for(int ll = NBVOLUMETET_TET_HCURL(oo-1);ll!=NBVOLUMETET_TET_HCURL(oo);ll++) {
1427  t_phi_v(i) = t_volume_interior(i);
1428  ++t_phi_v;
1429  ++t_volume_interior;
1430  ++cc;
1431  t_diff_phi_v(i,j) = t_diff_volume_interior(i,j);
1432  ++t_diff_phi_v;
1433  ++t_diff_volume_interior;
1434  }
1435  }
1436  // check consistency
1437  const int nb_base_fun_on_face = NBVOLUMETET_HCURL(p);
1438  if(cc!=nb_base_fun_on_face) {
1439  SETERRQ2(
1440  PETSC_COMM_SELF,
1442  "Wrong number of base functions %d != %d",
1443  cc,nb_base_fun_on_face
1444  );
1445  }
1446  }
1447  } else {
1448  for(int ii = 0;ii!=nb_integration_pts;ii++) {
1449  int cc = 0;
1450  for(int oo = 0;oo<=p;oo++) {
1451  for(int ll = NBVOLUMETET_FACE_HCURL(oo-1);ll!=NBVOLUMETET_FACE_HCURL(oo);ll++) {
1452  t_phi_v(i) = t_face_interior(i);
1453  ++t_phi_v;
1454  ++t_face_interior;
1455  ++cc;
1456  t_diff_phi_v(i,j) = t_diff_face_interior(i,j);
1457  ++t_diff_phi_v;
1458  ++t_diff_face_interior;
1459  }
1460  }
1461  // check consistency
1462  const int nb_base_fun_on_face = NBVOLUMETET_HCURL(p);
1463  if(cc!=nb_base_fun_on_face) {
1464  SETERRQ2(
1465  PETSC_COMM_SELF,
1467  "Wrong number of base functions %d != %d",
1468  cc,nb_base_fun_on_face
1469  );
1470  }
1471  }
1472  }
1473 
1474  PetscFunctionReturn(0);
1475 }
#define NBVOLUMETET_TET_HCURL(P)
PetscErrorCode Hcurl_VolumeInteriorFunctions_MBTET(int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Volume interior function.
Definition: Hcurl.cpp:868
CHKERRQ(ierr)
PetscErrorCode Hcurl_FaceInteriorFunctions_MBTET(int *faces_nodes, int p, double *N, double *diffN, double *phi_v, double *diff_phi_v, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face base interior function.
Definition: Hcurl.cpp:735
#define NBVOLUMETET_HCURL(P)
ublas::vector< double, DoubleAllacator > VectorDouble
Definition: Common.hpp:135
#define NBVOLUMETET_FACE_HCURL(P)

◆ Hcurl_VolumeInteriorFunctions_MBTET()

PetscErrorCode MoFEM::Hcurl_VolumeInteriorFunctions_MBTET ( int  p,
double N,
double diffN,
double phi_v,
double diff_phi_v,
int  nb_integration_pts,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Volume interior function.

On volume have (P-3)*(P-2)*(P-1)/2.

See NBVOLUMETET_TET_HCURL

Parameters
papproximation order
Nnodal shape functions
diffNderivatives of nodal shape functions
phi_vcalculated shape functions
diff_phi_vderivatives of shape functions
nb_integration_ptsnumber of shape functions
base_polynomialspolynomial base function (f.e. Legendre of Lobatto)
Returns
error code

Definition at line 868 of file Hcurl.cpp.

871  {
872 
873  PetscFunctionBegin;
874 
875  if(NBVOLUMETET_TET_HCURL(p)==0) PetscFunctionReturn(0);
876 
882 
883  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
884  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
885  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
886  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
887  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
888  };
889 
890  double diff_ksi0i[3],diff_ksi0j[3],diff_ksi0k[3];
891  FTensor::Tensor1<double*,3> t_diff_ksi0i(diff_ksi0i,&diff_ksi0i[1],&diff_ksi0i[2]);
892  FTensor::Tensor1<double*,3> t_diff_ksi0j(diff_ksi0j,&diff_ksi0j[1],&diff_ksi0j[2]);
893  FTensor::Tensor1<double*,3> t_diff_ksi0k(diff_ksi0k,&diff_ksi0k[1],&diff_ksi0k[2]);
894  t_diff_ksi0i(i) = t_node_diff_ksi[1](i)-t_node_diff_ksi[0](i);
895  t_diff_ksi0j(i) = t_node_diff_ksi[2](i)-t_node_diff_ksi[0](i);
896  t_diff_ksi0k(i) = t_node_diff_ksi[3](i)-t_node_diff_ksi[0](i);
897 
898  std::vector<double> v_psi_l_0i(p+1),v_diff_psi_l_0i(3*p+3);
899  std::vector<double> v_psi_l_0j(p+1),v_diff_psi_l_0j(3*p+3);
900  std::vector<double> v_psi_l_0k(p+1),v_diff_psi_l_0k(3*p+3);
901  double *psi_l_0i = &*v_psi_l_0i.begin();
902  double *diff_psi_l_0i = &*v_diff_psi_l_0i.begin();
903  double *psi_l_0j = &*v_psi_l_0j.begin();
904  double *diff_psi_l_0j = &*v_diff_psi_l_0j.begin();
905  double *psi_l_0k = &*v_psi_l_0k.begin();
906  double *diff_psi_l_0k = &*v_diff_psi_l_0k.begin();
907 
908  FTensor::Tensor1<double*,3> t_phi_v(&phi_v[0],&phi_v[1],&phi_v[2],3);
909  FTensor::Tensor2<double*,3,3> t_diff_phi_v(
910  &diff_phi_v[0],&diff_phi_v[3],&diff_phi_v[6],
911  &diff_phi_v[1],&diff_phi_v[4],&diff_phi_v[7],
912  &diff_phi_v[2],&diff_phi_v[5],&diff_phi_v[8],9
913  );
915 
916  for(int ii = 0;ii!=nb_integration_pts;ii++) {
917 
918  const int node_shift = ii*4;
919  const int n0 = node_shift+0;
920  const int n1 = node_shift+1;
921  const int n2 = node_shift+2;
922  const int n3 = node_shift+3;
923 
924  const double beta_v = N[n0]*N[n1]*N[n2]*N[n3];
925 
926  const double ksi_0i = N[n1]-N[n0];
927  ierr = base_polynomials(p,ksi_0i,diff_ksi0i,psi_l_0i,diff_psi_l_0i,3); CHKERRQ(ierr);
928  const double ksi_0j = N[n2]-N[n0];
929  ierr = base_polynomials(p,ksi_0j,diff_ksi0j,psi_l_0j,diff_psi_l_0j,3); CHKERRQ(ierr);
930  const double ksi_0k = N[n3]-N[n0];
931  ierr = base_polynomials(p,ksi_0k,diff_ksi0k,psi_l_0k,diff_psi_l_0k,3); CHKERRQ(ierr);
932 
933  FTensor::Tensor1<double,3> t_diff_beta_v;
934  t_diff_beta_v(j) =
935  t_node_diff_ksi[0](j)*N[n1]*N[n2]*N[n3]+
936  N[n0]*t_node_diff_ksi[1](j)*N[n2]*N[n3]+
937  N[n0]*N[n1]*t_node_diff_ksi[2](j)*N[n3]+
938  N[n0]*N[n1]*N[n2]*t_node_diff_ksi[3](j);
939 
940  int cc = 0;
941  for(int oo = 0;oo<=(p-4);oo++) {
942  FTensor::Tensor1<double*,3> t_diff_psi_l_0i(
943  &diff_psi_l_0i[0],
944  &diff_psi_l_0i[p+1],
945  &diff_psi_l_0i[2*p+2],1
946  );
947  for(int pp0 = 0;pp0<=oo;pp0++) {
948  FTensor::Tensor1<double*,3> t_diff_psi_l_0j(
949  &diff_psi_l_0j[0],
950  &diff_psi_l_0j[p+1],
951  &diff_psi_l_0j[2*p+2],1
952  );
953  for(int pp1 = 0;(pp0+pp1)<=oo;pp1++) {
954  const int pp2 = oo-pp0-pp1;
955  if(pp2>=0) {
956  FTensor::Tensor1<double*,3> t_diff_psi_l_0k(
957  &diff_psi_l_0k[0+pp2],
958  &diff_psi_l_0k[p+1+pp2],
959  &diff_psi_l_0k[2*p+2+pp2],1
960  );
961  const double t = psi_l_0i[pp0]*psi_l_0j[pp1]*psi_l_0k[pp2];
962  const double a = beta_v*t;
963  t_phi_v(0) = a;
964  t_phi_v(1) = 0;
965  t_phi_v(2) = 0;
966  ++t_phi_v;
967  ++cc;
968  t_phi_v(0) = 0;
969  t_phi_v(1) = a;
970  t_phi_v(2) = 0;
971  ++t_phi_v;
972  ++cc;
973  t_phi_v(0) = 0;
974  t_phi_v(1) = 0;
975  t_phi_v(2) = a;
976  ++t_phi_v;
977  ++cc;
978  t_b(j) =
979  t_diff_beta_v(j)*t+
980  beta_v*(
981  t_diff_psi_l_0i(j)*psi_l_0j[pp1]*psi_l_0k[pp2]+
982  psi_l_0i[pp0]*t_diff_psi_l_0j(j)*psi_l_0k[pp2]+
983  psi_l_0i[pp0]*psi_l_0j[pp1]*t_diff_psi_l_0k(j)
984  );
985  t_diff_phi_v(N0,j) = t_b(j);
986  t_diff_phi_v(N1,j) = 0;
987  t_diff_phi_v(N2,j) = 0;
988  ++t_diff_phi_v;
989  t_diff_phi_v(N0,j) = 0;
990  t_diff_phi_v(N1,j) = t_b(j);
991  t_diff_phi_v(N2,j) = 0;
992  ++t_diff_phi_v;
993  t_diff_phi_v(N0,j) = 0;
994  t_diff_phi_v(N1,j) = 0;
995  t_diff_phi_v(N2,j) = t_b(j);
996  ++t_diff_phi_v;
997  }
998  ++t_diff_psi_l_0j;
999  }
1000  ++t_diff_psi_l_0i;
1001  }
1002  }
1003 
1004  const int nb_base_fun_on_face = NBVOLUMETET_TET_HCURL(p);
1005  if(cc!=nb_base_fun_on_face) {
1006  SETERRQ2(
1007  PETSC_COMM_SELF,
1009  "Wrong number of base functions %d != %d",
1010  cc,nb_base_fun_on_face
1011  );
1012  }
1013 
1014  }
1015  PetscFunctionReturn(0);
1016 }
#define NBVOLUMETET_TET_HCURL(P)
CHKERRQ(ierr)

◆ Hdiv_Ainsworth_EdgeBasedVolumeShapeFunctions_MBTET()

PetscErrorCode MoFEM::Hdiv_Ainsworth_EdgeBasedVolumeShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double phi_v_e[6],
double diff_phi_v_e[6],
int  gdim,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Hdiv base function, Edge-based interior (volume) functions by Ainsworth [1].

Parameters
pvolume approx. order
Nshape functions
diffNderivatives of shape functions
phi_v_ebase functions (return)
diff_phi_v_ederivatives of base functions (returned)
gdimnumber of integration points
base_polynomialsbase function (Legendre/Lobbatto-Gauss)
Returns
error code

Definition at line 273 of file Hdiv.cpp.

283  {
284 
285  const int edges_nodes[6][2] = {
286  {0,1}, {1,2}, {2,0}, {0,3}, {1,3}, {2,3}
287  };
288 
289  PetscFunctionBegin;
290  if(p<2) PetscFunctionReturn(0);
291  if(diffN == NULL) {
292  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
293  }
294 
295  FTensor::Tensor1<double,3> t_coords[4] = {
300  };
301  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
302  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
303  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
304  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
305  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
306  };
307 
310 
312  FTensor::Tensor1<double,3> t_diff_ksi0i;
313  FTensor::Tensor1<double,3> t_diff_beta_e;
314 
315  double psi_l[p+1];
316  double diff_psi_l[3*(p+1)];
317 
318  for(int ee = 0;ee!=6;ee++) {
319  t_tou_e(i) =
320  t_coords[edges_nodes[ee][1]](i)-t_coords[edges_nodes[ee][0]](i);
321  t_diff_ksi0i(i) =
322  t_node_diff_ksi[edges_nodes[ee][1]](i)-t_node_diff_ksi[edges_nodes[ee][0]](i);
323  FTensor::Tensor1<double*,3> t_psi_v_e(
324  &phi_v_e[ee][0],&phi_v_e[ee][1],&phi_v_e[ee][2],3
325  );
326  FTensor::Tensor2<double*,3,3> t_diff_phi_v_e(
327  &diff_phi_v_e[ee][HDIV0_0],&diff_phi_v_e[ee][HDIV0_1],&diff_phi_v_e[ee][HDIV0_2],
328  &diff_phi_v_e[ee][HDIV1_0],&diff_phi_v_e[ee][HDIV1_1],&diff_phi_v_e[ee][HDIV1_2],
329  &diff_phi_v_e[ee][HDIV2_0],&diff_phi_v_e[ee][HDIV2_1],&diff_phi_v_e[ee][HDIV2_2],9
330  );
331  for(int ii = 0;ii!=gdim;ii++) {
332  const int node_shift = ii*4;
333  const double ni = N[node_shift+edges_nodes[ee][1]];
334  const double n0 = N[node_shift+edges_nodes[ee][0]];
335  const double beta_e = ni*n0;
336  const double ksi0i = ni-n0;
337  if(diff_phi_v_e) {
338  t_diff_beta_e(i) =
339  ni*t_node_diff_ksi[edges_nodes[ee][0]](i)+
340  t_node_diff_ksi[edges_nodes[ee][1]](i)*n0;
341  ierr = base_polynomials(p,ksi0i,&t_diff_ksi0i(0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
342  } else {
343  ierr = base_polynomials(p,ksi0i,NULL,psi_l,NULL,3); CHKERRQ(ierr);
344  }
345  FTensor::Tensor0<double*> t_psi_l(&psi_l[0]);
346  FTensor::Tensor1<double*,3> t_diff_psi_l(
347  &diff_psi_l[0],&diff_psi_l[p+1],&diff_psi_l[2*p+2],1
348  );
349  for(int l = 0;l<=p-2;l++) {
350  t_psi_v_e(i) = (beta_e*t_psi_l)*t_tou_e(i);
351  ++t_psi_v_e;
352  if(diff_phi_v_e) {
353  t_diff_phi_v_e(i,j) =
354  (t_diff_beta_e(j)*t_psi_l+beta_e*t_diff_psi_l(j))*
355  t_tou_e(i);
356  ++t_diff_phi_v_e;
357  ++t_diff_psi_l;
358  }
359  ++t_psi_l;
360  }
361  }
362  }
363 
364  PetscFunctionReturn(0);
365 }
CHKERRQ(ierr)

◆ Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET()

PetscErrorCode MoFEM::Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET ( int faces_nodes,
int p,
double N,
double diffN,
double phi_f_e[4][3],
double diff_phi_f_e[4][3],
int  gdim,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Hdiv base functions, Edge-based face functions by Ainsworth [1].

Parameters
faces_nodesFace nodes on tetrahedral
pApproximation order on faces
NShape functions
diffNDerivatives of shape functions
phi_f_eBase functions (returned)
diff_phi_f_eDerivatives of base functions (returned)
gdimNumber of integration pts
base_polynomialsbase function (Legendre/Lobbatto-Gauss)
Returns
error code

Definition at line 19 of file Hdiv.cpp.

22  {
23 
24  PetscFunctionBegin;
25  for(int ff = 0;ff<4;ff++) {
26  if(diff_phi_f_e!=NULL) {
28  &faces_nodes[3*ff],p[ff],N,diffN,phi_f_e[ff],diff_phi_f_e[ff],gdim,4,base_polynomials
29  ); CHKERRQ(ierr);
30  } else {
32  &faces_nodes[3*ff],p[ff],N,diffN,phi_f_e[ff],NULL,gdim,4,base_polynomials
33  ); CHKERRQ(ierr);
34  }
35  }
36  PetscFunctionReturn(0);
37 }
CHKERRQ(ierr)
PetscErrorCode Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Hdiv base functions, Edge-based face functions by Ainsworth .
Definition: Hdiv.cpp:39

◆ Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f_e[3],
double diff_phi_f_e[3],
int  gdim,
int  nb,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Hdiv base functions, Edge-based face functions by Ainsworth [1].

Parameters
faces_nodesFace nodes on face
pApproximation order on faces
NShape functions
diffNDerivatives of shape functions
phi_f_eBase functions (returned)
diff_phi_f_eDerivatives of base functions (returned)
gdimNumber of integration pts
base_polynomialsbase function (Legendre/Lobbatto-Gauss)
nbNumber of nodes on entity (4 if tet, 3 if triangle)
Returns
error code

Definition at line 39 of file Hdiv.cpp.

42  {
43 
44  const int face_edges_nodes[3][2] = { {0,1}, {1,2}, {2,0} };
45  const int face_oposite_edges_node[] = { 2, 0, 1 };
48 
49  PetscFunctionBegin;
50  if(p<1) PetscFunctionReturn(0);
51 
52  FTensor::Tensor1<double,3> t_edge_cross[3];
53  FTensor::Tensor1<double,3> t_node_diff_ksi[4];
54  FTensor::Tensor1<double,3> t_diff_ksi0i[3];
55  if(diffN!=NULL) {
56  t_node_diff_ksi[0] = FTensor::Tensor1<double,3>(diffN[0],diffN[ 1],diffN[ 2]);
57  t_node_diff_ksi[1] = FTensor::Tensor1<double,3>(diffN[3],diffN[ 4],diffN[ 5]);
58  t_node_diff_ksi[2] = FTensor::Tensor1<double,3>(diffN[6],diffN[ 7],diffN[ 8]);
59  t_node_diff_ksi[3] = FTensor::Tensor1<double,3>(diffN[9],diffN[10],diffN[11]);
60  for(int ee = 0;ee<3;ee++) {
61  const int n0 = faces_nodes[face_edges_nodes[ee][0]];
62  const int n1 = faces_nodes[face_edges_nodes[ee][1]];
63  t_diff_ksi0i[ee](i) =
64  t_node_diff_ksi[n1](i)-t_node_diff_ksi[n0](i);
65  t_edge_cross[ee](0) =
66  t_node_diff_ksi[n0](1)*t_node_diff_ksi[n1](2)-
67  t_node_diff_ksi[n0](2)*t_node_diff_ksi[n1](1);
68  t_edge_cross[ee](1) =
69  t_node_diff_ksi[n0](2)*t_node_diff_ksi[n1](0)-
70  t_node_diff_ksi[n0](0)*t_node_diff_ksi[n1](2);
71  t_edge_cross[ee](2) =
72  t_node_diff_ksi[n0](0)*t_node_diff_ksi[n1](1)-
73  t_node_diff_ksi[n0](1)*t_node_diff_ksi[n1](0);
74  }
75  } else {
76  for(int ee = 0;ee<3;ee++) {
77  t_edge_cross[ee](0) = 1;
78  t_edge_cross[ee](1) = 0;
79  t_edge_cross[ee](2) = 0;
80  }
81  }
82  double psi_l[p+1],diff_psi_l[3*(p+1)];
83  boost::shared_ptr<FTensor::Tensor2<double*,3,3> > t_diff_phi_f_e_ptr;
84 
85  for(int ee = 0;ee!=3;ee++) {
86  const int i0 = faces_nodes[face_edges_nodes[ee][0]];
87  const int i1 = faces_nodes[face_edges_nodes[ee][1]];
88  const int iO = faces_nodes[face_oposite_edges_node[ee]];
90  &phi_f_e[ee][0],&phi_f_e[ee][1],&phi_f_e[ee][2],3
91  );
92  if(diff_phi_f_e) {
93  t_diff_phi_f_e_ptr = boost::shared_ptr<FTensor::Tensor2<double*,3,3> >(
95  &diff_phi_f_e[ee][HDIV0_0],&diff_phi_f_e[ee][HDIV0_1],&diff_phi_f_e[ee][HDIV0_2],
96  &diff_phi_f_e[ee][HDIV1_0],&diff_phi_f_e[ee][HDIV1_1],&diff_phi_f_e[ee][HDIV1_2],
97  &diff_phi_f_e[ee][HDIV2_0],&diff_phi_f_e[ee][HDIV2_1],&diff_phi_f_e[ee][HDIV2_2],9
98  )
99  );
100  }
101  for(int ii = 0;ii!=gdim;ii++) {
102  const int node_shift = ii*nb;
103  const double n0 = N[node_shift+i0];
104  const double n1 = N[node_shift+i1];
105  const double lambda = N[node_shift+iO];
106  const double ksi0i = n1 - n0;
107  if(diff_phi_f_e) {
108  ierr = base_polynomials(p,ksi0i,&t_diff_ksi0i[ee](0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
109  } else {
110  ierr = base_polynomials(p,ksi0i,NULL,psi_l,NULL,3); CHKERRQ(ierr);
111  }
112  FTensor::Tensor0<double*> t_psi_l(&psi_l[0]);
113  FTensor::Tensor1<double*,3> t_diff_psi_l(
114  &diff_psi_l[0],&diff_psi_l[p+1],&diff_psi_l[2*p+2],1
115  );
116  for(int l = 0;l<=p-1;l++) {
117  t_psi_f_e(i) = lambda*t_psi_l*t_edge_cross[ee](i);
118  if(t_diff_phi_f_e_ptr) {
119  (*t_diff_phi_f_e_ptr)(i,j) =
120  (t_node_diff_ksi[iO](j)*t_psi_l+lambda*t_diff_psi_l(j))*t_edge_cross[ee](i);
121  ++t_diff_psi_l;
122  ++(*t_diff_phi_f_e_ptr);
123  }
124  ++t_psi_f_e;
125  ++t_psi_l;
126  }
127  }
128  }
129  PetscFunctionReturn(0);
130 }
CHKERRQ(ierr)

◆ Hdiv_Ainsworth_FaceBasedVolumeShapeFunctions_MBTET()

PetscErrorCode MoFEM::Hdiv_Ainsworth_FaceBasedVolumeShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double phi_v_f[],
double diff_phi_v_f[],
int  gdim,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Hdiv Face-based interior functions by Ainsworth [1]

Parameters
pApproximation order on face
NShape functions on face
diffNDerivatives of shape functions of face
phi_v_fBase functions (returned)
diff_phi_v_fDerivatives of base functions (returned)
gdimNumber of integration points
base_polynomialsBase function (Legendre/Lobbatto-Gauss)
Returns
Error code

Definition at line 367 of file Hdiv.cpp.

370  {
371 
372  const int faces_nodes[4][3] = { {0,1,3}, {1,2,3}, {0,2,3}, {0,1,2} };
373 
374  PetscFunctionBegin;
375  if(p<3) PetscFunctionReturn(0);
376 
377  FTensor::Tensor1<double,3> t_coords[4] = {
382  };
383 
384  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
385  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
386  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
387  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
388  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
389  };
390 
393 
394  FTensor::Tensor1<double,3> t_tau0i[4],t_tau0j[4];
395  FTensor::Tensor1<double,3> t_diff_ksi0i[4],t_diff_ksi0j[4];
396  for(int ff = 0;ff!=4;ff++) {
397  const int v0 = faces_nodes[ff][0];
398  const int vi = faces_nodes[ff][1];
399  const int vj = faces_nodes[ff][2];
400  t_tau0i[ff](i) = t_coords[vi](i)-t_coords[v0](i);
401  t_tau0j[ff](i) = t_coords[vj](i)-t_coords[v0](i);
402  t_diff_ksi0i[ff](i) = t_node_diff_ksi[vi](i)-t_node_diff_ksi[v0](i);
403  t_diff_ksi0j[ff](i) = t_node_diff_ksi[vj](i)-t_node_diff_ksi[v0](i);
404  }
405 
406  double psi_l[p+1],psi_m[p+1];
407  double diff_psi_l[3*(p+1)],diff_psi_m[3*(p+1)];
408  for(int ff = 0;ff!=4;ff++) {
409  const int v0 = faces_nodes[ff][0];
410  const int vi = faces_nodes[ff][1];
411  const int vj = faces_nodes[ff][2];
412  FTensor::Tensor1<double*,3> t_phi_v_f(
413  &phi_v_f[ff][HDIV0],&phi_v_f[ff][HDIV1],&phi_v_f[ff][HDIV2],3
414  );
415  FTensor::Tensor2<double*,3,3> t_diff_phi_v_f(
416  &diff_phi_v_f[ff][HDIV0_0],&diff_phi_v_f[ff][HDIV0_1],&diff_phi_v_f[ff][HDIV0_2],
417  &diff_phi_v_f[ff][HDIV1_0],&diff_phi_v_f[ff][HDIV1_1],&diff_phi_v_f[ff][HDIV1_2],
418  &diff_phi_v_f[ff][HDIV2_0],&diff_phi_v_f[ff][HDIV2_1],&diff_phi_v_f[ff][HDIV2_2],9
419  );
420  for(int ii = 0;ii<gdim;ii++) {
421  const int node_shift = 4*ii;
422  const double n0 = N[node_shift+v0];
423  const double ni = N[node_shift+vi];
424  const double nj = N[node_shift+vj];
425  const double beta_f = n0*ni*nj;
426  FTensor::Tensor1<double,3> t_diff_beta_f;
427  t_diff_beta_f(i) =
428  (ni*nj)*t_node_diff_ksi[v0](i)+
429  (n0*nj)*t_node_diff_ksi[vi](i)+
430  (n0*ni)*t_node_diff_ksi[vj](i);
431  const double ksi0i = ni-n0;
432  const double ksi0j = nj-n0;
433  ierr = base_polynomials(p,ksi0i,&t_diff_ksi0i[ff](0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
434  ierr = base_polynomials(p,ksi0j,&t_diff_ksi0j[ff](0),psi_m,diff_psi_m,3); CHKERRQ(ierr);
436  int jj = 0;
437  for(int oo = 0;oo<=p-3;oo++) {
438  FTensor::Tensor0<double*> t_psi_l(&psi_l[0],1);
439  FTensor::Tensor1<double*,3> t_diff_psi_l(
440  diff_psi_l,&diff_psi_l[p+1],&diff_psi_l[2*p+2],1
441  );
442  for(int l = 0;l<=oo;l++) {
443  int m = oo - l;
444  if(m>=0) {
445 
446  FTensor::Tensor1<double*,3> t_diff_psi_m(
447  &diff_psi_m[m],&diff_psi_m[p+1+m],&diff_psi_m[2*p+2+m],1
448  );
449  const double a = beta_f*t_psi_l*psi_m[m];
450  t_phi_v_f(i) = a*t_tau0i[ff](i);
451  ++t_phi_v_f;
452  ++jj;
453  t_phi_v_f(i) = a*t_tau0j[ff](i);
454  ++t_phi_v_f;
455  ++jj;
456 
457  t_diff_a(j) =
458  (t_psi_l*psi_m[m])*t_diff_beta_f(j)+
459  (beta_f*psi_m[m])*t_diff_psi_l(j)+
460  (beta_f*t_psi_l)*t_diff_psi_m(j);
461  t_diff_phi_v_f(i,j) = t_diff_a(j)*t_tau0i[ff](i);
462  ++t_diff_phi_v_f;
463  t_diff_phi_v_f(i,j) = t_diff_a(j)*t_tau0j[ff](i);
464  ++t_diff_phi_v_f;
465 
466  ++t_psi_l;
467  ++t_diff_psi_l;
468 
469  }
470  }
471  }
473  SETERRQ2(
474  PETSC_COMM_SELF,
476  "wrong order %d != %d",
478  }
479  }
480  }
481  PetscFunctionReturn(0);
482 }
#define NBVOLUMETET_AINSOWRTH_FACE_HDIV(P)
CHKERRQ(ierr)

◆ Hdiv_Ainsworth_FaceBubbleShapeFunctions()

PetscErrorCode MoFEM::Hdiv_Ainsworth_FaceBubbleShapeFunctions ( int faces_nodes,
int p,
double N,
double diffN,
double phi_f[],
double diff_phi_f[],
int  gdim,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face bubble functions by Ainsworth [1].

Parameters
faces_nodesFace nodes on tetrahedral
pApprox. order on faces
NShape function
diffNDerivatives of shape functions
phi_fBase functions
diff_phi_fDerivatives of base functions
gdimNumber of integration pts
base_polynomialsBase function (Legendre/Lobbatto-Gauss)
Returns
error code

Definition at line 132 of file Hdiv.cpp.

135  {
136 
137  PetscFunctionBegin;
138  for(int ff = 0;ff<4;ff++) {
139  double *diff;
140  if(diff_phi_f!=NULL) {
141  diff = diff_phi_f[ff];
142  } else {
143  diff = NULL;
144  }
146  &faces_nodes[3*ff],p[ff],N,diffN,phi_f[ff],diff,gdim,4,base_polynomials
147  ); CHKERRQ(ierr);
148  }
149  PetscFunctionReturn(0);
150 }
PetscErrorCode Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face bubble functions by Ainsworth .
Definition: Hdiv.cpp:152
CHKERRQ(ierr)

◆ Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE()

PetscErrorCode MoFEM::Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f,
double diff_phi_f,
int  gdim,
int  nb,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Face bubble functions by Ainsworth [1].

Parameters
faces_nodesFace nodes on face
pApprox. order on face
NShape function
diffNDerivatives of shape functions
phi_fBase functions
diff_phi_fDerivatives of base functions
gdimNumber of integration pts
nbNumber of nodes on entity (4 if tet, 3 if triangle)
base_polynomialsBase function (Legendre/Lobbatto-Gauss)
Returns
error code

Definition at line 152 of file Hdiv.cpp.

155  {
158 
159  PetscFunctionBegin;
160  if(p<3) PetscFunctionReturn(0);
161 
162  const int vert_i = face_nodes[1];
163  const int vert_j = face_nodes[2];
164  const int i0 = face_nodes[0];
166  FTensor::Tensor1<double,3> t_node_diff_ksi[4];
167  FTensor::Tensor1<double,3> t_diff_ksi0i;
168  FTensor::Tensor1<double,3> t_diff_ksi0j;
169 
170  if(diffN) {
171  t_node_diff_ksi[0] = FTensor::Tensor1<double,3>(diffN[0],diffN[ 1],diffN[ 2]);
172  t_node_diff_ksi[1] = FTensor::Tensor1<double,3>(diffN[3],diffN[ 4],diffN[ 5]);
173  t_node_diff_ksi[2] = FTensor::Tensor1<double,3>(diffN[6],diffN[ 7],diffN[ 8]);
174  t_node_diff_ksi[3] = FTensor::Tensor1<double,3>(diffN[9],diffN[10],diffN[11]);
175  t_diff_ksi0i(i) =
176  t_node_diff_ksi[vert_i](i)-t_node_diff_ksi[i0](i);
177  t_diff_ksi0j(i) =
178  t_node_diff_ksi[vert_j](i)-t_node_diff_ksi[i0](i);
179  t_cross(0) =
180  t_node_diff_ksi[vert_i](1)*t_node_diff_ksi[vert_j](2)-
181  t_node_diff_ksi[vert_i](2)*t_node_diff_ksi[vert_j](1);
182  t_cross(1) =
183  t_node_diff_ksi[vert_i](2)*t_node_diff_ksi[vert_j](0)-
184  t_node_diff_ksi[vert_i](0)*t_node_diff_ksi[vert_j](2);
185  t_cross(2) =
186  t_node_diff_ksi[vert_i](0)*t_node_diff_ksi[vert_j](1)-
187  t_node_diff_ksi[vert_i](1)*t_node_diff_ksi[vert_j](0);
188  } else {
189  t_cross(0) = 1;
190  t_cross(1) = 0;
191  t_cross(2) = 0;
192  }
193 
194  double psi_l[p+1],diff_psi_l[3*(p+1)];
195  double psi_m[p+1],diff_psi_m[3*(p+1)];
196  FTensor::Tensor1<double,3> t_diff_beta_0ij;
197 
199  &phi_f[HDIV0],&phi_f[HDIV1],&phi_f[HDIV2],3
200  );
201 
202  boost::shared_ptr<FTensor::Tensor2<double*,3,3> > t_diff_phi_f_ptr;
203  if(diff_phi_f) {
204  t_diff_phi_f_ptr = boost::shared_ptr<FTensor::Tensor2<double*,3,3> >(
206  &diff_phi_f[HDIV0_0],&diff_phi_f[HDIV0_1],&diff_phi_f[HDIV0_2],
207  &diff_phi_f[HDIV1_0],&diff_phi_f[HDIV1_1],&diff_phi_f[HDIV1_2],
208  &diff_phi_f[HDIV2_0],&diff_phi_f[HDIV2_1],&diff_phi_f[HDIV2_2],9
209  )
210  );
211  }
212 
213  for(int ii = 0;ii<gdim;ii++) {
214 
215  int node_shift = ii*nb;
216  const double ni = N[ node_shift+vert_i ];
217  const double nj = N[ node_shift+vert_j ];
218  const double n0 = N[ node_shift+i0 ];
219  const double ksi0i = ni - n0;
220  const double ksi0j = nj - n0;
221  double beta_0ij = n0*ni*nj;
222  if(diff_phi_f) {
223  t_diff_beta_0ij(i) =
224  (ni*nj)*t_node_diff_ksi[i0](i)+
225  (n0*nj)*t_node_diff_ksi[vert_i](i)+
226  (n0*ni)*t_node_diff_ksi[vert_j](i);
227  ierr = base_polynomials(p,ksi0i,&t_diff_ksi0i(0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
228  ierr = base_polynomials(p,ksi0j,&t_diff_ksi0j(0),psi_m,diff_psi_m,3); CHKERRQ(ierr);
229  } else {
230  ierr = base_polynomials(p,ksi0i,NULL,psi_l,NULL,3); CHKERRQ(ierr);
231  ierr = base_polynomials(p,ksi0j,NULL,psi_m,NULL,3); CHKERRQ(ierr);
232  }
233 
234  int jj = 0;
235  int oo = 0;
236  for(;oo<=p-3;oo++) {
237  FTensor::Tensor0<double*> t_psi_l(&psi_l[0]);
238  FTensor::Tensor1<double*,3> t_diff_psi_l(
239  diff_psi_l,&diff_psi_l[p+1],&diff_psi_l[2*p+2],1
240  );
241  for(int l = 0;l<=oo;l++) {
242  int m = oo - l;
243  if(m>=0) {
244  FTensor::Tensor1<double,3> t_diff_psi_m(
245  diff_psi_m[m],diff_psi_m[p+1+m],diff_psi_m[2*p+2+m]
246  );
247  t_psi_f(i) = (beta_0ij*t_psi_l*psi_m[m])*t_cross(i);
248  ++t_psi_f;
249  if(diff_phi_f) {
250  (*t_diff_phi_f_ptr)(i,j) = (
251  (t_psi_l*psi_m[m])*t_diff_beta_0ij(j)+
252  (beta_0ij*psi_m[m])*t_diff_psi_l(j)+
253  (beta_0ij*t_psi_l)*t_diff_psi_m(j)
254  )*t_cross(i);
255  ++(*t_diff_phi_f_ptr);
256  }
257  }
258  ++t_psi_l;
259  ++t_diff_psi_l;
260  ++jj;
261  }
262  }
263  if(jj!=NBFACETRI_AINSWORTH_FACE_HDIV(p)) {
264  SETERRQ2(
265  PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"wrong order %d != %d",
267  );
268  }
269  }
270  PetscFunctionReturn(0);
271 }
CHKERRQ(ierr)
#define NBFACETRI_AINSWORTH_FACE_HDIV(P)

◆ Hdiv_Ainsworth_VolumeBubbleShapeFunctions_MBTET()

PetscErrorCode MoFEM::Hdiv_Ainsworth_VolumeBubbleShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double phi_v,
double diff_phi_v,
int  gdim,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Interior bubble functions by Ainsworth [1].

Parameters
pVolume order
NShape functions
diffNDerivatives of shape functions
phi_vBase functions
diff_phi_vDerivatives of shape functions
gdimNumber of integration points
base_polynomialsBase function (Legendre/Lobbatto-Gauss)
Returns
Error code

Definition at line 484 of file Hdiv.cpp.

487  {
488 
489  PetscFunctionBegin;
490  if(p<4) PetscFunctionReturn(0);
491 
492  FTensor::Tensor1<double*,3> t_node_diff_ksi[4] = {
493  FTensor::Tensor1<double*,3>(&diffN[0],&diffN[ 1],&diffN[ 2]),
494  FTensor::Tensor1<double*,3>(&diffN[3],&diffN[ 4],&diffN[ 5]),
495  FTensor::Tensor1<double*,3>(&diffN[6],&diffN[ 7],&diffN[ 8]),
496  FTensor::Tensor1<double*,3>(&diffN[9],&diffN[10],&diffN[11])
497  };
498 
504 
505  FTensor::Tensor1<double,3> t_diff_ksi0i;
506  FTensor::Tensor1<double,3> t_diff_ksi0j;
507  FTensor::Tensor1<double,3> t_diff_ksi0k;
508 
509  t_diff_ksi0i(i) = t_node_diff_ksi[1](i)-t_node_diff_ksi[0](i);
510  t_diff_ksi0j(i) = t_node_diff_ksi[2](i)-t_node_diff_ksi[0](i);
511  t_diff_ksi0k(i) = t_node_diff_ksi[3](i)-t_node_diff_ksi[0](i);
512 
513  double psi_l[p+1];
514  double diff_psi_l[3*(p+1)];
515  double psi_m[p+1];
516  double diff_psi_m[3*(p+1)];
517  double psi_n[p+1];
518  double diff_psi_n[3*(p+1)];
519 
520  FTensor::Tensor1<double*,3> t_phi_v(phi_v,&phi_v[HDIV1],&phi_v[HDIV2],3);
521  FTensor::Tensor2<double*,3,3> t_diff_phi_v(
522  diff_phi_v,&diff_phi_v[HDIV0_1],&diff_phi_v[HDIV0_2],
523  &diff_phi_v[HDIV1_0],&diff_phi_v[HDIV1_1],&diff_phi_v[HDIV1_2],
524  &diff_phi_v[HDIV2_0],&diff_phi_v[HDIV2_1],&diff_phi_v[HDIV2_2],9
525  );
526 
527  FTensor::Tensor1<double,3> t_diff_beta_v;
528  for(int ii = 0;ii<gdim;ii++) {
529  const int node_shift = ii*4;
530  const double n0 = N[node_shift+0];
531  const double ni = N[node_shift+1];
532  const double nj = N[node_shift+2];
533  const double nk = N[node_shift+3];
534  const double ksi0i = ni-n0;
535  const double ksi0j = nj-n0;
536  const double ksi0k = nk-n0;
537  const double beta_v = n0*ni*nj*nk;
538  t_diff_beta_v(i) =
539  (ni*nj*nk)*t_node_diff_ksi[0](i)+
540  (n0*nj*nk)*t_node_diff_ksi[1](i)+
541  (n0*ni*nk)*t_node_diff_ksi[2](i)+
542  (n0*ni*nj)*t_node_diff_ksi[3](i);
543  ierr = base_polynomials(p,ksi0i,&t_diff_ksi0i(0),psi_l,diff_psi_l,3); CHKERRQ(ierr);
544  ierr = base_polynomials(p,ksi0j,&t_diff_ksi0j(0),psi_m,diff_psi_m,3); CHKERRQ(ierr);
545  ierr = base_polynomials(p,ksi0k,&t_diff_ksi0k(0),psi_n,diff_psi_n,3); CHKERRQ(ierr);
546 
548 
549  int jj = 0;
550  for(int oo = 0;oo<=p-4;oo++) {
551  FTensor::Tensor0<double*> t_psi_l(&psi_l[0]);
552  FTensor::Tensor1<double*,3> t_diff_psi_l(
553  diff_psi_l,&diff_psi_l[p+1],&diff_psi_l[2*p+2],1
554  );
555  for(int l = 0;l<=oo;l++) {
556  FTensor::Tensor0<double*> t_psi_m(&psi_m[0]);
557  FTensor::Tensor1<double*,3> t_diff_psi_m(
558  diff_psi_m,&diff_psi_m[p+1],&diff_psi_m[2*p+2],1
559  );
560  for(int m = 0;(l+m)<=oo;m++) {
561  int n = oo - l - m;
562  if(n>=0) {
563  FTensor::Tensor1<double,3> t_diff_psi_n(
564  diff_psi_n[n],diff_psi_n[p+1+n],diff_psi_n[2*p+2+n]
565  );
566  const double a = beta_v*t_psi_l*t_psi_m*psi_n[n];
567  t_phi_v(0) = a;
568  t_phi_v(1) = 0;
569  t_phi_v(2) = 0;
570  ++t_phi_v;
571  t_phi_v(0) = 0;
572  t_phi_v(1) = a;
573  t_phi_v(2) = 0;
574  ++t_phi_v;
575  t_phi_v(0) = 0;
576  t_phi_v(1) = 0;
577  t_phi_v(2) = a;
578  ++t_phi_v;
579  t_diff_a(j) =
580  (t_psi_l*t_psi_m*psi_n[n])*t_diff_beta_v(j)+
581  (beta_v*t_psi_m*psi_n[n])*t_diff_psi_l(j)+
582  (beta_v*t_psi_l*psi_n[n])*t_diff_psi_m(j)+
583  (beta_v*t_psi_l*t_psi_m)*t_diff_psi_n(j);
584  t_diff_phi_v(N0,j) = t_diff_a(j);
585  t_diff_phi_v(N1,j) = 0;
586  t_diff_phi_v(N2,j) = 0;
587  ++t_diff_phi_v;
588  t_diff_phi_v(N0,j) = 0;
589  t_diff_phi_v(N1,j) = t_diff_a(j);
590  t_diff_phi_v(N2,j) = 0;
591  ++t_diff_phi_v;
592  t_diff_phi_v(N0,j) = 0;
593  t_diff_phi_v(N1,j) = 0;
594  t_diff_phi_v(N2,j) = t_diff_a(j);
595  ++t_diff_phi_v;
596  ++jj;
597  }
598  ++t_psi_m;
599  ++t_diff_psi_m;
600  }
601  ++t_psi_l;
602  ++t_diff_psi_l;
603  }
604  }
605 
606  if(3*jj!=NBVOLUMETET_AINSWORTH_VOLUME_HDIV(p)) {
607  SETERRQ2(
608  PETSC_COMM_SELF,
610  "wrong order %d != %d",jj,
612  );
613  }
614  }
615 
616  PetscFunctionReturn(0);
617 }
#define NBVOLUMETET_AINSWORTH_VOLUME_HDIV(P)
CHKERRQ(ierr)

◆ Hdiv_Demkowicz_Face_MBTET_ON_FACE()

PetscErrorCode MoFEM::Hdiv_Demkowicz_Face_MBTET_ON_FACE ( int faces_nodes,
int  p,
double N,
double diffN,
double phi_f,
double diff_phi_f,
int  gdim,
int  nb 
)

HDiv base finuctions on triangle by Demkowicz [22]

Parameters
faces_nodesNodes on face
pApprox. order
NShape functions
diffNDerivatives of shape functions
phi_fReturned base functions on face
diff_phi_fReturned derivatives of base functions on face
gdimNumber of integration points
nbnb is 4 for face on tetrahedral and 3 for face
Returns
error code

Definition at line 619 of file Hdiv.cpp.

625  {
626  const int face_edges_nodes[3][2] = { {0,1}, {1,2}, {2,0} };
627  const int face_oposite_edges_node[] = { 2, 0, 1 };
628 
629  PetscFunctionBegin;
630 
633 
634  FTensor::Tensor1<double,3> t_cross[3];
635  FTensor::Tensor2<double,3,3> t_diff_cross;
636  FTensor::Tensor1<double,3> t_node_diff_ksi[4];
637  FTensor::Tensor1<double,3> t_node_diff_sum_n0_n1;
638 
639  const int i0 = faces_nodes[0];
640  const int i1 = faces_nodes[1];
641  const int i2 = faces_nodes[2];
642  const int o[] = {
643  faces_nodes[face_oposite_edges_node[0]],
644  faces_nodes[face_oposite_edges_node[1]],
645  faces_nodes[face_oposite_edges_node[2]]
646  };
647 
648  FTensor::Tensor1<double,3> t_diff_n0_p_n1;
649  FTensor::Tensor1<double,3> t_diff_n0_p_n1_p_n2;
650 
651  if(diff_phi_f) {
652  t_node_diff_ksi[0] = FTensor::Tensor1<double,3>(diffN[0],diffN[ 1],diffN[ 2]);
653  t_node_diff_ksi[1] = FTensor::Tensor1<double,3>(diffN[3],diffN[ 4],diffN[ 5]);
654  t_node_diff_ksi[2] = FTensor::Tensor1<double,3>(diffN[6],diffN[ 7],diffN[ 8]);
655  t_node_diff_ksi[3] = FTensor::Tensor1<double,3>(diffN[9],diffN[10],diffN[11]);
656  t_diff_cross(i,j) = 0;
657  for(int ee = 0;ee!=3;ee++) {
658  int ei0 = faces_nodes[face_edges_nodes[ee][0]];
659  int ei1 = faces_nodes[face_edges_nodes[ee][1]];
660  t_cross[ee](0) =
661  t_node_diff_ksi[ei0](1)*t_node_diff_ksi[ei1](2)-
662  t_node_diff_ksi[ei0](2)*t_node_diff_ksi[ei1](1);
663  t_cross[ee](1) =
664  t_node_diff_ksi[ei0](2)*t_node_diff_ksi[ei1](0)-
665  t_node_diff_ksi[ei0](0)*t_node_diff_ksi[ei1](2);
666  t_cross[ee](2) =
667  t_node_diff_ksi[ei0](0)*t_node_diff_ksi[ei1](1)-
668  t_node_diff_ksi[ei0](1)*t_node_diff_ksi[ei1](0);
670  diffN[3*o[ee]+0],diffN[3*o[ee]+1],diffN[3*o[ee]+2]
671  );
672  t_diff_cross(i,j) += t_cross[ee](i)*t_diff_o(j);
673  // cerr << t_cross[ee](0) << " " << t_cross[ee](1) << " " << t_cross[ee](2) << endl;
674  }
675  // cerr << endl << endl;
676  t_diff_n0_p_n1(i) = t_node_diff_ksi[i0](i)+t_node_diff_ksi[i1](i);
677  t_diff_n0_p_n1_p_n2(i) = t_diff_n0_p_n1(i)+t_node_diff_ksi[i2](i);
678  } else {
679  for(int ee = 0;ee!=3;ee++) {
680  t_cross[ee](0) = 1;
681  t_cross[ee](1) = 0;
682  t_cross[ee](2) = 0;
683  }
684  }
685 
686  FTensor::Tensor1<double*,3> t_phi(&phi_f[HDIV0],&phi_f[HDIV1],&phi_f[HDIV2],3);
687  boost::shared_ptr<FTensor::Tensor2<double*,3,3> > t_diff_phi_ptr;
688  if(diff_phi_f) {
689  t_diff_phi_ptr = boost::shared_ptr<FTensor::Tensor2<double*,3,3> >(
691  &diff_phi_f[HDIV0_0],&diff_phi_f[HDIV0_1],&diff_phi_f[HDIV0_2],
692  &diff_phi_f[HDIV1_0],&diff_phi_f[HDIV1_1],&diff_phi_f[HDIV1_2],
693  &diff_phi_f[HDIV2_0],&diff_phi_f[HDIV2_1],&diff_phi_f[HDIV2_2],9
694  )
695  );
696  }
697 
698  double fi[p+1],diff_fi[3*p+3];
699  double fj[p+1],diff_fj[3*p+3];
700  double tmp_fj[p+1],tmp_diff_fj[3*p+3];
701  for(int ii = 0;ii!=gdim;ii++) {
702  const int shift = ii*nb;
703  double n0 = N[shift+i0];
704  double n1 = N[shift+i1];
705  double n2 = N[shift+i2];
706  double *diff_n1 = (diff_phi_f) ? &t_node_diff_ksi[i1](0) : NULL;
707  double *diff_n0_p_n1 = (diff_phi_f) ? &t_diff_n0_p_n1(0) : NULL;
709  p,0,n1,n0+n1,diff_n1,diff_n0_p_n1,fi,diff_phi_f?diff_fi:NULL,3
710  ); CHKERRQ(ierr);
711  for(int pp = 0;pp<=p;pp++) {
712  double *diff_n2 = (diff_phi_f) ? &t_node_diff_ksi[i2](0) : NULL;
713  double *diff_n0_p_n1_p_n2 = (diff_phi_f) ? &t_diff_n0_p_n1_p_n2(0) : NULL;
715  pp,2*pp+1,n2,n0+n1+n2,diff_n2,diff_n0_p_n1_p_n2,tmp_fj,diff_phi_f?tmp_diff_fj:NULL,3
716  ); CHKERRQ(ierr);
717  fj[pp] = tmp_fj[pp];
718  if(diff_phi_f) {
719  diff_fj[0*(p+1)+pp] = tmp_diff_fj[0*(pp+1)+pp];
720  diff_fj[1*(p+1)+pp] = tmp_diff_fj[1*(pp+1)+pp];
721  diff_fj[2*(p+1)+pp] = tmp_diff_fj[2*(pp+1)+pp];
722  }
723  }
724  double no0 = N[shift+o[0]];
725  double no1 = N[shift+o[1]];
726  double no2 = N[shift+o[2]];
728  base0(i) = no0*t_cross[0](i)+no1*t_cross[1](i)+no2*t_cross[2](i);
729  int jj = 0;
730  for(int oo = 0;oo<p;oo++) {
731  FTensor::Tensor0<double*> t_fi(fi);
732  FTensor::Tensor1<double*,3> t_diff_fi(
733  &diff_fi[0],&diff_fi[p+1],&diff_fi[2*p+2],1
734  );
735  for(int ll = 0;ll<=oo;ll++) {
736  const int mm = oo-ll;
737  if(mm>=0) {
738  const double a = t_fi*fj[mm];
739  // cerr << ll << " " << mm << " " << a << endl;
740  t_phi(i) = a*base0(i);
741  if(diff_phi_f) {
742  FTensor::Tensor1<double*,3> t_diff_fj(
743  &diff_fj[0+mm],&diff_fj[p+1+mm],&diff_fj[2*p+2+mm],1
744  );
745  (*t_diff_phi_ptr)(i,j) =
746  a*t_diff_cross(i,j)+
747  (
748  t_diff_fi(j)*fj[mm]+t_fi*t_diff_fj(j)
749  )*base0(i);
750  ++(*t_diff_phi_ptr);
751  ++t_diff_fi;
752  }
753  ++t_fi;
754  ++t_phi;
755  ++jj;
756  }
757  }
758  }
759  if(jj!=NBFACETRI_DEMKOWICZ_HDIV(p)) {