v0.14.0
Loading...
Searching...
No Matches
Classes | Public Member Functions | Private Types | Private Member Functions | Static Private Member Functions | Private Attributes | Static Private Attributes | Friends | List of all members
Example Struct Reference

[Example] More...

Collaboration diagram for Example:
[legend]

Classes

struct  BoundaryOp
 
struct  CommonData
 [Example] More...
 
struct  DynamicFirstOrderConsConstantTimeScale
 
struct  DynamicFirstOrderConsSinusTimeScale
 
struct  OpCalcSurfaceAverageTemperature
 
struct  OpError
 
struct  OpError< 1 >
 
struct  OpFirst
 
struct  OpFluxRhs
 
struct  OpRadiationLhs
 
struct  OpRadiationRhs
 
struct  OpRhs
 
struct  OpSecond
 [Operator] More...
 
struct  OpZero
 [Common data] More...
 
struct  ScaledTimeScale
 

Public Member Functions

 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 [Run problem] More...
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 
 Example (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProblem ()
 

Private Types

enum  BoundingBox {
  CENTER_X = 0 , CENTER_Y , MAX_X , MAX_Y ,
  MIN_X , MIN_Y , LAST_BB
}
 

Private Member Functions

MoFEMErrorCode setupProblem ()
 [Run problem] More...
 
MoFEMErrorCode createCommonData ()
 [Set up problem] More...
 
MoFEMErrorCode bC ()
 [Create common data] More...
 
MoFEMErrorCode OPs ()
 [Boundary condition] More...
 
MoFEMErrorCode tsSolve ()
 
MoFEMErrorCode readMesh ()
 [Run problem] More...
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 [Set up problem] More...
 
MoFEMErrorCode assembleSystem ()
 [Push operators to pipeline] More...
 
MoFEMErrorCode solveSystem ()
 [Solve] More...
 
MoFEMErrorCode outputResults ()
 [Solve] More...
 
MoFEMErrorCode checkResults ()
 [Postprocess results] More...
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode setUp ()
 [Run all] More...
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode setFieldValues ()
 [Create common data] More...
 
MoFEMErrorCode pushOperators ()
 [Set density distribution] More...
 
MoFEMErrorCode integrateElements ()
 [Push operators to pipeline] More...
 
MoFEMErrorCode postProcess ()
 [Integrate] More...
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode setIntegrationRules ()
 [Set up problem] More...
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode assembleSystemIntensity ()
 [Calculate flux on boundary] More...
 
MoFEMErrorCode assembleSystemFlux ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode calculateFlux (double &calc_flux)
 [Set up problem] More...
 
MoFEMErrorCode outputResults (const int i)
 [Solve] More...
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode setIntegrationRules ()
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode bC ()
 
MoFEMErrorCode OPs ()
 
MoFEMErrorCode kspSolve ()
 [Push operators to pipeline] More...
 
MoFEMErrorCode postProcess ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode setIntegrationRules ()
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode addMatBlockOps (boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pipeline, std::string field_name, std::string block_name, boost::shared_ptr< MatrixDouble > mat_D_Ptr, Sev sev)
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode gettingNorms ()
 [Solve] More...
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 
MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode createCommonData ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 
MoFEMErrorCode checkResults ()
 

Static Private Member Functions

static std::pair< int, int > getCoordsInImage (double x, double y)
 
static double rhsSource (const double x, const double y, const double)
 
static double lhsFlux (const double x, const double y, const double)
 
static int integrationRule (int, int, int p_data)
 

Private Attributes

MoFEM::InterfacemField
 
boost::shared_ptr< DomainElereactionFe
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
 
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
 
SimplesimpleInterface
 
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
 
boost::shared_ptr< CommonDatacommonDataPtr
 
FieldApproximationBase base
 
FieldSpace space
 
boost::shared_ptr< VectorDouble > approxVals
 
boost::shared_ptr< MatrixDouble > approxGradVals
 
Range pinchNodes
 
PetscBool doEvalField
 
std::array< double, SPACE_DIMfieldEvalCoords
 
boost::shared_ptr< FieldEvaluatorInterface::SetPtsDatafieldEvalData
 
boost::shared_ptr< MatrixDouble > vectorFieldPtr
 
boost::shared_ptr< MatrixDouble > matDPtr
 
SmartPetscObj< Mat > M
 
SmartPetscObj< Mat > K
 
SmartPetscObj< EPS > ePS
 
std::array< SmartPetscObj< Vec >, 6 > rigidBodyMotion
 
boost::shared_ptr< FEMethoddomianLhsFEPtr
 
boost::shared_ptr< FEMethoddomianRhsFEPtr
 

Static Private Attributes

static std::vector< doublerZ
 
static std::vector< MatrixIntiI
 
static std::array< double, LAST_BBaveMaxMin
 
static int focalIndex
 
static int savitzkyGolayNormalisation
 
static const int * savitzkyGolayWeights
 
static ApproxFieldFunction< FIELD_DIMapproxFunction
 

Friends

struct TSPrePostProc
 

Detailed Description

[Example]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 226 of file plastic.cpp.

Member Enumeration Documentation

◆ BoundingBox

enum Example::BoundingBox
private
Enumerator
CENTER_X 
CENTER_Y 
MAX_X 
MAX_Y 
MIN_X 
MIN_Y 
LAST_BB 
Examples
phase.cpp.

Definition at line 98 of file phase.cpp.

98 {
99 CENTER_X = 0,
100 CENTER_Y,
101 MAX_X,
102 MAX_Y,
103 MIN_X,
104 MIN_Y,
105 LAST_BB
106 };
@ MAX_X
Definition: phase.cpp:101
@ MIN_X
Definition: phase.cpp:103
@ MIN_Y
Definition: phase.cpp:104
@ CENTER_X
Definition: phase.cpp:99
@ MAX_Y
Definition: phase.cpp:102
@ CENTER_Y
Definition: phase.cpp:100
@ LAST_BB
Definition: phase.cpp:105

Constructor & Destructor Documentation

◆ Example() [1/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 228 of file plastic.cpp.

228: mField(m_field) {}
MoFEM::Interface & mField
Definition: plastic.cpp:233

◆ Example() [2/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 421 of file dynamic_first_order_con_law.cpp.

421: mField(m_field) {}

◆ Example() [3/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 35 of file helmholtz.cpp.

35: mField(m_field) {}

◆ Example() [4/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 25 of file integration.cpp.

25: mField(m_field) {}

◆ Example() [5/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 50 of file plot_base.cpp.

50: mField(m_field) {}

◆ Example() [6/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 81 of file phase.cpp.

81: mField(m_field) {}

◆ Example() [7/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 53 of file approximaton.cpp.

53: mField(m_field) {}

◆ Example() [8/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 45 of file radiation.cpp.

45: mField(m_field) {}

◆ Example() [9/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 46 of file heat_method.cpp.

46: mField(m_field) {}

◆ Example() [10/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 87 of file elastic.cpp.

87: mField(m_field) {}

◆ Example() [11/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 49 of file eigen_elastic.cpp.

49: mField(m_field) {}

◆ Example() [12/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 59 of file nonlinear_elastic.cpp.

59: mField(m_field) {}

◆ Example() [13/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 55 of file nonlinear_dynamic_elastic.cpp.

55: mField(m_field) {}

◆ Example() [14/14]

Example::Example ( MoFEM::Interface m_field)
inline

Definition at line 339 of file shallow_wave.cpp.

339: mField(m_field) {}

Member Function Documentation

◆ addMatBlockOps()

MoFEMErrorCode Example::addMatBlockOps ( boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &  pipeline,
std::string  field_name,
std::string  block_name,
boost::shared_ptr< MatrixDouble >  mat_D_Ptr,
Sev  sev 
)
private

[Calculate elasticity tensor]

[Calculate elasticity tensor]

Definition at line 113 of file elastic.cpp.

116 {
118
119 struct OpMatBlocks : public DomainEleOp {
120 OpMatBlocks(std::string field_name, boost::shared_ptr<MatrixDouble> m,
121 double bulk_modulus_K, double shear_modulus_G,
122 MoFEM::Interface &m_field, Sev sev,
123 std::vector<const CubitMeshSets *> meshset_vec_ptr)
124 : DomainEleOp(field_name, DomainEleOp::OPROW), matDPtr(m),
125 bulkModulusKDefault(bulk_modulus_K),
126 shearModulusGDefault(shear_modulus_G) {
127 std::fill(&(doEntities[MBEDGE]), &(doEntities[MBMAXTYPE]), false);
128 CHK_THROW_MESSAGE(extractBlockData(m_field, meshset_vec_ptr, sev),
129 "Can not get data from block");
130 }
131
132 MoFEMErrorCode doWork(int side, EntityType type,
135
136 for (auto &b : blockData) {
137
138 if (b.blockEnts.find(getFEEntityHandle()) != b.blockEnts.end()) {
139 CHKERR getMatDPtr(matDPtr, b.bulkModulusK, b.shearModulusG);
141 }
142 }
143
144 CHKERR getMatDPtr(matDPtr, bulkModulusKDefault, shearModulusGDefault);
146 }
147
148 private:
149 boost::shared_ptr<MatrixDouble> matDPtr;
150
151 struct BlockData {
152 double bulkModulusK;
153 double shearModulusG;
154 Range blockEnts;
155 };
156
157 double bulkModulusKDefault;
158 double shearModulusGDefault;
159 std::vector<BlockData> blockData;
160
162 extractBlockData(MoFEM::Interface &m_field,
163 std::vector<const CubitMeshSets *> meshset_vec_ptr,
164 Sev sev) {
166
167 for (auto m : meshset_vec_ptr) {
168 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock") << *m;
169 std::vector<double> block_data;
170 CHKERR m->getAttributes(block_data);
171 if (block_data.size() < 2) {
172 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
173 "Expected that block has two attributes");
174 }
175 auto get_block_ents = [&]() {
176 Range ents;
177 CHKERR
178 m_field.get_moab().get_entities_by_handle(m->meshset, ents, true);
179 return ents;
180 };
181
182 double young_modulus = block_data[0];
183 double poisson_ratio = block_data[1];
184 double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
185 double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
186
187 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock")
188 << "E = " << young_modulus << " nu = " << poisson_ratio;
189
190 blockData.push_back(
191 {bulk_modulus_K, shear_modulus_G, get_block_ents()});
192 }
193 MOFEM_LOG_CHANNEL("WORLD");
195 }
196
197 MoFEMErrorCode getMatDPtr(boost::shared_ptr<MatrixDouble> mat_D_ptr,
198 double bulk_modulus_K, double shear_modulus_G) {
200 //! [Calculate elasticity tensor]
201 auto set_material_stiffness = [&]() {
207 double A = 1.;
208 if (SPACE_DIM == 2 && !is_plane_strain) {
209 A = 2 * shear_modulus_G /
210 (bulk_modulus_K + (4. / 3.) * shear_modulus_G);
211 }
212 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*mat_D_ptr);
213 t_D(i, j, k, l) =
214 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
215 A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) * t_kd(i, j) *
216 t_kd(k, l);
217 };
218 //! [Calculate elasticity tensor]
219 constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
220 mat_D_ptr->resize(size_symm * size_symm, 1);
221 set_material_stiffness();
223 }
224 };
225
226 pipeline.push_back(new OpMatBlocks(
228
229 // Get blockset using regular expression
231
232 (boost::format("%s(.*)") % block_name).str()
233
234 ))
235
236 ));
237
239}
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
constexpr int SPACE_DIM
Kronecker Delta class symmetric.
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
double bulk_modulus_K
double shear_modulus_G
PetscBool is_plane_strain
Definition: elastic.cpp:83
FTensor::Index< 'm', SPACE_DIM > m
constexpr auto t_kd
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
constexpr AssemblyType A
double young_modulus
Young modulus.
Definition: plastic.cpp:172
constexpr auto size_symm
Definition: plastic.cpp:42
constexpr auto field_name
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ assembleSystem() [1/10]

MoFEMErrorCode Example::assembleSystem ( )
private

[Push operators to pipeline]

[Boundary condition]

[Applying essential BC]

[Push operators to pipeline]

[Integration rule]

[Integration rule]

[Push domain stiffness matrix to pipeline]

[Push domain stiffness matrix to pipeline]

[Push Internal forces]

[Push Internal forces]

[Push Body forces]

[Push Body forces]

[Push natural boundary conditions]

[Push natural boundary conditions]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 647 of file dynamic_first_order_con_law.cpp.

647 {
649 auto *simple = mField.getInterface<Simple>();
650 auto *pipeline_mng = mField.getInterface<PipelineManager>();
651
652 auto get_body_force = [this](const double, const double, const double) {
655 t_source(i) = 0.;
656 t_source(0) = 0.1;
657 t_source(1) = 1.;
658 return t_source;
659 };
660
661 // specific time scaling
662 auto get_time_scale = [this](const double time) {
663 return sin(time * omega * M_PI);
664 };
665
666 auto apply_rhs = [&](auto &pip) {
668
670 "GEOMETRY");
671
672 // Calculate Gradient of velocity
673 auto mat_v_grad_ptr = boost::make_shared<MatrixDouble>();
675 "V", mat_v_grad_ptr));
676
677 auto gravity_vector_ptr = boost::make_shared<MatrixDouble>();
678 gravity_vector_ptr->resize(SPACE_DIM, 1);
679 auto set_body_force = [&]() {
682 auto t_force = getFTensor1FromMat<SPACE_DIM, 0>(*gravity_vector_ptr);
683 double unit_weight = 0.;
684 CHKERR PetscOptionsGetReal(PETSC_NULL, "", "-unit_weight", &unit_weight,
685 PETSC_NULL);
686 t_force(i) = 0;
687 if (SPACE_DIM == 2) {
688 t_force(1) = -unit_weight;
689 } else if (SPACE_DIM == 3) {
690 t_force(2) = unit_weight;
691 }
693 };
694
695 CHKERR set_body_force();
696 pip.push_back(new OpBodyForce("V", gravity_vector_ptr,
697 [](double, double, double) { return 1.; }));
698
699 // Calculate unknown F
700 auto mat_H_tensor_ptr = boost::make_shared<MatrixDouble>();
702 "F", mat_H_tensor_ptr));
703
704 // // Calculate F
705 double tau = 0.2;
706 CHKERR PetscOptionsGetReal(PETSC_NULL, "", "-tau", &tau, PETSC_NULL);
707
708 double xi = 0.;
709 CHKERR PetscOptionsGetReal(PETSC_NULL, "", "-xi", &xi, PETSC_NULL);
710
711 // Calculate P stab
712 auto one = [&](const double, const double, const double) {
713 return 3. * bulk_modulus_K;
714 };
715 auto minus_one = [](const double, const double, const double) {
716 return -1.;
717 };
718
719 auto mat_dot_F_tensor_ptr = boost::make_shared<MatrixDouble>();
721 "F_dot", mat_dot_F_tensor_ptr));
722
723 // Calculate Gradient of Spatial Positions
724 auto mat_x_grad_ptr = boost::make_shared<MatrixDouble>();
726 "x_2", mat_x_grad_ptr));
727
728 auto mat_F_tensor_ptr = boost::make_shared<MatrixDouble>();
730 mat_F_tensor_ptr, mat_H_tensor_ptr));
731
732 auto mat_F_stab_ptr = boost::make_shared<MatrixDouble>();
734 mat_F_tensor_ptr, mat_F_stab_ptr, mat_dot_F_tensor_ptr, tau, xi,
735 mat_x_grad_ptr, mat_v_grad_ptr));
736
737 PetscBool is_linear_elasticity = PETSC_TRUE;
738 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-is_linear_elasticity",
739 &is_linear_elasticity, PETSC_NULL);
740
741 auto mat_P_stab_ptr = boost::make_shared<MatrixDouble>();
742 if (is_linear_elasticity) {
745 mat_F_stab_ptr));
746 } else {
747 auto inv_F = boost::make_shared<MatrixDouble>();
748 auto det_ptr = boost::make_shared<VectorDouble>();
749
750 pip.push_back(
751 new OpInvertMatrix<SPACE_DIM>(mat_F_stab_ptr, det_ptr, inv_F));
752
753 // OpCalculatePiolaIncompressibleNH
756 mat_F_stab_ptr, inv_F, det_ptr));
757 }
758
759 pip.push_back(new OpGradTimesTensor2("V", mat_P_stab_ptr, minus_one));
760 pip.push_back(new OpRhsTestPiola("F", mat_v_grad_ptr, one));
761
763 };
764
765 CHKERR apply_rhs(pipeline_mng->getOpDomainExplicitRhsPipeline());
766
767 auto integration_rule = [](int, int, int approx_order) {
768 return 2 * approx_order;
769 };
770 CHKERR pipeline_mng->setDomainExplicitRhsIntegrationRule(integration_rule);
771
773}
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
@ H1
continuous field
Definition: definitions.h:85
double lamme_lambda
FormsIntegrators< DomainEleOp >::Assembly< AssemblyType::PETSC >::LinearForm< IntegrationType::GAUSS >::OpBaseTimesVector< 1, SPACE_DIM *SPACE_DIM, 1 > OpRhsTestPiola
FormsIntegrators< DomainEleOp >::Assembly< AssemblyType::PETSC >::LinearForm< IntegrationType::GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpGradTimesTensor2
auto integration_rule
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpSource< 1, 3 > OpBodyForce
PetscErrorCode PetscOptionsGetReal(PetscOptions *, const char pre[], const char name[], PetscReal *dval, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
static constexpr int approx_order
constexpr double mu
constexpr double omega
FTensor::Index< 'i', 3 > i
Add operators pushing bases from local to physical configuration.
Get values at integration pts for tensor filed rank 2, i.e. matrix field.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
PipelineManager interface.
Simple interface for fast problem set-up.
Definition: Simple.hpp:27

◆ assembleSystem() [2/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [3/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [4/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [5/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [6/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [7/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [8/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [9/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystem() [10/10]

MoFEMErrorCode Example::assembleSystem ( )
private

◆ assembleSystemFlux()

MoFEMErrorCode Example::assembleSystemFlux ( )
private
Examples
phase.cpp.

◆ assembleSystemIntensity()

MoFEMErrorCode Example::assembleSystemIntensity ( )
private

[Calculate flux on boundary]

[Push operators to pipeline]

Examples
phase.cpp.

Definition at line 433 of file phase.cpp.

433 {
435
436 auto *pipeline_mng = mField.getInterface<PipelineManager>();
437
438 pipeline_mng->getDomainLhsFE().reset();
439 pipeline_mng->getDomainRhsFE().reset();
440 pipeline_mng->getBoundaryRhsFE().reset();
441
442 auto rule_vol = [](int, int, int order) { return 2 * (order + 1); };
443 pipeline_mng->setDomainLhsIntegrationRule(rule_vol);
444 pipeline_mng->setDomainRhsIntegrationRule(rule_vol);
445
446 auto det_ptr = boost::make_shared<VectorDouble>();
447 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
448 auto jac_ptr = boost::make_shared<MatrixDouble>();
449 pipeline_mng->getOpDomainLhsPipeline().push_back(new OpSetHOWeightsOnFace());
450 pipeline_mng->getOpDomainLhsPipeline().push_back(
451 new OpCalculateHOJacForFace(jac_ptr));
452 pipeline_mng->getOpDomainLhsPipeline().push_back(
453 new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
454 pipeline_mng->getOpDomainLhsPipeline().push_back(new OpMakeHdivFromHcurl());
455 pipeline_mng->getOpDomainLhsPipeline().push_back(
457 pipeline_mng->getOpDomainLhsPipeline().push_back(
458 new OpSetInvJacHcurlFace(inv_jac_ptr));
459 pipeline_mng->getOpDomainLhsPipeline().push_back(new OpSetHOWeightsOnFace());
460
461 pipeline_mng->getOpDomainLhsPipeline().push_back(
462 new OpHdivHdiv("S", "S", lhsFlux));
463 auto unity = []() { return 1; };
464 pipeline_mng->getOpDomainLhsPipeline().push_back(
465 new OpHdivU("S", "PHI", unity, true));
466
467 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpSetHOWeightsOnFace());
468 pipeline_mng->getOpDomainRhsPipeline().push_back(
469 new OpDomainSource("PHI", rhsSource));
470
472}
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
constexpr int order
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMixDivTimesScalar< 2 > OpHdivU
OpSetInvJacHcurlFaceImpl< 2 > OpSetInvJacHcurlFace
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
OpSetContravariantPiolaTransformOnFace2DImpl< 2 > OpSetContravariantPiolaTransformOnFace2D
static double rhsSource(const double x, const double y, const double)
Definition: phase.cpp:150
static double lhsFlux(const double x, const double y, const double)
Definition: phase.cpp:165
Make Hdiv space from Hcurl space in 2d.
Modify integration weights on face to take in account higher-order geometry.
boost::shared_ptr< FEMethod > & getDomainLhsFE()
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 3, SPACE_DIM > OpHdivHdiv
[Linear elastic problem]

◆ bC() [1/2]

MoFEMErrorCode Example::bC ( )
private

[Create common data]

[Boundary condition]

Examples
plastic.cpp.

Definition at line 548 of file plastic.cpp.

548 {
550
552 auto bc_mng = mField.getInterface<BcManager>();
553 auto prb_mng = mField.getInterface<ProblemsManager>();
554
555 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_X",
556 "U", 0, 0);
557 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Y",
558 "U", 1, 1);
559 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Z",
560 "U", 2, 2);
561 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(),
562 "REMOVE_ALL", "U", 0, 3);
563
564#ifdef ADD_CONTACT
565 for (auto b : {"FIX_X", "REMOVE_X"})
566 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), b,
567 "SIGMA", 0, 0, false, true);
568 for (auto b : {"FIX_Y", "REMOVE_Y"})
569 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), b,
570 "SIGMA", 1, 1, false, true);
571 for (auto b : {"FIX_Z", "REMOVE_Z"})
572 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), b,
573 "SIGMA", 2, 2, false, true);
574 for (auto b : {"FIX_ALL", "REMOVE_ALL"})
575 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), b,
576 "SIGMA", 0, 3, false, true);
577 CHKERR bc_mng->removeBlockDOFsOnEntities(
578 simple->getProblemName(), "NO_CONTACT", "SIGMA", 0, 3, false, true);
579#endif
580
581 CHKERR bc_mng->pushMarkDOFsOnEntities<DisplacementCubitBcData>(
582 simple->getProblemName(), "U");
583
584 auto &bc_map = bc_mng->getBcMapByBlockName();
585 for (auto bc : bc_map)
586 MOFEM_LOG("PLASTICITY", Sev::verbose) << "Marker " << bc.first;
587
589}
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
Simple interface for fast problem set-up.
Definition: BcManager.hpp:25
Definition of the displacement bc data structure.
Definition: BCData.hpp:72
Problem manager is used to build and partition problems.

