v0.15.0
Loading...
Searching...
No Matches
PoissonExample::CreateFiniteElements Struct Reference

Create finite elements instances. More...

#include "tutorials/cor-2to5/src/PoissonOperators.hpp"

Collaboration diagram for PoissonExample::CreateFiniteElements:
[legend]

Public Member Functions

 CreateFiniteElements (MoFEM::Interface &m_field)
 
MoFEMErrorCode createFEToAssembleMatrixAndVector (boost::function< double(const double, const double, const double)> f_u, boost::function< double(const double, const double, const double)> f_source, boost::shared_ptr< ForcesAndSourcesCore > &domain_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &domain_rhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_rhs_fe, bool trans=true) const
 Create finite element to calculate matrix and vectors.
 
MoFEMErrorCode createFEToEvaluateError (boost::function< double(const double, const double, const double)> f_u, boost::function< FTensor::Tensor1< double, 3 >(const double, const double, const double)> g_u, Vec global_error, boost::shared_ptr< ForcesAndSourcesCore > &domain_error) const
 Create finite element to calculate error.
 
MoFEMErrorCode creatFEToPostProcessResults (boost::shared_ptr< PostProcFE > &post_proc_volume) const
 Create finite element to post-process results.
 
MoFEMErrorCode createFEToAssembleMatrixAndVectorForNonlinearProblem (boost::function< double(const double, const double, const double)> f_u, boost::function< double(const double, const double, const double)> f_source, boost::function< double(const double)> a, boost::function< double(const double)> diff_a, boost::shared_ptr< ForcesAndSourcesCore > &domain_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &domain_rhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_rhs_fe, ForcesAndSourcesCore::RuleHookFun vol_rule, ForcesAndSourcesCore::RuleHookFun face_rule=FaceRule(), bool trans=true) const
 Create finite element to calculate matrix and vectors.
 

Private Attributes

MoFEM::InterfacemField
 

Detailed Description

Create finite elements instances.

Create finite element instances and add operators to finite elements.

Examples
analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, and analytical_poisson_field_split.cpp.

Definition at line 861 of file PoissonOperators.hpp.

Constructor & Destructor Documentation

◆ CreateFiniteElements()

PoissonExample::CreateFiniteElements::CreateFiniteElements ( MoFEM::Interface & m_field)
inline
Examples
PoissonOperators.hpp.

Definition at line 863 of file PoissonOperators.hpp.

863: mField(m_field) {}

Member Function Documentation

◆ createFEToAssembleMatrixAndVector()

MoFEMErrorCode PoissonExample::CreateFiniteElements::createFEToAssembleMatrixAndVector ( boost::function< double(const double, const double, const double)> f_u,
boost::function< double(const double, const double, const double)> f_source,
boost::shared_ptr< ForcesAndSourcesCore > & domain_lhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & boundary_lhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & domain_rhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & boundary_rhs_fe,
bool trans = true ) const
inline

Create finite element to calculate matrix and vectors.

Examples
PoissonOperators.hpp, analytical_poisson.cpp, and analytical_poisson_field_split.cpp.

Definition at line 868 of file PoissonOperators.hpp.

876 {
878
879 // Create elements element instances
880 domain_lhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
882 boundary_lhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
884 domain_rhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
886 boundary_rhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
888
889 // Set integration rule to elements instances
890 domain_lhs_fe->getRuleHook = VolRule();
891 domain_rhs_fe->getRuleHook = VolRule();
892 boundary_lhs_fe->getRuleHook = FaceRule();
893 boundary_rhs_fe->getRuleHook = FaceRule();
894
895 // Add operators to element instances
896 // Add operator grad-grad for calculate matrix
897 domain_lhs_fe->getOpPtrVector().push_back(new OpK());
898 // Add operator to calculate source terms
899 domain_rhs_fe->getOpPtrVector().push_back(new OpF(f_source));
900 // Add operator calculating constrains matrix
901 boundary_lhs_fe->getOpPtrVector().push_back(new OpC(trans));
902 // Add operator calculating constrains vector
903 boundary_rhs_fe->getOpPtrVector().push_back(new Op_g(f_u));
904
906 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradSymTensorGrad< 1, SPACE_DIM, SPACE_DIM, 0 > OpK
Set integration rule to boundary elements.
Set integration rule to volume elements.

◆ createFEToAssembleMatrixAndVectorForNonlinearProblem()

MoFEMErrorCode PoissonExample::CreateFiniteElements::createFEToAssembleMatrixAndVectorForNonlinearProblem ( boost::function< double(const double, const double, const double)> f_u,
boost::function< double(const double, const double, const double)> f_source,
boost::function< double(const double)> a,
boost::function< double(const double)> diff_a,
boost::shared_ptr< ForcesAndSourcesCore > & domain_lhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & boundary_lhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & domain_rhs_fe,
boost::shared_ptr< ForcesAndSourcesCore > & boundary_rhs_fe,
ForcesAndSourcesCore::RuleHookFun vol_rule,
ForcesAndSourcesCore::RuleHookFun face_rule = FaceRule(),
bool trans = true ) const
inline

Create finite element to calculate matrix and vectors.

Examples
analytical_nonlinear_poisson.cpp.

Definition at line 1013 of file PoissonOperators.hpp.

1025 {
1027
1028 // Create elements element instances
1029 domain_lhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
1031 boundary_lhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
1033 domain_rhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
1035 boundary_rhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
1037
1038 // Set integration rule to elements instances
1039 domain_lhs_fe->getRuleHook = vol_rule;
1040 domain_rhs_fe->getRuleHook = vol_rule;
1041 boundary_lhs_fe->getRuleHook = face_rule;
1042 boundary_rhs_fe->getRuleHook = face_rule;
1043
1044 // Set integration rule
1045 // Crate shared vector storing values of field "u" on integration points on
1046 // element. element is local and is used to exchange data between operators.
1047 boost::shared_ptr<VectorDouble> values_at_integration_ptr =
1048 boost::make_shared<VectorDouble>();
1049 // Storing gradients of field
1050 boost::shared_ptr<MatrixDouble> grad_at_integration_ptr =
1051 boost::make_shared<MatrixDouble>();
1052 // multipliers values
1053 boost::shared_ptr<VectorDouble> multiplier_at_integration_ptr =
1054 boost::make_shared<VectorDouble>();
1055
1056 // Add operators to element instances
1057 // Add default operator to calculate field values at integration points
1058 domain_lhs_fe->getOpPtrVector().push_back(
1059 new OpCalculateScalarFieldValues("U", values_at_integration_ptr));
1060 // Add default operator to calculate field gradient at integration points
1061 domain_lhs_fe->getOpPtrVector().push_back(
1062 new OpCalculateScalarFieldGradient<3>("U", grad_at_integration_ptr));
1063 // Add operator grad-(1+u^2)grad for calculate matrix
1064 domain_lhs_fe->getOpPtrVector().push_back(new OpKt(
1065 a, diff_a, values_at_integration_ptr, grad_at_integration_ptr));
1066
1067 // Add default operator to calculate field values at integration points
1068 domain_rhs_fe->getOpPtrVector().push_back(
1069 new OpCalculateScalarFieldValues("U", values_at_integration_ptr));
1070 // Add default operator to calculate field gradient at integration points
1071 domain_rhs_fe->getOpPtrVector().push_back(
1072 new OpCalculateScalarFieldGradient<3>("U", grad_at_integration_ptr));
1073 // Add operator to calculate source terms
1074 domain_rhs_fe->getOpPtrVector().push_back(new OpResF_Domain(
1075 f_source, a, values_at_integration_ptr, grad_at_integration_ptr));
1076
1077 // Add operator calculating constrains matrix
1078 boundary_lhs_fe->getOpPtrVector().push_back(new OpC(trans));
1079
1080 // Add default operator to calculate field values at integration points
1081 boundary_rhs_fe->getOpPtrVector().push_back(
1082 new OpCalculateScalarFieldValues("U", values_at_integration_ptr));
1083 // Add default operator to calculate values of Lagrange multipliers
1084 boundary_rhs_fe->getOpPtrVector().push_back(
1085 new OpCalculateScalarFieldValues("L", multiplier_at_integration_ptr));
1086 // Add operator calculating constrains vector
1087 boundary_rhs_fe->getOpPtrVector().push_back(
1088 new OpRes_g(f_u, values_at_integration_ptr));
1089 boundary_rhs_fe->getOpPtrVector().push_back(
1090 new OpResF_Boundary(multiplier_at_integration_ptr));
1091
1093 }
constexpr double a
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Get value at integration points for scalar field.

◆ createFEToEvaluateError()

MoFEMErrorCode PoissonExample::CreateFiniteElements::createFEToEvaluateError ( boost::function< double(const double, const double, const double)> f_u,
boost::function< FTensor::Tensor1< double, 3 >(const double, const double, const double)> g_u,
Vec global_error,
boost::shared_ptr< ForcesAndSourcesCore > & domain_error ) const
inline

Create finite element to calculate error.

Examples
PoissonOperators.hpp, analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, and analytical_poisson_field_split.cpp.

Definition at line 911 of file PoissonOperators.hpp.

917 {
919 // Create finite element instance to calculate error
920 domain_error = boost::shared_ptr<ForcesAndSourcesCore>(
922 domain_error->getRuleHook = VolRule();
923 // Set integration rule
924 // Crate shared vector storing values of field "u" on integration points on
925 // element. element is local and is used to exchange data between operators.
926 boost::shared_ptr<VectorDouble> values_at_integration_ptr =
927 boost::make_shared<VectorDouble>();
928 // Storing gradients of field
929 boost::shared_ptr<MatrixDouble> grad_at_integration_ptr =
930 boost::make_shared<MatrixDouble>();
931 // Add default operator to calculate field values at integration points
932 domain_error->getOpPtrVector().push_back(
933 new OpCalculateScalarFieldValues("U", values_at_integration_ptr));
934 // Add default operator to calculate field gradient at integration points
935 domain_error->getOpPtrVector().push_back(
936 new OpCalculateScalarFieldGradient<3>("U", grad_at_integration_ptr));
937 // Add operator to integrate error element by element.
938 domain_error->getOpPtrVector().push_back(
939 new OpError(f_u, g_u, values_at_integration_ptr,
940 grad_at_integration_ptr, global_error));
942 }

◆ creatFEToPostProcessResults()

MoFEMErrorCode PoissonExample::CreateFiniteElements::creatFEToPostProcessResults ( boost::shared_ptr< PostProcFE > & post_proc_volume) const
inline

Create finite element to post-process results.

Examples
PoissonOperators.hpp, analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, and analytical_poisson_field_split.cpp.

Definition at line 947 of file PoissonOperators.hpp.

948 {
949
951
952 // Note that user can stack together arbitrary number of operators to
953 // compose complex PDEs.
954
955 // Post-process results. This is standard element, with functionality
956 // enabling refining mesh for post-processing. In addition in
957 // PostProcOnRefMesh.hpp are implanted set of users operators to
958 // post-processing fields. Here using simplified mechanism for
959 // post-processing finite element, we add operators to save data from field
960 // on mesh tags for ParaView visualization.
961 post_proc_volume =
962 boost::shared_ptr<PostProcFE>(new PostProcFE(mField));
963 // Add operators to the elements, starting with some generic
964
965 auto det_ptr = boost::make_shared<VectorDouble>();
966 auto jac_ptr = boost::make_shared<MatrixDouble>();
967 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
968 post_proc_volume->getOpPtrVector().push_back(
969 new OpCalculateHOJac<3>(jac_ptr));
970 post_proc_volume->getOpPtrVector().push_back(
971 new OpInvertMatrix<3>(jac_ptr, det_ptr, inv_jac_ptr));
972 post_proc_volume->getOpPtrVector().push_back(
973 new OpSetHOInvJacToScalarBases<3>(H1, inv_jac_ptr));
974
975 auto u_ptr = boost::make_shared<VectorDouble>();
976 auto grad_ptr = boost::make_shared<MatrixDouble>();
977 auto e_ptr = boost::make_shared<VectorDouble>();
978
979 post_proc_volume->getOpPtrVector().push_back(
980 new OpCalculateScalarFieldValues("U", u_ptr));
981 post_proc_volume->getOpPtrVector().push_back(
982 new OpCalculateScalarFieldGradient<3>("U", grad_ptr));
983 post_proc_volume->getOpPtrVector().push_back(
984 new OpCalculateScalarFieldValues("ERROR", e_ptr));
985
987
988 post_proc_volume->getOpPtrVector().push_back(
989
990 new OpPPMap(
991
992 post_proc_volume->getPostProcMesh(),
993 post_proc_volume->getMapGaussPts(),
994
995 {{"U", u_ptr}, {"ERROR", e_ptr}},
996
997 {{"GRAD", grad_ptr}},
998
999 {},
1000
1001 {}
1002
1003 )
1004
1005 );
1006
1008 }
@ H1
continuous field
Definition definitions.h:85
PostProcBrokenMeshInMoab< VolumeElementForcesAndSourcesCore > PostProcFE
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Post post-proc data at points from hash maps.
Set inverse jacobian to base functions.

Member Data Documentation

◆ mField

MoFEM::Interface& PoissonExample::CreateFiniteElements::mField
private
Examples
PoissonOperators.hpp.

Definition at line 1096 of file PoissonOperators.hpp.


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