◆ bC() [2/2]

MoFEMErrorCode Example::bC ( )
private

◆ boundaryCondition() [1/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

[Set up problem]

[Create common data]

[Boundary condition]

[Applying essential BC]

[Define gravity vector]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 536 of file dynamic_first_order_con_law.cpp.

536 {
538
540 auto bc_mng = mField.getInterface<BcManager>();
541 auto *pipeline_mng = mField.getInterface<PipelineManager>();
542 auto time_scale = boost::make_shared<TimeScale>();
543
544 PetscBool sin_time_function = PETSC_FALSE;
545 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-sin_time_function",
546 &sin_time_function, PETSC_NULL);
547
548 if (sin_time_function)
549 time_scale = boost::make_shared<DynamicFirstOrderConsSinusTimeScale>();
550 else
551 time_scale = boost::make_shared<DynamicFirstOrderConsConstantTimeScale>();
552
553 pipeline_mng->getBoundaryExplicitRhsFE().reset();
555 pipeline_mng->getOpBoundaryExplicitRhsPipeline(), {NOSPACE}, "GEOMETRY");
556
558 pipeline_mng->getOpBoundaryExplicitRhsPipeline(), mField, "V",
559 {time_scale}, "FORCE", Sev::inform);
560
561 auto integration_rule = [](int, int, int approx_order) {
562 return 2 * approx_order;
563 };
564
565 CHKERR pipeline_mng->setBoundaryExplicitRhsIntegrationRule(integration_rule);
566 CHKERR pipeline_mng->setDomainExplicitRhsIntegrationRule(integration_rule);
567
568 CHKERR bc_mng->removeBlockDOFsOnEntities<DisplacementCubitBcData>(
569 simple->getProblemName(), "V");
570
571 auto get_pre_proc_hook = [&]() {
573 mField, pipeline_mng->getDomainExplicitRhsFE(), {time_scale});
574 };
575 pipeline_mng->getDomainExplicitRhsFE()->preProcessHook = get_pre_proc_hook();
576
578}
Class (Function) to enforce essential constrains.
Definition: Essential.hpp:39

◆ boundaryCondition() [2/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [3/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [4/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [5/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [6/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [7/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [8/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [9/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ boundaryCondition() [10/10]

MoFEMErrorCode Example::boundaryCondition ( )
private

◆ calculateFlux()

MoFEMErrorCode Example::calculateFlux ( double calc_flux)
private

[Set up problem]

[Calculate flux on boundary]

Examples
phase.cpp.

Definition at line 402 of file phase.cpp.

402 {
404 auto pipeline_mng = mField.getInterface<PipelineManager>();
405
406 pipeline_mng->getDomainLhsFE().reset();
407 pipeline_mng->getDomainRhsFE().reset();
408 pipeline_mng->getBoundaryRhsFE().reset();
409
410 auto rule_vol = [](int, int, int order) { return 2 * (order + 1); };
411 pipeline_mng->setBoundaryRhsIntegrationRule(rule_vol);
412
413 auto flux_ptr = boost::make_shared<MatrixDouble>();
414 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
416 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
417 new OpCalculateHVecVectorField<3>("S", flux_ptr));
418 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
419 new BoundaryOp(flux_ptr, calc_flux));
420
421 calc_flux = 0;
422 CHKERR pipeline_mng->loopFiniteElements();
423 double global_flux_assembeld = 0;
424 MPI_Allreduce(&calc_flux, &global_flux_assembeld, 1, MPI_DOUBLE, MPI_SUM,
425 mField.get_comm());
426 calc_flux = global_flux_assembeld;
427
429}
virtual MPI_Comm & get_comm() const =0
Get vector field for H-div approximation.

◆ checkResults() [1/12]

MoFEMErrorCode Example::checkResults ( )
private

[Postprocess results]

[Postprocessing results]

[Print results]

[Check]

[Check results]

[Test example]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 1205 of file dynamic_first_order_con_law.cpp.

1205 {
1208}

◆ checkResults() [2/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [3/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [4/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [5/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [6/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [7/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [8/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [9/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [10/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [11/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ checkResults() [12/12]

MoFEMErrorCode Example::checkResults ( )
private

◆ createCommonData() [1/8]

MoFEMErrorCode Example::createCommonData ( )
private

[Set up problem]

[Set integration rule]

[Create common data]

< true if tau order is set

< true if tau order is set

Examples
eigen_elastic.cpp, heat_method.cpp, plastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 440 of file plastic.cpp.

440 {
442
443 auto get_command_line_parameters = [&]() {
445
446 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-scale", &scale, PETSC_NULL);
447 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-young_modulus",
448 &young_modulus, PETSC_NULL);
449 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-poisson_ratio",
450 &poisson_ratio, PETSC_NULL);
451 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-hardening", &H, PETSC_NULL);
452 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-hardening_viscous", &visH,
453 PETSC_NULL);
454 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-yield_stress", &sigmaY,
455 PETSC_NULL);
456 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-cn0", &cn0, PETSC_NULL);
457 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-cn1", &cn1, PETSC_NULL);
458 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-zeta", &zeta, PETSC_NULL);
459 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-Qinf", &Qinf, PETSC_NULL);
460 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-b_iso", &b_iso, PETSC_NULL);
461 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-C1_k", &C1_k, PETSC_NULL);
462 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-large_strains",
463 &is_large_strains, PETSC_NULL);
464 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-set_timer", &set_timer,
465 PETSC_NULL);
466
467 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
468 PetscBool tau_order_is_set; ///< true if tau order is set
469 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-tau_order", &tau_order,
470 &tau_order_is_set);
471 PetscBool ep_order_is_set; ///< true if tau order is set
472 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-ep_order", &ep_order,
473 &ep_order_is_set);
474 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-geom_order", &geom_order,
475 PETSC_NULL);
476
477 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-rho", &rho, PETSC_NULL);
478 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-alpha_damping",
479 &alpha_damping, PETSC_NULL);
480
481 MOFEM_LOG("PLASTICITY", Sev::inform) << "Young modulus " << young_modulus;
482 MOFEM_LOG("PLASTICITY", Sev::inform) << "Poisson ratio " << poisson_ratio;
483 MOFEM_LOG("PLASTICITY", Sev::inform) << "Yield stress " << sigmaY;
484 MOFEM_LOG("PLASTICITY", Sev::inform) << "Hardening " << H;
485 MOFEM_LOG("PLASTICITY", Sev::inform) << "Viscous hardening " << visH;
486 MOFEM_LOG("PLASTICITY", Sev::inform) << "Saturation yield stress " << Qinf;
487 MOFEM_LOG("PLASTICITY", Sev::inform) << "Saturation exponent " << b_iso;
488 MOFEM_LOG("PLASTICITY", Sev::inform) << "Kinematic hardening " << C1_k;
489 MOFEM_LOG("PLASTICITY", Sev::inform) << "cn0 " << cn0;
490 MOFEM_LOG("PLASTICITY", Sev::inform) << "cn1 " << cn1;
491 MOFEM_LOG("PLASTICITY", Sev::inform) << "zeta " << zeta;
492
493 if (tau_order_is_set == PETSC_FALSE)
494 tau_order = order - 2;
495 if (ep_order_is_set == PETSC_FALSE)
496 ep_order = order - 1;
497
498 MOFEM_LOG("PLASTICITY", Sev::inform) << "Approximation order " << order;
499 MOFEM_LOG("PLASTICITY", Sev::inform)
500 << "Ep approximation order " << ep_order;
501 MOFEM_LOG("PLASTICITY", Sev::inform)
502 << "Tau approximation order " << tau_order;
503 MOFEM_LOG("PLASTICITY", Sev::inform)
504 << "Geometry approximation order " << geom_order;
505
506 MOFEM_LOG("PLASTICITY", Sev::inform) << "Density " << rho;
507 MOFEM_LOG("PLASTICITY", Sev::inform) << "alpha_damping " << alpha_damping;
508
509 PetscBool is_scale = PETSC_TRUE;
510 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-is_scale", &is_scale,
511 PETSC_NULL);
512 if (is_scale) {
514 }
515
516 MOFEM_LOG("PLASTICITY", Sev::inform) << "Scale " << scale;
517
518#ifdef ADD_CONTACT
519 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-cn_contact",
520 &ContactOps::cn_contact, PETSC_NULL);
521 MOFEM_LOG("CONTACT", Sev::inform)
522 << "cn_contact " << ContactOps::cn_contact;
523#endif // ADD_CONTACT
524
525 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-quasi_static",
526 &is_quasi_static, PETSC_NULL);
527 MOFEM_LOG("PLASTICITY", Sev::inform)
528 << "Is quasi static: " << (is_quasi_static ? "true" : "false");
529
531 };
532
533 CHKERR get_command_line_parameters();
534
535#ifdef ADD_CONTACT
536#ifdef PYTHON_SFD
537 sdfPythonPtr = boost::make_shared<ContactOps::SDFPython>();
538 CHKERR sdfPythonPtr->sdfInit("sdf.py");
539 ContactOps::sdfPythonWeakPtr = sdfPythonPtr;
540#endif
541#endif // ADD_CONTACT
542
544}
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
double C1_k
Kinematic hardening.
Definition: plastic.cpp:180
double Qinf
Saturation yield stress.
Definition: plastic.cpp:178
double rho
Definition: plastic.cpp:191
PetscBool is_quasi_static
Definition: plastic.cpp:190
double alpha_damping
Definition: plastic.cpp:192
double visH
Viscous hardening.
Definition: plastic.cpp:176
PetscBool set_timer
Set timer.
Definition: plastic.cpp:168
double scale
Definition: plastic.cpp:170
double zeta
Viscous hardening.
Definition: plastic.cpp:177
double H
Hardening.
Definition: plastic.cpp:175
int tau_order
Order of tau files.
Definition: plastic.cpp:186
double cn0
Definition: plastic.cpp:182
double b_iso
Saturation exponent.
Definition: plastic.cpp:179
PetscBool is_large_strains
Large strains.
Definition: plastic.cpp:167
int geom_order
Order if fixed.
Definition: plastic.cpp:188
double sigmaY
Yield stress.
Definition: plastic.cpp:174
int ep_order
Order of ep files.
Definition: plastic.cpp:187
double cn1
Definition: plastic.cpp:183

◆ createCommonData() [2/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [3/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [4/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [5/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [6/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [7/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ createCommonData() [8/8]

MoFEMErrorCode Example::createCommonData ( )
private

◆ getCoordsInImage()

std::pair< int, int > Example::getCoordsInImage ( double  x,
double  y 
)
staticprivate
Examples
phase.cpp.

Definition at line 131 of file phase.cpp.

131 {
132
133 auto &m = iI[focalIndex];
134 x -= aveMaxMin[MIN_X];
135 y -= aveMaxMin[MIN_Y];
136 x *= (m.size1() - 1) / (aveMaxMin[MAX_X] - aveMaxMin[MIN_X]);
137 y *= (m.size2() - 1) / (aveMaxMin[MAX_Y] - aveMaxMin[MIN_Y]);
138 const auto p = std::make_pair<int, int>(std::round(x), std::round(y));
139
140#ifndef NDEBUG
141 if (p.first < 0 && p.first >= m.size1())
142 THROW_MESSAGE("Wrong index");
143 if (p.second < 0 && p.second >= m.size2())
144 THROW_MESSAGE("Wrong index");
145#endif
146
147 return p;
148}
static Index< 'p', 3 > p
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
FTensor::Index< 'm', 3 > m
static std::array< double, LAST_BB > aveMaxMin
Definition: phase.cpp:110
static int focalIndex
Definition: phase.cpp:112
static std::vector< MatrixInt > iI
Definition: phase.cpp:109

◆ gettingNorms()

MoFEMErrorCode Example::gettingNorms ( )
private

[Solve]

[Getting norms]

Examples
nonlinear_elastic.cpp.

Definition at line 368 of file nonlinear_elastic.cpp.

368 {
370
372 auto dm = simple->getDM();
373
374 auto T = createDMVector(simple->getDM());
375 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
376 SCATTER_FORWARD);
377 double nrm2;
378 CHKERR VecNorm(T, NORM_2, &nrm2);
379 MOFEM_LOG("EXAMPLE", Sev::inform) << "Solution norm " << nrm2;
380
381 auto post_proc_norm_fe = boost::make_shared<DomainEle>(mField);
382
383 auto post_proc_norm_rule_hook = [](int, int, int p) -> int { return 2 * p; };
384 post_proc_norm_fe->getRuleHook = post_proc_norm_rule_hook;
385
387 post_proc_norm_fe->getOpPtrVector(), {H1});
388
389 enum NORMS { U_NORM_L2 = 0, PIOLA_NORM, LAST_NORM };
390 auto norms_vec =
392 (mField.get_comm_rank() == 0) ? LAST_NORM : 0, LAST_NORM);
393 CHKERR VecZeroEntries(norms_vec);
394
395 auto u_ptr = boost::make_shared<MatrixDouble>();
396 post_proc_norm_fe->getOpPtrVector().push_back(
398
399 post_proc_norm_fe->getOpPtrVector().push_back(
400 new OpCalcNormL2Tensor1<SPACE_DIM>(u_ptr, norms_vec, U_NORM_L2));
401
402 auto common_ptr = commonDataFactory<SPACE_DIM, GAUSS, DomainEleOp>(
403 mField, post_proc_norm_fe->getOpPtrVector(), "U", "MAT_ELASTIC",
404 Sev::inform);
405
406 post_proc_norm_fe->getOpPtrVector().push_back(
408 common_ptr->getMatFirstPiolaStress(), norms_vec, PIOLA_NORM));
409
410 CHKERR DMoFEMLoopFiniteElements(dm, simple->getDomainFEName(),
411 post_proc_norm_fe);
412
413 CHKERR VecAssemblyBegin(norms_vec);
414 CHKERR VecAssemblyEnd(norms_vec);
415
416 MOFEM_LOG_CHANNEL("SELF"); // Clear channel from old tags
417 if (mField.get_comm_rank() == 0) {
418 const double *norms;
419 CHKERR VecGetArrayRead(norms_vec, &norms);
420 MOFEM_TAG_AND_LOG("SELF", Sev::inform, "example")
421 << "norm_u: " << std::scientific << std::sqrt(norms[U_NORM_L2]);
422 MOFEM_TAG_AND_LOG("SELF", Sev::inform, "example")
423 << "norm_piola: " << std::scientific << std::sqrt(norms[PIOLA_NORM]);
424 CHKERR VecRestoreArrayRead(norms_vec, &norms);
425 }
426
428}
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:509
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:572
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1003
const double T
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
virtual int get_comm_rank() const =0
Get norm of input MatrixDouble for Tensor1.
Get norm of input MatrixDouble for Tensor2.
Get values at integration pts for tensor filed rank 1, i.e. vector field.

◆ integrateElements()

MoFEMErrorCode Example::integrateElements ( )
private

[Push operators to pipeline]

[Integrate]

Definition at line 218 of file integration.cpp.

218 {
220 // Zero global vector
221 CHKERR VecZeroEntries(commonDataPtr->petscVec);
222
223 // Integrate elements by executing operators in the pipeline
225 CHKERR pipeline_mng->loopFiniteElements();
226
227 // Assemble MPI vector
228 CHKERR VecAssemblyBegin(commonDataPtr->petscVec);
229 CHKERR VecAssemblyEnd(commonDataPtr->petscVec);
231}
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::shared_ptr< CommonData > commonDataPtr
Definition: integration.cpp:42

◆ integrationRule()

static int Example::integrationRule ( int  ,
int  ,
int  p_data 
)
inlinestaticprivate

Definition at line 52 of file radiation.cpp.

52{ return 2 * p_data; };

◆ kspSolve()

MoFEMErrorCode Example::kspSolve ( )
private

[Push operators to pipeline]

[Solve]

Definition at line 230 of file radiation.cpp.

230 {
234 auto ts = pipeline_mng->createTSIM();
235
236 double ftime = 1;
237 CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
238 CHKERR TSSetFromOptions(ts);
239 CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
240
241 auto T = createDMVector(simple->getDM());
242 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
243 SCATTER_FORWARD);
244
245 CHKERR TSSolve(ts, T);
246 CHKERR TSGetTime(ts, &ftime);
247
248 PetscInt steps, snesfails, rejects, nonlinits, linits;
249 CHKERR TSGetTimeStepNumber(ts, &steps);
250 CHKERR TSGetSNESFailures(ts, &snesfails);
251 CHKERR TSGetStepRejections(ts, &rejects);
252 CHKERR TSGetSNESIterations(ts, &nonlinits);
253 CHKERR TSGetKSPIterations(ts, &linits);
254 MOFEM_LOG_C("EXAMPLE", Sev::inform,
255 "steps %d (%d rejected, %d SNES fails), ftime %g, nonlinits "
256 "%d, linits %d",
257 steps, rejects, snesfails, ftime, nonlinits, linits);
258
260}
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
SmartPetscObj< TS > createTSIM(SmartPetscObj< DM > dm=nullptr)
Create TS (time) implicit solver.

◆ lhsFlux()

double Example::lhsFlux ( const double  x,
const double  y,
const double   
)
staticprivate
Examples
phase.cpp.

Definition at line 165 of file phase.cpp.

165 {
166 const auto idx = getCoordsInImage(x, y);
167 const auto &m = iI[focalIndex];
168 return 1. / m(idx.first, idx.second);
169}
static std::pair< int, int > getCoordsInImage(double x, double y)
Definition: phase.cpp:131

◆ OPs() [1/2]

MoFEMErrorCode Example::OPs ( )
private

[Boundary condition]

[Push operators to pipeline]

[Only used for dynamics]

[Only used for dynamics]

[Only used for dynamics]

[Only used for dynamics]

Examples
plastic.cpp.

Definition at line 593 of file plastic.cpp.

593 {
595 auto pip_mng = mField.getInterface<PipelineManager>();
597 auto bc_mng = mField.getInterface<BcManager>();
598
599 auto integration_rule_bc = [](int, int, int ao) { return 2 * ao; };
600
601 auto vol_rule = [](int, int, int ao) { return 2 * ao + geom_order - 1; };
602
603 auto add_boundary_ops_lhs_mechanical = [&](auto &pip) {
605
607 "GEOMETRY");
608 pip.push_back(new OpSetHOWeightsOnSubDim<SPACE_DIM>());
609
610 // Add Natural BCs to LHS
612 pip, mField, "U", Sev::inform);
613
614#ifdef ADD_CONTACT
615 CHKERR
616 ContactOps::opFactoryBoundaryLhs<SPACE_DIM, AT, GAUSS, BoundaryEleOp>(
617 pip, "SIGMA", "U");
618 CHKERR
619 ContactOps::opFactoryBoundaryToDomainLhs<SPACE_DIM, AT, IT, DomainEle>(
620 mField, pip, simple->getDomainFEName(), "SIGMA", "U", "GEOMETRY",
621 vol_rule);
622#endif // ADD_CONTACT
623
625 };
626
627 auto add_boundary_ops_rhs_mechanical = [&](auto &pip) {
629
631 "GEOMETRY");
632 pip.push_back(new OpSetHOWeightsOnSubDim<SPACE_DIM>());
633
634 // Add Natural BCs to RHS
636 pip, mField, "U", {boost::make_shared<ScaledTimeScale>()}, Sev::inform);
637
638#ifdef ADD_CONTACT
639 CHKERR ContactOps::opFactoryBoundaryRhs<SPACE_DIM, AT, IT, BoundaryEleOp>(
640 pip, "SIGMA", "U");
641#endif // ADD_CONTACT
642
644 };
645
646 auto add_domain_ops_lhs = [this](auto &pip) {
649 "GEOMETRY");
650
651 if (is_quasi_static == PETSC_FALSE) {
652
653 //! [Only used for dynamics]
656 //! [Only used for dynamics]
657
658 auto get_inertia_and_mass_damping = [this](const double, const double,
659 const double) {
660 auto *pip = mField.getInterface<PipelineManager>();
661 auto &fe_domain_lhs = pip->getDomainLhsFE();
662 return (rho / scale) * fe_domain_lhs->ts_aa +
663 (alpha_damping / scale) * fe_domain_lhs->ts_a;
664 };
665 pip.push_back(new OpMass("U", "U", get_inertia_and_mass_damping));
666 }
667
668 CHKERR PlasticOps::opFactoryDomainLhs<SPACE_DIM, AT, IT, DomainEleOp>(
669 mField, "MAT_PLASTIC", pip, "U", "EP", "TAU");
670
672 };
673
674 auto add_domain_ops_rhs = [this](auto &pip) {
676
678 "GEOMETRY");
679
681 pip, mField, "U",
682 {boost::make_shared<ScaledTimeScale>("body_force_hist.txt")},
683 Sev::inform);
684
685 // only in case of dynamics
686 if (is_quasi_static == PETSC_FALSE) {
687
688 //! [Only used for dynamics]
691 //! [Only used for dynamics]
692
693 auto mat_acceleration = boost::make_shared<MatrixDouble>();
695 "U", mat_acceleration));
696 pip.push_back(
697 new OpInertiaForce("U", mat_acceleration, [](double, double, double) {
698 return rho / scale;
699 }));
700 if (alpha_damping > 0) {
701 auto mat_velocity = boost::make_shared<MatrixDouble>();
702 pip.push_back(
703 new OpCalculateVectorFieldValuesDot<SPACE_DIM>("U", mat_velocity));
704 pip.push_back(
705 new OpInertiaForce("U", mat_velocity, [](double, double, double) {
706 return alpha_damping / scale;
707 }));
708 }
709 }
710
711 CHKERR PlasticOps::opFactoryDomainRhs<SPACE_DIM, AT, IT, DomainEleOp>(
712 mField, "MAT_PLASTIC", pip, "U", "EP", "TAU");
713
714#ifdef ADD_CONTACT
715 CHKERR ContactOps::opFactoryDomainRhs<SPACE_DIM, AT, IT, DomainEleOp>(
716 pip, "SIGMA", "U");
717#endif // ADD_CONTACT
718
720 };
721
722 CHKERR add_domain_ops_lhs(pip_mng->getOpDomainLhsPipeline());
723 CHKERR add_domain_ops_rhs(pip_mng->getOpDomainRhsPipeline());
724
725 // Boundary
726 CHKERR add_boundary_ops_lhs_mechanical(pip_mng->getOpBoundaryLhsPipeline());
727 CHKERR add_boundary_ops_rhs_mechanical(pip_mng->getOpBoundaryRhsPipeline());
728
729 CHKERR pip_mng->setDomainRhsIntegrationRule(vol_rule);
730 CHKERR pip_mng->setDomainLhsIntegrationRule(vol_rule);
731
732 CHKERR pip_mng->setBoundaryLhsIntegrationRule(integration_rule_bc);
733 CHKERR pip_mng->setBoundaryRhsIntegrationRule(integration_rule_bc);
734
735 auto create_reaction_pipeline = [&](auto &pip) {
738 "GEOMETRY");
739 CHKERR PlasticOps::opFactoryDomainReactions<SPACE_DIM, AT, IT, DomainEleOp>(
740 mField, "MAT_PLASTIC", pip, "U", "EP", "TAU");
742 };
743
744 reactionFe = boost::make_shared<DomainEle>(mField);
745 reactionFe->getRuleHook = vol_rule;
746 CHKERR create_reaction_pipeline(reactionFe->getOpPtrVector());
747 reactionFe->postProcessHook =
749
751}
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpBaseTimesVector< 1, 3, 1 > OpInertiaForce
constexpr AssemblyType AT
Definition: plastic.cpp:44
boost::shared_ptr< DomainEle > reactionFe
Definition: plastic.cpp:241
Class (Function) to calculate residual side diagonal.
Definition: Essential.hpp:63
Approximate field values for given petsc vector.

◆ OPs() [2/2]

MoFEMErrorCode Example::OPs ( )
private

◆ outputResults() [1/11]

MoFEMErrorCode Example::outputResults ( )
private

[Solve]

[Getting norms]

[Postprocess results]

[Postprocessing results]

[Postprocess clean]

[Postprocess clean]

[Postprocess initialise]

[Postprocess initialise]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 1183 of file dynamic_first_order_con_law.cpp.

1183 {
1185 PetscBool test_flg = PETSC_FALSE;
1186 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-test", &test_flg, PETSC_NULL);
1187 if (test_flg) {
1188 auto *simple = mField.getInterface<Simple>();
1189 auto T = createDMVector(simple->getDM());
1190 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
1191 SCATTER_FORWARD);
1192 double nrm2;
1193 CHKERR VecNorm(T, NORM_2, &nrm2);
1194 MOFEM_LOG("EXAMPLE", Sev::inform) << "Regression norm " << nrm2;
1195 constexpr double regression_value = 0.0194561;
1196 if (fabs(nrm2 - regression_value) > 1e-2)
1197 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
1198 "Regression test failed; wrong norm value.");
1199 }
1201}
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40

◆ outputResults() [2/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [3/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [4/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [5/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [6/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [7/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [8/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [9/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [10/11]

MoFEMErrorCode Example::outputResults ( )
private

◆ outputResults() [11/11]

MoFEMErrorCode Example::outputResults ( const int  i)
private

[Solve]

[Postprocess results]

Definition at line 504 of file phase.cpp.

504 {
507 pipeline_mng->getDomainLhsFE().reset();
508 pipeline_mng->getDomainRhsFE().reset();
509 pipeline_mng->getBoundaryRhsFE().reset();
510
512
513 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
514 auto jac_ptr = boost::make_shared<MatrixDouble>();
515 post_proc_fe->getOpPtrVector().push_back(
516 new OpCalculateHOJacForFace(jac_ptr));
517 post_proc_fe->getOpPtrVector().push_back(new OpMakeHdivFromHcurl());
518 post_proc_fe->getOpPtrVector().push_back(
520
521 auto s_ptr = boost::make_shared<VectorDouble>();
522 auto phi_ptr = boost::make_shared<MatrixDouble>();
523 post_proc_fe->getOpPtrVector().push_back(
524 new OpCalculateScalarFieldValues("S", s_ptr));
525 post_proc_fe->getOpPtrVector().push_back(
526 new OpCalculateHVecVectorField<3>("PHI", phi_ptr));
527
529
530 post_proc_fe->getOpPtrVector().push_back(
531
532 new OpPPMap(post_proc_fe->getPostProcMesh(),
533 post_proc_fe->getMapGaussPts(),
534
535 OpPPMap::DataMapVec{{"S", s_ptr}},
536
537 OpPPMap::DataMapMat{{"PHI", phi_ptr}},
538
540
542
543 )
544
545 );
546
547 // post_proc_fe->addFieldValuesPostProc("S");
548 // post_proc_fe->addFieldValuesPostProc("PHI");
549
550
551 pipeline_mng->getDomainRhsFE() = post_proc_fe;
552 CHKERR pipeline_mng->loopFiniteElements();
553 CHKERR post_proc_fe->writeFile("out_" + boost::lexical_cast<std::string>(i) +
554 ".h5m");
556}
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Get value at integration points for scalar field.
Post post-proc data at points from hash maps.
std::map< std::string, boost::shared_ptr< VectorDouble > > DataMapVec
std::map< std::string, boost::shared_ptr< MatrixDouble > > DataMapMat
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()

◆ postProcess() [1/2]

MoFEMErrorCode Example::postProcess ( )
private

[Integrate]

[Solve]

[Print results]

[Postprocess results]

Definition at line 235 of file integration.cpp.

235 {
237 const double *array;
238 CHKERR VecGetArrayRead(commonDataPtr->petscVec, &array);
239 if (mField.get_comm_rank() == 0) {
240 MOFEM_LOG_C("SELF", Sev::inform, "Mass %6.4e", array[CommonData::ZERO]);
241 MOFEM_LOG_C("SELF", Sev::inform,
242 "First moment of inertia [ %6.4e, %6.4e, %6.4e ]",
244 array[CommonData::FIRST_Z]);
245 MOFEM_LOG_C("SELF", Sev::inform,
246 "Second moment of inertia [ %6.4e, %6.4e, %6.4e; %6.4e %6.4e; "
247 "%6.4e ]",
251 }
252 CHKERR VecRestoreArrayRead(commonDataPtr->petscVec, &array);
254}

◆ postProcess() [2/2]

MoFEMErrorCode Example::postProcess ( )
private

◆ pushOperators()

MoFEMErrorCode Example::pushOperators ( )
private

[Set density distribution]

[Push operators to pipeline]

Definition at line 188 of file integration.cpp.

188 {
191
192 // Push an operator which calculates values of density at integration points
193 pipeline_mng->getOpDomainRhsPipeline().push_back(
195 "rho", commonDataPtr->getRhoAtIntegrationPtsPtr()));
196
197 // Push an operator to pipeline to calculate zero moment of inertia (mass)
198 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpZero(commonDataPtr));
199
200 // Push an operator to the pipeline to calculate first moment of inertaia
201 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpFirst(commonDataPtr));
202
203 // Push an operator to the pipeline to calculate second moment of inertaia
204 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpSecond(commonDataPtr));
205
206 // Set integration rule. Integration rule is equal to the polynomial order of
207 // the density field plus 2, since under the integral of the second moment of
208 // inertia term x*x is present
209 auto integration_rule = [](int, int, int p_data) { return p_data + 2; };
210
211 // Add integration rule to the element
214}
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)

◆ readMesh() [1/11]

MoFEMErrorCode Example::readMesh ( )
private

[Run problem]

[Run programme]

[run problem]

[Read mesh]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 463 of file dynamic_first_order_con_law.cpp.

463 {
466
467 CHKERR simple->getOptions();
468 CHKERR simple->loadFile();
470}

◆ readMesh() [2/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [3/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [4/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [5/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [6/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [7/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [8/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [9/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [10/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ readMesh() [11/11]

MoFEMErrorCode Example::readMesh ( )
private

◆ rhsSource()

double Example::rhsSource ( const double  x,
const double  y,
const double   
)
staticprivate
Examples
phase.cpp.

Definition at line 150 of file phase.cpp.

150 {
151 const auto idx = getCoordsInImage(x, y);
152
153 double v = 0;
154 for (auto w = 0; w != window_savitzky_golay; ++w) {
155 const auto i = focalIndex - (window_savitzky_golay - 1) / 2 + w;
156 const auto &intensity = iI[i];
157 v += intensity(idx.first, idx.second) * savitzkyGolayWeights[w];
158 }
159 v = static_cast<double>(v) / savitzkyGolayNormalisation;
160
161 const auto dz = rZ[focalIndex + 1] - rZ[focalIndex - 1];
162 return -k * v / dz;
163}
const double v
phase velocity of light in medium (cm/ns)
double w(const double g, const double t)
static int window_savitzky_golay
Definition: phase.cpp:77
static int savitzkyGolayNormalisation
Definition: phase.cpp:113
static std::vector< double > rZ
Definition: phase.cpp:108
static const int * savitzkyGolayWeights
Definition: phase.cpp:114

◆ runProblem() [1/14]

MoFEMErrorCode Example::runProblem ( )

[Run problem]

[Create common data]

[Run programme]

[Operator]

[run problem]

[Run all]

[Run problem]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 262 of file plastic.cpp.

262 {
266 CHKERR bC();
267 CHKERR OPs();
268 CHKERR tsSolve();
270}
MoFEMErrorCode tsSolve()
Definition: plastic.cpp:778
MoFEMErrorCode createCommonData()
[Set up problem]
Definition: plastic.cpp:440
MoFEMErrorCode OPs()
[Boundary condition]
Definition: plastic.cpp:593
MoFEMErrorCode setupProblem()
[Run problem]
Definition: plastic.cpp:274
MoFEMErrorCode bC()
[Create common data]
Definition: plastic.cpp:548

◆ runProblem() [2/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [3/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [4/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [5/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [6/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [7/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [8/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [9/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [10/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [11/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [12/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [13/14]

MoFEMErrorCode Example::runProblem ( )

◆ runProblem() [14/14]

MoFEMErrorCode Example::runProblem ( )

◆ setFieldValues()

MoFEMErrorCode Example::setFieldValues ( )
private

[Create common data]

[Set density distribution]

Definition at line 172 of file integration.cpp.

172 {
174 auto set_density = [&](VectorAdaptor &&field_data, double *xcoord,
175 double *ycoord, double *zcoord) {
177 field_data[0] = 1;
179 };
180 FieldBlas *field_blas;
181 CHKERR mField.getInterface(field_blas);
182 CHKERR field_blas->setVertexDofs(set_density, "rho");
184}
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:115
Basic algebra on fields.
Definition: FieldBlas.hpp:21
MoFEMErrorCode setVertexDofs(VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
Set DOFs on vertices using user function.
Definition: FieldBlas.cpp:318

◆ setIntegrationRules() [1/3]

MoFEMErrorCode Example::setIntegrationRules ( )
private

[Set up problem]

[Set integration rule]

Examples
heat_method.cpp, and plot_base.cpp.

Definition at line 203 of file plot_base.cpp.

203 {
206}

◆ setIntegrationRules() [2/3]

MoFEMErrorCode Example::setIntegrationRules ( )
private

◆ setIntegrationRules() [3/3]

MoFEMErrorCode Example::setIntegrationRules ( )
private

◆ setUp()

MoFEMErrorCode Example::setUp ( )
private

[Run all]

[Set up problem]

Definition at line 137 of file integration.cpp.

137 {
140 CHKERR simple->getOptions();
141 CHKERR simple->loadFile();
142 // Add field
143 CHKERR simple->addDomainField("rho", H1, AINSWORTH_LEGENDRE_BASE, 1);
144 constexpr int order = 1;
145 CHKERR simple->setFieldOrder("rho", order);
146 CHKERR simple->setUp();
148}
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60

◆ setupProblem() [1/13]

MoFEMErrorCode Example::setupProblem ( )
private

[Run problem]

[Read mesh]

[Set up problem]

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plastic.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 274 of file plastic.cpp.

274 {
277
278 Range domain_ents;
279 CHKERR mField.get_moab().get_entities_by_dimension(0, SPACE_DIM, domain_ents,
280 true);
281 auto get_ents_by_dim = [&](const auto dim) {
282 if (dim == SPACE_DIM) {
283 return domain_ents;
284 } else {
285 Range ents;
286 if (dim == 0)
287 CHKERR mField.get_moab().get_connectivity(domain_ents, ents, true);
288 else
289 CHKERR mField.get_moab().get_entities_by_dimension(0, dim, ents, true);
290 return ents;
291 }
292 };
293
294 auto get_base = [&]() {
295 auto domain_ents = get_ents_by_dim(SPACE_DIM);
296 if (domain_ents.empty())
297 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Empty mesh");
298 const auto type = type_from_handle(domain_ents[0]);
299 switch (type) {
300 case MBQUAD:
302 case MBHEX:
304 case MBTRI:
306 case MBTET:
308 default:
309 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Element type not handled");
310 }
311 return NOBASE;
312 };
313
314 const auto base = get_base();
315 MOFEM_LOG("PLASTICITY", Sev::inform)
316 << "Base " << ApproximationBaseNames[base];
317
318 CHKERR simple->addDomainField("U", H1, base, SPACE_DIM);
319 CHKERR simple->addDomainField("EP", L2, base, size_symm);
320 CHKERR simple->addDomainField("TAU", L2, base, 1);
321 CHKERR simple->addBoundaryField("U", H1, base, SPACE_DIM);
322
323 CHKERR simple->addDataField("GEOMETRY", H1, base, SPACE_DIM);
324
325 PetscBool order_edge = PETSC_FALSE;
326 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-order_edge", &order_edge,
327 PETSC_NULL);
328 PetscBool order_face = PETSC_FALSE;
329 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-order_face", &order_face,
330 PETSC_NULL);
331 PetscBool order_volume = PETSC_FALSE;
332 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-order_volume", &order_volume,
333 PETSC_NULL);
334
335 if (order_edge || order_face || order_volume) {
336
337 MOFEM_LOG("PLASTICITY", Sev::inform) << "Order edge " << order_edge
338 ? "true"
339 : "false";
340 MOFEM_LOG("PLASTICITY", Sev::inform) << "Order face " << order_face
341 ? "true"
342 : "false";
343 MOFEM_LOG("PLASTICITY", Sev::inform) << "Order volume " << order_volume
344 ? "true"
345 : "false";
346
347 auto ents = get_ents_by_dim(0);
348 if (order_edge)
349 ents.merge(get_ents_by_dim(1));
350 if (order_face)
351 ents.merge(get_ents_by_dim(2));
352 if (order_volume)
353 ents.merge(get_ents_by_dim(3));
354 CHKERR simple->setFieldOrder("U", order, &ents);
355 } else {
356 CHKERR simple->setFieldOrder("U", order);
357 }
358 CHKERR simple->setFieldOrder("EP", ep_order);
359 CHKERR simple->setFieldOrder("TAU", tau_order);
360
361 CHKERR simple->setFieldOrder("GEOMETRY", geom_order);
362
363#ifdef ADD_CONTACT
364 CHKERR simple->addDomainField("SIGMA", CONTACT_SPACE, DEMKOWICZ_JACOBI_BASE,
365 SPACE_DIM);
366 CHKERR simple->addBoundaryField("SIGMA", CONTACT_SPACE, DEMKOWICZ_JACOBI_BASE,
367 SPACE_DIM);
368
369 auto get_skin = [&]() {
370 Range body_ents;
371 CHKERR mField.get_moab().get_entities_by_dimension(0, SPACE_DIM, body_ents);
372 Skinner skin(&mField.get_moab());
373 Range skin_ents;
374 CHKERR skin.find_skin(0, body_ents, false, skin_ents);
375 return skin_ents;
376 };
377
378 auto filter_blocks = [&](auto skin) {
379 bool is_contact_block = false;
380 Range contact_range;
381 for (auto m :
383
384 (boost::format("%s(.*)") % "CONTACT").str()
385
386 ))
387
388 ) {
389 is_contact_block =
390 true; ///< bloks interation is collectibe, so that is set irrespective
391 ///< if there are enerities in given rank or not in the block
392 MOFEM_LOG("CONTACT", Sev::inform)
393 << "Find contact block set: " << m->getName();
394 auto meshset = m->getMeshset();
395 Range contact_meshset_range;
396 CHKERR mField.get_moab().get_entities_by_dimension(
397 meshset, SPACE_DIM - 1, contact_meshset_range, true);
398
399 CHKERR mField.getInterface<CommInterface>()->synchroniseEntities(
400 contact_meshset_range);
401 contact_range.merge(contact_meshset_range);
402 }
403 if (is_contact_block) {
404 MOFEM_LOG("SYNC", Sev::inform)
405 << "Nb entities in contact surface: " << contact_range.size();
407 skin = intersect(skin, contact_range);
408 }
409 return skin;
410 };
411
412 auto filter_true_skin = [&](auto skin) {
413 Range boundary_ents;
414 ParallelComm *pcomm =
415 ParallelComm::get_pcomm(&mField.get_moab(), MYPCOMM_INDEX);
416 CHKERR pcomm->filter_pstatus(skin, PSTATUS_SHARED | PSTATUS_MULTISHARED,
417 PSTATUS_NOT, -1, &boundary_ents);
418 return boundary_ents;
419 };
420
421 auto boundary_ents = filter_true_skin(filter_blocks(get_skin()));
422 CHKERR simple->setFieldOrder("SIGMA", 0);
423 CHKERR simple->setFieldOrder("SIGMA", order - 1, &boundary_ents);
424#endif
425
426 CHKERR simple->setUp();
427 CHKERR simple->addFieldToEmptyFieldBlocks("U", "TAU");
428
429 auto project_ho_geometry = [&]() {
430 Projection10NodeCoordsOnField ent_method(mField, "GEOMETRY");
431 return mField.loop_dofs("GEOMETRY", ent_method);
432 };
433 CHKERR project_ho_geometry();
434
436}
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:345
@ NOBASE
Definition: definitions.h:59
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
@ L2
field with C-1 continuity
Definition: definitions.h:88
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
@ MOFEM_NOT_FOUND
Definition: definitions.h:33
static const char *const ApproximationBaseNames[]
Definition: definitions.h:72
const int dim
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1918
constexpr FieldSpace CONTACT_SPACE
Definition: plastic.cpp:52
FieldApproximationBase base
Definition: plot_base.cpp:68
Managing BitRefLevels.
Projection of edge entities with one mid-node on hierarchical basis.

◆ setupProblem() [2/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [3/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [4/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [5/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [6/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [7/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [8/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [9/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [10/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [11/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [12/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ setupProblem() [13/13]

MoFEMErrorCode Example::setupProblem ( )
private

◆ solveSystem() [1/11]

MoFEMErrorCode Example::solveSystem ( )
private

[Solve]

[Push operators to pipeline]

[Solve]

< Mass matrix

< Linear solver

Examples
dynamic_first_order_con_law.cpp, eigen_elastic.cpp, heat_method.cpp, helmholtz.cpp, nonlinear_elastic.cpp, phase.cpp, plot_base.cpp, and shallow_wave.cpp.

Definition at line 893 of file dynamic_first_order_con_law.cpp.

893 {
895 auto *simple = mField.getInterface<Simple>();
896 auto *pipeline_mng = mField.getInterface<PipelineManager>();
897
898 auto dm = simple->getDM();
899
900 auto calculate_stress_ops = [&](auto &pip) {
902
903 auto v_ptr = boost::make_shared<MatrixDouble>();
904 pip.push_back(new OpCalculateVectorFieldValues<SPACE_DIM>("V", v_ptr));
905 auto X_ptr = boost::make_shared<MatrixDouble>();
906 pip.push_back(
907 new OpCalculateVectorFieldValues<SPACE_DIM>("GEOMETRY", X_ptr));
908
909 auto x_ptr = boost::make_shared<MatrixDouble>();
910 pip.push_back(new OpCalculateVectorFieldValues<SPACE_DIM>("x_1", x_ptr));
911
912 // Calculate unknown F
913 auto mat_H_tensor_ptr = boost::make_shared<MatrixDouble>();
915 "F", mat_H_tensor_ptr));
916
917 auto u_ptr = boost::make_shared<MatrixDouble>();
918 pip.push_back(new OpCalculateDisplacement<SPACE_DIM>(x_ptr, X_ptr, u_ptr));
919 // Calculate P
920
921 auto mat_F_ptr = boost::make_shared<MatrixDouble>();
923 mat_F_ptr, mat_H_tensor_ptr));
924
925 PetscBool is_linear_elasticity = PETSC_TRUE;
926 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-is_linear_elasticity",
927 &is_linear_elasticity, PETSC_NULL);
928
929 auto mat_P_ptr = boost::make_shared<MatrixDouble>();
930 if (is_linear_elasticity) {
933 mat_F_ptr));
934 } else {
935 auto inv_F = boost::make_shared<MatrixDouble>();
936 auto det_ptr = boost::make_shared<VectorDouble>();
937
938 pip.push_back(new OpInvertMatrix<SPACE_DIM>(mat_F_ptr, det_ptr, inv_F));
939
942 mat_F_ptr, inv_F, det_ptr));
943 }
944
945 auto mat_v_grad_ptr = boost::make_shared<MatrixDouble>();
947 "V", mat_v_grad_ptr));
948
949 return boost::make_tuple(v_ptr, X_ptr, x_ptr, mat_P_ptr, mat_F_ptr, u_ptr);
950 };
951
952 auto post_proc_boundary = [&]() {
953 auto boundary_post_proc_fe = boost::make_shared<PostProcFaceEle>(mField);
954
956 boundary_post_proc_fe->getOpPtrVector(), {}, "GEOMETRY");
957 auto op_loop_side =
958 new OpLoopSide<SideEle>(mField, simple->getDomainFEName(), SPACE_DIM);
959 // push ops to side element, through op_loop_side operator
960 auto [boundary_v_ptr, boundary_X_ptr, boundary_x_ptr, boundary_mat_P_ptr,
961 boundary_mat_F_ptr, boundary_u_ptr] =
962 calculate_stress_ops(op_loop_side->getOpPtrVector());
963 boundary_post_proc_fe->getOpPtrVector().push_back(op_loop_side);
964
966
967 boundary_post_proc_fe->getOpPtrVector().push_back(
968
969 new OpPPMap(
970
971 boundary_post_proc_fe->getPostProcMesh(),
972 boundary_post_proc_fe->getMapGaussPts(),
973
975
976 OpPPMap::DataMapMat{{"V", boundary_v_ptr},
977 {"GEOMETRY", boundary_X_ptr},
978 {"x", boundary_x_ptr},
979 {"U", boundary_u_ptr}},
980
981 OpPPMap::DataMapMat{{"FIRST_PIOLA", boundary_mat_P_ptr},
982 {"F", boundary_mat_F_ptr}},
983
985
986 )
987
988 );
989 return boundary_post_proc_fe;
990 };
991
992 // Add monitor to time solver
993
994 double rho = 1.;
995 CHKERR PetscOptionsGetReal(PETSC_NULL, "", "-density", &rho, PETSC_NULL);
996 auto get_rho = [rho](const double, const double, const double) {
997 return rho;
998 };
999
1000 SmartPetscObj<Mat> M; ///< Mass matrix
1001 SmartPetscObj<KSP> ksp; ///< Linear solver
1002
1003 auto ts_pre_post_proc = boost::make_shared<TSPrePostProc>();
1004 tsPrePostProc = ts_pre_post_proc;
1005
1007 CHKERR MatZeroEntries(M);
1008
1009 boost::shared_ptr<DomainEle> vol_mass_ele(new DomainEle(mField));
1010
1011 vol_mass_ele->B = M;
1012
1013 auto integration_rule = [](int, int, int approx_order) {
1014 return 2 * approx_order;
1015 };
1016
1017 vol_mass_ele->getRuleHook = integration_rule;
1018
1019 vol_mass_ele->getOpPtrVector().push_back(new OpMassV("V", "V", get_rho));
1020 vol_mass_ele->getOpPtrVector().push_back(new OpMassF("F", "F"));
1021
1022 CHKERR DMoFEMLoopFiniteElements(dm, simple->getDomainFEName(), vol_mass_ele);
1023 CHKERR MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY);
1024 CHKERR MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY);
1025
1026 auto lumpVec = createDMVector(simple->getDM());
1027 CHKERR MatGetRowSum(M, lumpVec);
1028
1029 CHKERR MatZeroEntries(M);
1030 CHKERR MatDiagonalSet(M, lumpVec, INSERT_VALUES);
1031
1032 // Create and septup KSP (linear solver), we need this to calculate g(t,u) =
1033 // M^-1G(t,u)
1034 ksp = createKSP(mField.get_comm());
1035 CHKERR KSPSetOperators(ksp, M, M);
1036 CHKERR KSPSetFromOptions(ksp);
1037 CHKERR KSPSetUp(ksp);
1038
1039 auto solve_boundary_for_g = [&]() {
1041 if (*(pipeline_mng->getBoundaryExplicitRhsFE()->vecAssembleSwitch)) {
1042
1043 CHKERR VecGhostUpdateBegin(pipeline_mng->getBoundaryExplicitRhsFE()->ts_F,
1044 ADD_VALUES, SCATTER_REVERSE);
1045 CHKERR VecGhostUpdateEnd(pipeline_mng->getBoundaryExplicitRhsFE()->ts_F,
1046 ADD_VALUES, SCATTER_REVERSE);
1047 CHKERR VecAssemblyBegin(pipeline_mng->getBoundaryExplicitRhsFE()->ts_F);
1048 CHKERR VecAssemblyEnd(pipeline_mng->getBoundaryExplicitRhsFE()->ts_F);
1049 *(pipeline_mng->getBoundaryExplicitRhsFE()->vecAssembleSwitch) = false;
1050
1051 auto D =
1052 smartVectorDuplicate(pipeline_mng->getBoundaryExplicitRhsFE()->ts_F);
1053 CHKERR KSPSolve(ksp, pipeline_mng->getBoundaryExplicitRhsFE()->ts_F, D);
1054 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
1055 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
1056 CHKERR VecCopy(D, pipeline_mng->getBoundaryExplicitRhsFE()->ts_F);
1057 }
1058
1060 };
1061
1062 pipeline_mng->getBoundaryExplicitRhsFE()->postProcessHook =
1063 solve_boundary_for_g;
1064
1066 ts = pipeline_mng->createTSEX(dm);
1067
1068 // Field eval
1069 PetscBool field_eval_flag = PETSC_TRUE;
1070 boost::shared_ptr<MatrixDouble> velocity_field_ptr;
1071 boost::shared_ptr<MatrixDouble> geometry_field_ptr;
1072 boost::shared_ptr<MatrixDouble> spatial_position_field_ptr;
1073 boost::shared_ptr<SetPtsData> field_eval_data;
1074
1075 std::array<double, 3> field_eval_coords = {0.5, 0.5, 5.};
1076 int dim = 3;
1077 CHKERR PetscOptionsGetRealArray(NULL, NULL, "-field_eval_coords",
1078 field_eval_coords.data(), &dim,
1079 &field_eval_flag);
1080
1081 if (field_eval_flag) {
1082 field_eval_data =
1083 mField.getInterface<FieldEvaluatorInterface>()->getData<DomainEle>();
1084 if (SPACE_DIM == 3) {
1086 field_eval_data, simple->getDomainFEName());
1087 } else {
1089 field_eval_data, simple->getDomainFEName());
1090 }
1091
1092 field_eval_data->setEvalPoints(field_eval_coords.data(), 1);
1093
1094 auto no_rule = [](int, int, int) { return -1; };
1095
1096 auto fe_ptr = field_eval_data->feMethodPtr.lock();
1097 fe_ptr->getRuleHook = no_rule;
1098 velocity_field_ptr = boost::make_shared<MatrixDouble>();
1099 geometry_field_ptr = boost::make_shared<MatrixDouble>();
1100 spatial_position_field_ptr = boost::make_shared<MatrixDouble>();
1101 fe_ptr->getOpPtrVector().push_back(
1102 new OpCalculateVectorFieldValues<SPACE_DIM>("V", velocity_field_ptr));
1103 fe_ptr->getOpPtrVector().push_back(
1105 geometry_field_ptr));
1106 fe_ptr->getOpPtrVector().push_back(
1108 "x_2", spatial_position_field_ptr));
1109 }
1110
1111 auto post_proc_domain = [&]() {
1112 auto post_proc_fe_vol = boost::make_shared<PostProcEle>(mField);
1113
1115
1116 auto [boundary_v_ptr, boundary_X_ptr, boundary_x_ptr, boundary_mat_P_ptr,
1117 boundary_mat_F_ptr, boundary_u_ptr] =
1118 calculate_stress_ops(post_proc_fe_vol->getOpPtrVector());
1119
1120 post_proc_fe_vol->getOpPtrVector().push_back(
1121
1122 new OpPPMap(
1123
1124 post_proc_fe_vol->getPostProcMesh(),
1125 post_proc_fe_vol->getMapGaussPts(),
1126
1127 {},
1128
1129 {{"V", boundary_v_ptr},
1130 {"GEOMETRY", boundary_X_ptr},
1131 {"x", boundary_x_ptr},
1132 {"U", boundary_u_ptr}},
1133
1134 {{"FIRST_PIOLA", boundary_mat_P_ptr}, {"F", boundary_mat_F_ptr}},
1135
1136 {}
1137
1138 )
1139
1140 );
1141 return post_proc_fe_vol;
1142 };
1143
1144 boost::shared_ptr<FEMethod> null_fe;
1145 auto monitor_ptr = boost::make_shared<Monitor>(
1146 SmartPetscObj<DM>(dm, true), mField, post_proc_domain(),
1147 post_proc_boundary(), velocity_field_ptr, spatial_position_field_ptr,
1148 geometry_field_ptr, field_eval_coords, field_eval_data);
1149
1150 CHKERR DMMoFEMTSSetMonitor(dm, ts, simple->getDomainFEName(), null_fe,
1151 null_fe, monitor_ptr);
1152
1153 double ftime = 1;
1154 // CHKERR TSSetMaxTime(ts, ftime);
1155 CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
1156
1157 auto T = createDMVector(simple->getDM());
1158 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
1159 SCATTER_FORWARD);
1160 CHKERR TSSetSolution(ts, T);
1161 CHKERR TSSetFromOptions(ts);
1162
1163 auto fb = mField.getInterface<FieldBlas>();
1164
1165 CHKERR TSSetPostStage(ts, TSPrePostProc::tsPostStage);
1166 CHKERR TSSetPostStep(ts, TSPrePostProc::tsPostStep);
1167 CHKERR TSSetPreStep(ts, TSPrePostProc::tsPreStep);
1168
1169 boost::shared_ptr<ForcesAndSourcesCore> null;
1170
1171 if (auto ptr = tsPrePostProc.lock()) {
1172 ptr->fsRawPtr = this;
1173 CHKERR TSSetUp(ts);
1174 CHKERR TSSolve(ts, NULL);
1175 CHKERR TSGetTime(ts, &ftime);
1176 }
1177
1179}
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMassV
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM *SPACE_DIM > OpMassF
static boost::weak_ptr< TSPrePostProc > tsPrePostProc
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMoFEM.cpp:1183
double D
auto createKSP(MPI_Comm comm)
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMoFEM.cpp:1042
DEPRECATED SmartPetscObj< Vec > smartVectorDuplicate(Vec vec)
PetscErrorCode PetscOptionsGetRealArray(PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
ElementsAndOps< FE_DIM >::DomainEle DomainEle
SmartPetscObj< Mat > M
Field evaluator interface.
Element used to execute operators on side of the element.
intrusive_ptr for managing petsc objects
static MoFEMErrorCode tsPostStep(TS ts)
static MoFEMErrorCode tsPreStep(TS ts)
static MoFEMErrorCode tsPostStage(TS ts, PetscReal stagetime, PetscInt stageindex, Vec *Y)
[Boundary condition]

◆ solveSystem() [2/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [3/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [4/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [5/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [6/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [7/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [8/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [9/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [10/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ solveSystem() [11/11]

MoFEMErrorCode Example::solveSystem ( )
private

◆ tsSolve()

MoFEMErrorCode Example::tsSolve ( )
private
Examples
plastic.cpp.

Definition at line 778 of file plastic.cpp.

778 {
780
783 ISManager *is_manager = mField.getInterface<ISManager>();
784
785 auto snes_ctx_ptr = getDMSnesCtx(simple->getDM());
786
787 auto set_section_monitor = [&](auto solver) {
789 SNES snes;
790 CHKERR TSGetSNES(solver, &snes);
791 CHKERR SNESMonitorSet(snes,
792 (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal,
794 (void *)(snes_ctx_ptr.get()), nullptr);
796 };
797
798 auto create_post_process_elements = [&]() {
799 auto pp_fe = boost::make_shared<PostProcEle>(mField);
800 auto &pip = pp_fe->getOpPtrVector();
801
802 auto push_vol_ops = [this](auto &pip) {
804 "GEOMETRY");
805
806 auto [common_plastic_ptr, common_henky_ptr] =
807 PlasticOps::createCommonPlasticOps<SPACE_DIM, IT, DomainEleOp>(
808 mField, "MAT_PLASTIC", pip, "U", "EP", "TAU", 1., Sev::inform);
809
810 if (common_henky_ptr) {
811 if (common_plastic_ptr->mGradPtr != common_henky_ptr->matGradPtr)
812 CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "Wrong pointer for grad");
813 }
814
815 return std::make_pair(common_plastic_ptr, common_henky_ptr);
816 };
817
818 auto push_vol_post_proc_ops = [this](auto &pp_fe, auto &&p) {
820
821 auto &pip = pp_fe->getOpPtrVector();
822
823 auto [common_plastic_ptr, common_henky_ptr] = p;
824
826
827 auto x_ptr = boost::make_shared<MatrixDouble>();
828 pip.push_back(
829 new OpCalculateVectorFieldValues<SPACE_DIM>("GEOMETRY", x_ptr));
830 auto u_ptr = boost::make_shared<MatrixDouble>();
831 pip.push_back(new OpCalculateVectorFieldValues<SPACE_DIM>("U", u_ptr));
832
833 if (is_large_strains) {
834
835 pip.push_back(
836
837 new OpPPMap(
838
839 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
840
841 {{"PLASTIC_SURFACE",
842 common_plastic_ptr->getPlasticSurfacePtr()},
843 {"PLASTIC_MULTIPLIER",
844 common_plastic_ptr->getPlasticTauPtr()}},
845
846 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
847
848 {{"GRAD", common_henky_ptr->matGradPtr},
849 {"FIRST_PIOLA", common_henky_ptr->getMatFirstPiolaStress()}},
850
851 {{"PLASTIC_STRAIN", common_plastic_ptr->getPlasticStrainPtr()},
852 {"PLASTIC_FLOW", common_plastic_ptr->getPlasticFlowPtr()}}
853
854 )
855
856 );
857
858 } else {
859
860 pip.push_back(
861
862 new OpPPMap(
863
864 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
865
866 {{"PLASTIC_SURFACE",
867 common_plastic_ptr->getPlasticSurfacePtr()},
868 {"PLASTIC_MULTIPLIER",
869 common_plastic_ptr->getPlasticTauPtr()}},
870
871 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
872
873 {},
874
875 {{"STRAIN", common_plastic_ptr->mStrainPtr},
876 {"STRESS", common_plastic_ptr->mStressPtr},
877 {"PLASTIC_STRAIN", common_plastic_ptr->getPlasticStrainPtr()},
878 {"PLASTIC_FLOW", common_plastic_ptr->getPlasticFlowPtr()}}
879
880 )
881
882 );
883 }
884
886 };
887
888 auto vol_post_proc = [this, push_vol_post_proc_ops, push_vol_ops]() {
889 PetscBool post_proc_vol = PETSC_FALSE;
890 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-post_proc_vol",
891 &post_proc_vol, PETSC_NULL);
892 if (post_proc_vol == PETSC_FALSE)
893 return boost::shared_ptr<PostProcEle>();
894 auto pp_fe = boost::make_shared<PostProcEle>(mField);
896 push_vol_post_proc_ops(pp_fe, push_vol_ops(pp_fe->getOpPtrVector())),
897 "push_vol_post_proc_ops");
898 return pp_fe;
899 };
900
901 auto skin_post_proc = [this, push_vol_post_proc_ops, push_vol_ops]() {
902 PetscBool post_proc_skin = PETSC_TRUE;
903 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-post_proc_skin",
904 &post_proc_skin, PETSC_NULL);
905 if (post_proc_skin == PETSC_FALSE)
906 return boost::shared_ptr<SkinPostProcEle>();
907
909 auto pp_fe = boost::make_shared<SkinPostProcEle>(mField);
910 auto op_side = new OpLoopSide<SideEle>(
911 mField, simple->getDomainFEName(), SPACE_DIM, Sev::verbose);
912 pp_fe->getOpPtrVector().push_back(op_side);
913 CHK_MOAB_THROW(push_vol_post_proc_ops(
914 pp_fe, push_vol_ops(op_side->getOpPtrVector())),
915 "push_vol_post_proc_ops");
916 return pp_fe;
917 };
918
919 return std::make_pair(vol_post_proc(), skin_post_proc());
920 };
921
922 auto scatter_create = [&](auto D, auto coeff) {
924 CHKERR is_manager->isCreateProblemFieldAndRank(simple->getProblemName(),
925 ROW, "U", coeff, coeff, is);
926 int loc_size;
927 CHKERR ISGetLocalSize(is, &loc_size);
928 Vec v;
929 CHKERR VecCreateMPI(mField.get_comm(), loc_size, PETSC_DETERMINE, &v);
930 VecScatter scatter;
931 CHKERR VecScatterCreate(D, is, v, PETSC_NULL, &scatter);
932 return std::make_tuple(SmartPetscObj<Vec>(v),
934 };
935
936 auto set_time_monitor = [&](auto dm, auto solver) {
938 boost::shared_ptr<Monitor> monitor_ptr(
939 new Monitor(dm, create_post_process_elements(), reactionFe, uXScatter,
941 boost::shared_ptr<ForcesAndSourcesCore> null;
942 CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
943 monitor_ptr, null, null);
945 };
946
947 auto set_schur_pc = [&](auto solver,
948 boost::shared_ptr<SetUpSchur> &schur_ptr) {
950
951 auto bc_mng = mField.getInterface<BcManager>();
952 auto name_prb = simple->getProblemName();
953
954 // create sub dm for Schur complement
955 auto create_sub_u_dm = [&](SmartPetscObj<DM> base_dm,
956 SmartPetscObj<DM> &dm_sub) {
958 dm_sub = createDM(mField.get_comm(), "DMMOFEM");
959 CHKERR DMMoFEMCreateSubDM(dm_sub, base_dm, "SUB_U");
960 CHKERR DMMoFEMSetSquareProblem(dm_sub, PETSC_TRUE);
961 CHKERR DMMoFEMAddElement(dm_sub, simple->getDomainFEName());
962 CHKERR DMMoFEMAddElement(dm_sub, simple->getBoundaryFEName());
963 for (auto f : {"U"}) {
966 }
967 CHKERR DMSetUp(dm_sub);
968
970 };
971
972 // Create nested (sub BC) Schur DM
973 if constexpr (AT == AssemblyType::SCHUR) {
974 SmartPetscObj<IS> is_epp;
975 CHKERR mField.getInterface<ISManager>()->isCreateProblemFieldAndRank(
976 simple->getProblemName(), ROW, "EP", 0, MAX_DOFS_ON_ENTITY, is_epp);
977 SmartPetscObj<IS> is_tau;
978 CHKERR mField.getInterface<ISManager>()->isCreateProblemFieldAndRank(
979 simple->getProblemName(), ROW, "TAU", 0, MAX_DOFS_ON_ENTITY, is_tau);
980
981 IS is_union_raw;
982 CHKERR ISExpand(is_epp, is_tau, &is_union_raw);
983 SmartPetscObj<IS> is_union(is_union_raw);
984
985#ifdef ADD_CONTACT
986 auto add_sigma_to_is = [&](auto is_union) {
987 SmartPetscObj<IS> is_union_sigma;
988 auto add_sigma_to_is_impl = [&]() {
990 SmartPetscObj<IS> is_sigma;
991 CHKERR mField.getInterface<ISManager>()->isCreateProblemFieldAndRank(
992 simple->getProblemName(), ROW, "SIGMA", 0, MAX_DOFS_ON_ENTITY,
993 is_sigma);
994 IS is_union_raw_sigma;
995 CHKERR ISExpand(is_union, is_sigma, &is_union_raw_sigma);
996 is_union_sigma = SmartPetscObj<IS>(is_union_raw_sigma);
998 };
999 CHK_THROW_MESSAGE(add_sigma_to_is_impl(), "Can not add sigma to IS");
1000 return is_union_sigma;
1001 };
1002 is_union = add_sigma_to_is(is_union);
1003#endif // ADD_CONTACT
1004
1005 SmartPetscObj<DM> dm_u_sub;
1006 CHKERR create_sub_u_dm(simple->getDM(), dm_u_sub);
1007
1008 // Indices has to be map fro very to level, while assembling Schur
1009 // complement.
1010 auto is_up = getDMSubData(dm_u_sub)->getSmartRowIs();
1011 auto ao_up = createAOMappingIS(is_up, PETSC_NULL);
1012 schur_ptr =
1013 SetUpSchur::createSetUpSchur(mField, dm_u_sub, is_union, ao_up);
1014 CHKERR schur_ptr->setUp(solver);
1015 }
1016
1018 };
1019
1020 auto dm = simple->getDM();
1021 auto D = createDMVector(dm);
1022 auto DD = vectorDuplicate(D);
1023 uXScatter = scatter_create(D, 0);
1024 uYScatter = scatter_create(D, 1);
1025 if constexpr (SPACE_DIM == 3)
1026 uZScatter = scatter_create(D, 2);
1027
1028 auto create_solver = [pip_mng]() {
1029 if (is_quasi_static == PETSC_TRUE)
1030 return pip_mng->createTSIM();
1031 else
1032 return pip_mng->createTSIM2();
1033 };
1034
1035 auto solver = create_solver();
1036
1037 auto active_pre_lhs = []() {
1039 std::fill(PlasticOps::CommonData::activityData.begin(),
1042 };
1043
1044 auto active_post_lhs = [&]() {
1046 auto get_iter = [&]() {
1047 SNES snes;
1048 CHK_THROW_MESSAGE(TSGetSNES(solver, &snes), "Can not get SNES");
1049 int iter;
1050 CHK_THROW_MESSAGE(SNESGetIterationNumber(snes, &iter),
1051 "Can not get iter");
1052 return iter;
1053 };
1054
1055 auto iter = get_iter();
1056 if (iter >= 0) {
1057
1058 std::array<int, 5> activity_data;
1059 std::fill(activity_data.begin(), activity_data.end(), 0);
1060 MPI_Allreduce(PlasticOps::CommonData::activityData.data(),
1061 activity_data.data(), activity_data.size(), MPI_INT,
1062 MPI_SUM, mField.get_comm());
1063
1064 int &active_points = activity_data[0];
1065 int &avtive_full_elems = activity_data[1];
1066 int &avtive_elems = activity_data[2];
1067 int &nb_points = activity_data[3];
1068 int &nb_elements = activity_data[4];
1069
1070 if (nb_points) {
1071
1072 double proc_nb_points =
1073 100 * static_cast<double>(active_points) / nb_points;
1074 double proc_nb_active =
1075 100 * static_cast<double>(avtive_elems) / nb_elements;
1076 double proc_nb_full_active = 100;
1077 if (avtive_elems)
1078 proc_nb_full_active =
1079 100 * static_cast<double>(avtive_full_elems) / avtive_elems;
1080
1081 MOFEM_LOG_C("PLASTICITY", Sev::inform,
1082 "Iter %d nb pts %d nb avtive pts %d (%3.3f\%) nb active "
1083 "elements %d "
1084 "(%3.3f\%) nb full active elems %d (%3.3f\%)",
1085 iter, nb_points, active_points, proc_nb_points,
1086 avtive_elems, proc_nb_active, avtive_full_elems,
1087 proc_nb_full_active, iter);
1088 }
1089 }
1090
1092 };
1093
1094 auto add_active_dofs_elem = [&](auto dm) {
1096 auto fe_pre_proc = boost::make_shared<FEMethod>();
1097 fe_pre_proc->preProcessHook = active_pre_lhs;
1098 auto fe_post_proc = boost::make_shared<FEMethod>();
1099 fe_post_proc->postProcessHook = active_post_lhs;
1100 auto ts_ctx_ptr = getDMTsCtx(dm);
1101 ts_ctx_ptr->getPreProcessIJacobian().push_front(fe_pre_proc);
1102 ts_ctx_ptr->getPostProcessIJacobian().push_back(fe_post_proc);
1104 };
1105
1106 auto set_essential_bc = [&](auto dm, auto solver) {
1108 // This is low level pushing finite elements (pipelines) to solver
1109
1110 auto pre_proc_ptr = boost::make_shared<FEMethod>();
1111 auto post_proc_rhs_ptr = boost::make_shared<FEMethod>();
1112 auto post_proc_lhs_ptr = boost::make_shared<FEMethod>();
1113
1114 // Add boundary condition scaling
1115 auto disp_time_scale = boost::make_shared<TimeScale>();
1116
1117 auto get_bc_hook_rhs = [this, pre_proc_ptr, disp_time_scale]() {
1119 {disp_time_scale}, false);
1120 return hook;
1121 };
1122 pre_proc_ptr->preProcessHook = get_bc_hook_rhs();
1123
1124 auto get_post_proc_hook_rhs = [this, post_proc_rhs_ptr]() {
1127 mField, post_proc_rhs_ptr, nullptr, Sev::verbose)();
1129 mField, post_proc_rhs_ptr, 1.)();
1131 };
1132 auto get_post_proc_hook_lhs = [this, post_proc_lhs_ptr]() {
1134 mField, post_proc_lhs_ptr, 1.);
1135 };
1136 post_proc_rhs_ptr->postProcessHook = get_post_proc_hook_rhs;
1137
1138 auto ts_ctx_ptr = getDMTsCtx(dm);
1139 ts_ctx_ptr->getPreProcessIFunction().push_front(pre_proc_ptr);
1140 ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_ptr);
1141 ts_ctx_ptr->getPostProcessIFunction().push_back(post_proc_rhs_ptr);
1142
1143 SNES snes;
1144 CHKERR TSGetSNES(solver, &snes);
1145 KSP ksp;
1146 CHKERR SNESGetKSP(snes, &ksp);
1147 PC pc;
1148 CHKERR KSPGetPC(ksp, &pc);
1149 PetscBool is_pcfs = PETSC_FALSE;
1150 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
1151
1152 if (is_pcfs == PETSC_FALSE) {
1153 post_proc_lhs_ptr->postProcessHook = get_post_proc_hook_lhs();
1154 ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_lhs_ptr);
1155 }
1157 };
1158
1159 if (is_quasi_static == PETSC_TRUE) {
1160 CHKERR TSSetSolution(solver, D);
1161 } else {
1162 CHKERR TS2SetSolution(solver, D, DD);
1163 }
1164
1165 CHKERR set_section_monitor(solver);
1166 CHKERR set_time_monitor(dm, solver);
1167 CHKERR TSSetFromOptions(solver);
1168 CHKERR TSSetUp(solver);
1169
1170 CHKERR add_active_dofs_elem(dm);
1171 boost::shared_ptr<SetUpSchur> schur_ptr;
1172 CHKERR set_schur_pc(solver, schur_ptr);
1173 CHKERR set_essential_bc(dm, solver);
1174
1175 MOFEM_LOG_CHANNEL("TIMER");
1176 MOFEM_LOG_TAG("TIMER", "timer");
1177 if (set_timer)
1178 BOOST_LOG_SCOPED_THREAD_ATTR("Timeline", attrs::timer());
1179 MOFEM_LOG("TIMER", Sev::verbose) << "TSSetUp";
1180 CHKERR TSSetUp(solver);
1181 MOFEM_LOG("TIMER", Sev::verbose) << "TSSetUp <= done";
1182 MOFEM_LOG("TIMER", Sev::verbose) << "TSSolve";
1183 CHKERR TSSolve(solver, NULL);
1184 MOFEM_LOG("TIMER", Sev::verbose) << "TSSolve <= done";
1185
1187}
@ ROW
Definition: definitions.h:123
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
Definition: definitions.h:236
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:576
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMoFEM.cpp:219
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:483
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMoFEM.cpp:442
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:242
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:275
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
const FTensor::Tensor2< T, Dim, Dim > Vec
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition: DMMoFEM.hpp:1045
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition: SnesCtx.cpp:226
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
auto getDMSubData(DM dm)
Get sub problem data structure.
Definition: DMMoFEM.hpp:1061
auto createAOMappingIS(IS isapp, IS ispetsc)
Creates an application mapping using two index sets.
auto getDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition: DMMoFEM.hpp:1031
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
Definition: plastic.cpp:245
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
Definition: plastic.cpp:243
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
Definition: plastic.cpp:244
Class (Function) to enforce essential constrains on the left hand side diagonal.
Definition: Essential.hpp:47
Class (Function) to enforce essential constrains on the right hand side diagonal.
Definition: Essential.hpp:55
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
static std::array< int, 5 > activityData
Definition: PlasticOps.hpp:107
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field, SmartPetscObj< DM > sub_dm, SmartPetscObj< IS > field_split_it, SmartPetscObj< AO > ao_map)
Create data structure for handling Schur complement.
Definition: plastic.cpp:1474

Friends And Related Function Documentation

◆ TSPrePostProc

friend struct TSPrePostProc
friend

Definition at line 435 of file dynamic_first_order_con_law.cpp.

Member Data Documentation

◆ approxFunction

ApproxFieldFunction< FIELD_DIM > Example::approxFunction
staticprivate
Initial value:

Definition at line 61 of file approximaton.cpp.

◆ approxGradVals

boost::shared_ptr<MatrixDouble> Example::approxGradVals
private

Definition at line 63 of file radiation.cpp.

◆ approxVals

boost::shared_ptr<VectorDouble> Example::approxVals
private

Definition at line 62 of file radiation.cpp.

◆ aveMaxMin

std::array< double, Example::BoundingBox::LAST_BB > Example::aveMaxMin
staticprivate
Examples
phase.cpp.

Definition at line 110 of file phase.cpp.

◆ base

FieldApproximationBase Example::base
private

◆ boundaryMarker

boost::shared_ptr< std::vector< unsigned char > > Example::boundaryMarker
private
Examples
helmholtz.cpp, and phase.cpp.

Definition at line 42 of file helmholtz.cpp.

◆ commonDataPtr

boost::shared_ptr< CommonData > Example::commonDataPtr
private

Definition at line 42 of file integration.cpp.

◆ doEvalField

PetscBool Example::doEvalField
private

Definition at line 94 of file elastic.cpp.

◆ domianLhsFEPtr

boost::shared_ptr<FEMethod> Example::domianLhsFEPtr
private
Examples
shallow_wave.cpp.

Definition at line 355 of file shallow_wave.cpp.

◆ domianRhsFEPtr

boost::shared_ptr<FEMethod> Example::domianRhsFEPtr
private
Examples
shallow_wave.cpp.

Definition at line 356 of file shallow_wave.cpp.

◆ ePS

SmartPetscObj<EPS> Example::ePS
private
Examples
eigen_elastic.cpp.

Definition at line 69 of file eigen_elastic.cpp.

◆ fieldEvalCoords

std::array<double, SPACE_DIM> Example::fieldEvalCoords
private

Definition at line 95 of file elastic.cpp.

◆ fieldEvalData

boost::shared_ptr<FieldEvaluatorInterface::SetPtsData> Example::fieldEvalData
private

Definition at line 96 of file elastic.cpp.

◆ focalIndex

int Example::focalIndex
staticprivate
Examples
phase.cpp.

Definition at line 112 of file phase.cpp.

◆ iI

std::vector< MatrixInt > Example::iI
staticprivate
Examples
phase.cpp.

Definition at line 109 of file phase.cpp.

◆ K

SmartPetscObj<Mat> Example::K
private
Examples
eigen_elastic.cpp.

Definition at line 68 of file eigen_elastic.cpp.

◆ M

SmartPetscObj<Mat> Example::M
private
Examples
eigen_elastic.cpp.

Definition at line 67 of file eigen_elastic.cpp.

◆ matDPtr

boost::shared_ptr<MatrixDouble> Example::matDPtr
private
Examples
eigen_elastic.cpp.

Definition at line 65 of file eigen_elastic.cpp.

◆ mField

MoFEM::Interface & Example::mField
private

◆ pinchNodes

Range Example::pinchNodes
private
Examples
heat_method.cpp.

Definition at line 54 of file heat_method.cpp.

◆ reactionFe

boost::shared_ptr<DomainEle> Example::reactionFe
private
Examples
plastic.cpp.

Definition at line 241 of file plastic.cpp.

◆ rigidBodyMotion

std::array<SmartPetscObj<Vec>, 6> Example::rigidBodyMotion
private
Examples
eigen_elastic.cpp.

Definition at line 71 of file eigen_elastic.cpp.

◆ rZ

std::vector< double > Example::rZ
staticprivate
Examples
phase.cpp.

Definition at line 108 of file phase.cpp.

◆ savitzkyGolayNormalisation

int Example::savitzkyGolayNormalisation
staticprivate
Examples
phase.cpp.

Definition at line 113 of file phase.cpp.

◆ savitzkyGolayWeights

const int * Example::savitzkyGolayWeights
staticprivate
Examples
phase.cpp.

Definition at line 114 of file phase.cpp.

◆ simpleInterface

Simple * Example::simpleInterface
private
Examples
heat_method.cpp, helmholtz.cpp, phase.cpp, and plot_base.cpp.

Definition at line 41 of file helmholtz.cpp.

◆ space

FieldSpace Example::space
private
Examples
plot_base.cpp.

Definition at line 69 of file plot_base.cpp.

◆ uXScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > Example::uXScatter
private
Examples
plastic.cpp.

Definition at line 243 of file plastic.cpp.

◆ uYScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > Example::uYScatter
private
Examples
plastic.cpp.

Definition at line 244 of file plastic.cpp.

◆ uZScatter

std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter> > Example::uZScatter
private
Examples
plastic.cpp.

Definition at line 245 of file plastic.cpp.

◆ vectorFieldPtr

boost::shared_ptr<MatrixDouble> Example::vectorFieldPtr
private

Definition at line 97 of file elastic.cpp.


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