v0.14.0
higher_derivatives.cpp
Go to the documentation of this file.
1 /**
2  * \example higher_derivatives.cpp
3  *
4  * Testing higher derivatives of base functions
5  *
6  */
7 
8 
9 
10 #include <MoFEM.hpp>
11 
12 using namespace MoFEM;
13 
14 static char help[] = "...\n\n";
15 
16 constexpr char FIELD_NAME[] = "U";
17 constexpr int BASE_DIM = 1;
18 constexpr int FIELD_DIM = 1;
19 constexpr int SPACE_DIM = 2;
20 
21 template <int DIM> struct ElementsAndOps {};
22 
23 template <> struct ElementsAndOps<2> {
26 };
27 
28 using DomainEle = ElementsAndOps<SPACE_DIM>::DomainEle; ///< Finite elenent type
29 using DomainEleOp =
30  DomainEle::UserDataOperator; ///< Finire element operator type
31 using EntData = EntitiesFieldData::EntData; ///< Data on entities
32 
33 /**
34  * @brief Function to approximate
35  *
36  */
37 auto fun = [](const double x, const double y, const double z) {
38  return x * x + y * y + x * y + pow(x, 3) + pow(y, 3) + pow(x, 4) + pow(y, 4);
39 };
40 
41 /**
42  * @brief Function derivative
43  *
44  */
45 auto diff_fun = [](const double x, const double y, const double z) {
47  2 * x + y + 3 * pow(x, 2) + 4 * pow(x, 3),
48  2 * y + x + 3 * pow(y, 2) + 4 * pow(y, 3)};
49 };
50 
51 /**
52  * @brief Function second derivative
53  *
54  */
55 auto diff2_fun = [](const double x, const double y, const double z) {
57  2 + 6 * x + 12 * pow(x, 2), 1.,
58 
59  1., 2 + 6 * y + 12 * pow(y, 2)};
60 };
61 
62 /**
63  * @brief OPerator to integrate mass matrix for least square approximation
64  *
65  */
68 
69 /**
70  * @brief Operator to integrate the right hand side matrix for the problem
71  *
72  */
74  PETSC>::LinearForm<GAUSS>::OpSource<BASE_DIM, FIELD_DIM>;
75 
76 struct AtomTest {
77 
78  AtomTest(MoFEM::Interface &m_field) : mField(m_field) {}
79 
80  MoFEMErrorCode runProblem();
81 
82 private:
83  MoFEM::Interface &mField;
84  Simple *simpleInterface;
85 
86  MoFEMErrorCode readMesh();
87  MoFEMErrorCode setupProblem();
88  MoFEMErrorCode assembleSystem();
89  MoFEMErrorCode solveSystem();
90  MoFEMErrorCode checkResults();
91 
92  /**
93  * @brief Collected data use d by operator to evaluate errors for the test
94  *
95  */
96  struct CommonData {
97  boost::shared_ptr<MatrixDouble> invJacPtr;
98  boost::shared_ptr<VectorDouble> approxVals;
99  boost::shared_ptr<MatrixDouble> approxGradVals;
100  boost::shared_ptr<MatrixDouble> approxHessianVals;
101  SmartPetscObj<Vec> L2Vec;
102  };
103 
104  /**
105  * @brief Operator to evaluate errors
106  *
107  */
108  struct OpError;
109 };
110 
111 /**
112  * @brief Operator to evaluate errors
113  *
114  */
115 struct AtomTest::OpError : public DomainEleOp {
116  boost::shared_ptr<CommonData> commonDataPtr;
117  OpError(boost::shared_ptr<CommonData> &common_data_ptr)
118  : DomainEleOp(FIELD_NAME, OPROW), commonDataPtr(common_data_ptr) {
119  std::fill(doEntities.begin(), doEntities.end(), false);
120  doEntities[MBVERTEX] = true;
121  }
122  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
124 
125  const int nb_integration_pts = getGaussPts().size2();
126  auto t_w = getFTensor0IntegrationWeight(); // ger integration weights
127  auto t_val = getFTensor0FromVec(*(
128  commonDataPtr->approxVals)); // get function approximation at gauss pts
129  auto t_grad_val = getFTensor1FromMat<SPACE_DIM>(
130  *(commonDataPtr
131  ->approxGradVals)); // get gradient of approximation at gauss pts
132  auto t_hessian_val = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(
133  *(commonDataPtr)->approxHessianVals); // get hessian of approximation
134  // at integration pts
135 
136  auto t_inv_jac = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(
137  *(commonDataPtr->invJacPtr)); // get inverse of element jacobian
138  auto t_coords = getFTensor1CoordsAtGaussPts(); // get coordinates of
139  // integration points
140 
141  // Indices used for tensor operations
146 
147  const double volume = getMeasure(); // get finite element area
148 
149 
150  std::array<double, 3> error = {0, 0,
151  0}; // array for storing operator errors
152 
153  for (int gg = 0; gg != nb_integration_pts; ++gg) {
154 
155  const double alpha = t_w * volume;
156 
157  // Calculate errors
158 
159  double diff = t_val - fun(t_coords(0), t_coords(1), t_coords(2));
160  error[0] += alpha * pow(diff, 2);
161  auto t_diff_grad = diff_fun(t_coords(0), t_coords(1), t_coords(2));
162  t_diff_grad(i) -= t_grad_val(i);
163 
164  error[1] += alpha * t_diff_grad(i) *
165  t_diff_grad(i); // note push forward derivatives
166 
168  MOFEM_LOG("SELF", Sev::noisy) << "t_hessian_val " << t_hessian_push;
169 
170  // hessian expected to have symmetry
171  if (std::abs(t_hessian_val(0, 1) - t_hessian_val(1, 0)) >
172  std::numeric_limits<float>::epsilon()) {
173  MOFEM_LOG("SELF", Sev::error) << "t_hessian_val " << t_hessian_val;
174  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
175  "Hessian should be symmetric");
176  }
177 
178  auto t_diff_hessian = diff2_fun(t_coords(0), t_coords(1), t_coords(2));
179  t_diff_hessian(i, j) -= t_hessian_val(i, j);
180  error[2] = t_diff_hessian(i, j) * t_diff_hessian(i, j);
181 
182  // move data to next integration point
183  ++t_w;
184  ++t_val;
185  ++t_grad_val;
186  ++t_hessian_val;
187  ++t_inv_jac;
188  ++t_coords;
189  }
190 
191  // assemble error vector
192  std::array<int, 3> index = {0, 1, 2};
193  CHKERR VecSetValues(commonDataPtr->L2Vec, 3, index.data(), error.data(),
194  ADD_VALUES);
195 
197  }
198 };
199 
200 //! [Run programme]
203  CHKERR readMesh();
204  CHKERR setupProblem();
205  CHKERR assembleSystem();
206  CHKERR solveSystem();
207  CHKERR checkResults();
209 }
210 //! [Run programme]
211 
212 //! [Read mesh]
215 
216  CHKERR mField.getInterface(simpleInterface);
217  CHKERR simpleInterface->getOptions();
218  CHKERR simpleInterface->loadFile();
219 
221 }
222 //! [Read mesh]
223 
224 //! [Set up problem]
227  // Add field
228  CHKERR simpleInterface->addDomainField(FIELD_NAME, H1,
230  constexpr int order = 4;
231  CHKERR simpleInterface->setFieldOrder(FIELD_NAME, order);
232  CHKERR simpleInterface->setUp();
233 
235 }
236 //! [Set up problem]
237 
238 //! [Push operators to pipeline]
241 
242  auto rule = [](int, int, int p) -> int { return 2 * p; };
243 
244  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
245  CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule);
246  CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
247 
248  auto beta = [](const double, const double, const double) { return 1; };
249  pipeline_mng->getOpDomainLhsPipeline().push_back(
250  new OpDomainMass(FIELD_NAME, FIELD_NAME, beta));
251  pipeline_mng->getOpDomainRhsPipeline().push_back(
253 
255 }
256 //! [Push operators to pipeline]
257 
258 //! [Solve]
261  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
262 
263  MOFEM_LOG("WORLD", Sev::inform) << "Solve problem";
264 
265  auto solver = pipeline_mng->createKSP();
266  CHKERR KSPSetFromOptions(solver);
267  CHKERR KSPSetUp(solver);
268 
269  auto dm = simpleInterface->getDM();
270  auto D = createDMVector(dm);
271  auto F = vectorDuplicate(D);
272 
273  CHKERR KSPSolve(solver, F, D);
274  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
275  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
276  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
278 }
279 
280 //! [Check results]
283  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
284  pipeline_mng->getDomainLhsFE().reset();
285  pipeline_mng->getDomainRhsFE().reset();
286  pipeline_mng->getOpDomainRhsPipeline().clear();
287 
288  auto rule = [](int, int, int p) -> int { return 2 * p; };
289  CHKERR pipeline_mng->setDomainRhsIntegrationRule(
290  rule); // set integration rule
291 
292  // create data structures for operator
293  auto common_data_ptr = boost::make_shared<CommonData>();
294  common_data_ptr->L2Vec = createVectorMPI(
295  mField.get_comm(), (!mField.get_comm_rank()) ? 3 : 0, 3);
296  common_data_ptr->approxVals = boost::make_shared<VectorDouble>();
297  common_data_ptr->approxGradVals = boost::make_shared<MatrixDouble>();
298  common_data_ptr->approxHessianVals = boost::make_shared<MatrixDouble>();
299 
300  // create data strutires for evaluation of higher order derivatives
301  auto base_mass = boost::make_shared<MatrixDouble>();
302  auto data_l2 = boost::make_shared<EntitiesFieldData>(MBENTITYSET);
303  auto jac_ptr = boost::make_shared<MatrixDouble>();
304  auto det_ptr = boost::make_shared<VectorDouble>();
305  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
306  common_data_ptr->invJacPtr = inv_jac_ptr;
307 
308  // calculate jacobian at integration points
309  pipeline_mng->getOpDomainRhsPipeline().push_back(
310  new OpCalculateHOJacForFace(jac_ptr));
311  // calculate incerse of jacobian at integration points
312  pipeline_mng->getOpDomainRhsPipeline().push_back(
313  new OpInvertMatrix<SPACE_DIM>(jac_ptr, det_ptr, inv_jac_ptr));
314 
315  // calculate mass matrix to project derivatives
316  pipeline_mng->getOpDomainRhsPipeline().push_back(
317  new OpBaseDerivativesMass<BASE_DIM>(base_mass, data_l2,
319  // calculate second derivative of base functions, i.e. hessian
320  pipeline_mng->getOpDomainRhsPipeline().push_back(
321  new OpBaseDerivativesNext<BASE_DIM>(BaseDerivatives::SecondDerivative,
322  base_mass, data_l2,
324  // calculate third derivative
325  pipeline_mng->getOpDomainRhsPipeline().push_back(
326  new OpBaseDerivativesNext<BASE_DIM>(BaseDerivatives::ThirdDerivative,
327  base_mass, data_l2,
329  // calculate forth derivative
330  pipeline_mng->getOpDomainRhsPipeline().push_back(
331  new OpBaseDerivativesNext<BASE_DIM>(BaseDerivatives::ForthDerivative,
332  base_mass, data_l2,
334 
335  // push first base derivatives tp physical element shape
336  pipeline_mng->getOpDomainRhsPipeline().push_back(
337  new OpSetInvJacH1ForFace<1>(inv_jac_ptr));
338  // push second base derivatives tp physical element shape
339  pipeline_mng->getOpDomainRhsPipeline().push_back(
340  new OpSetInvJacH1ForFace<2>(inv_jac_ptr));
341 
342  // calculate value of function at integration points
343  pipeline_mng->getOpDomainRhsPipeline().push_back(
345  common_data_ptr->approxVals));
346 
347  // calculate gradient at integration points
348  pipeline_mng->getOpDomainRhsPipeline().push_back(
350  FIELD_NAME, common_data_ptr->approxGradVals));
351 
352 
353 
354  // calculate hessian at integration points
355  pipeline_mng->getOpDomainRhsPipeline().push_back(
357  FIELD_NAME, common_data_ptr->approxHessianVals));
358 
359  //FIXME: Note third and forth derivative is calculated but not properly tested
360 
361  // push operator to evaluate errors
362  pipeline_mng->getOpDomainRhsPipeline().push_back(
363  new OpError(common_data_ptr));
364 
365  // zero error vector and iterate over all elements on the mesh
366  CHKERR VecZeroEntries(common_data_ptr->L2Vec);
367  CHKERR pipeline_mng->loopFiniteElements();
368 
369  // assemble error vector
370  CHKERR VecAssemblyBegin(common_data_ptr->L2Vec);
371  CHKERR VecAssemblyEnd(common_data_ptr->L2Vec);
372 
373  // check if errors are small and print results
374 
375  constexpr double eps = 1e-8;
376 
377  const double *array;
378  CHKERR VecGetArrayRead(common_data_ptr->L2Vec, &array);
379  MOFEM_LOG_C("WORLD", Sev::inform,
380  "Error %6.4e Diff Error %6.4e Diff2 Error %6.4e\n",
381  std::sqrt(array[0]), std::sqrt(array[1]), std::sqrt(array[2]));
382  if (std::sqrt(array[0]) > eps)
383  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Wrong function value");
384  if (std::sqrt(array[1]) > eps)
385  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
386  "Wrong function first direcative");
387  if (std::sqrt(array[2]) > eps)
388  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
389  "Wrong function second direcative");
390 
391  CHKERR VecRestoreArrayRead(common_data_ptr->L2Vec, &array);
392 
394 }
395 //! [Check results]
396 
397 int main(int argc, char *argv[]) {
398 
399  // Initialisation of MoFEM/PETSc and MOAB data structures
400  MoFEM::Core::Initialize(&argc, &argv, NULL, help);
401 
402  try {
403 
404  //! [Register MoFEM discrete manager in PETSc]
405  DMType dm_name = "DMMOFEM";
406  CHKERR DMRegister_MoFEM(dm_name);
407  //! [Register MoFEM discrete manager in PETSc
408 
409  //! [Create MoAB]
410  moab::Core mb_instance; ///< mesh database
411  moab::Interface &moab = mb_instance; ///< mesh database interface
412  //! [Create MoAB]
413 
414  //! [Create MoFEM]
415  MoFEM::Core core(moab); ///< finite element database
416  MoFEM::Interface &m_field = core; ///< finite element database insterface
417  //! [Create MoFEM]
418 
419  //! [AtomTest]
420  AtomTest ex(m_field);
421  CHKERR ex.runProblem();
422  //! [AtomTest]
423  }
424  CATCH_ERRORS;
425 
427 }
AtomTest
Definition: child_and_parent.cpp:57
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
MoFEM::PipelineManager::getDomainRhsFE
boost::shared_ptr< FEMethod > & getDomainRhsFE()
Definition: PipelineManager.hpp:405
FIELD_DIM
constexpr int FIELD_DIM
Definition: higher_derivatives.cpp:18
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
OpError
Definition: initial_diffusion.cpp:61
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
MoFEM::PipelineManager::loopFiniteElements
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
Definition: PipelineManager.cpp:19
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::PipelineManager::getOpDomainRhsPipeline
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
Definition: PipelineManager.hpp:773
MoFEM::PETSC
@ PETSC
Definition: FormsIntegrators.hpp:105
MoFEM::PipelineManager
PipelineManager interface.
Definition: PipelineManager.hpp:24
MoFEM.hpp
MoFEM::DMoFEMMeshToLocalVector
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:523
MoFEM::OpBaseDerivativesNext
Definition: BaseDerivativesDataOperators.hpp:67
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
AtomTest::assembleSystem
MoFEMErrorCode assembleSystem()
[Push operators to pipeline]
Definition: child_and_parent.cpp:290
MoFEM::Simple
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
OpDomainMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, FIELD_DIM > OpDomainMass
Definition: child_and_parent.cpp:53
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
FTensor::Tensor2
Definition: Tensor2_value.hpp:16
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::OpCalculateScalarFieldGradient
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Definition: UserDataOperators.hpp:1293
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::OpCalculateScalarFieldHessian
Evaluate field gradient values for scalar field, i.e. gradient is tensor rank 1 (vector),...
Definition: UserDataOperators.hpp:1306
AtomTest::OpError
Operator to evaluate errors.
Definition: child_and_parent.cpp:82
AtomTest::OpError::OpError
OpError(boost::shared_ptr< CommonData > &common_data_ptr)
Definition: higher_derivatives.cpp:117
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::OpBaseDerivativesMass
Definition: BaseDerivativesDataOperators.hpp:35
MoFEM::createDMVector
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1099
diff2_fun
auto diff2_fun
Function second derivative.
Definition: higher_derivatives.cpp:55
BASE_DIM
constexpr int BASE_DIM
Definition: higher_derivatives.cpp:17
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator
default operator for TRI element
Definition: FaceElementForcesAndSourcesCore.hpp:94
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
double
convert.type
type
Definition: convert.py:64
MoFEM::PipelineManager::FaceEle
MoFEM::FaceElementForcesAndSourcesCore FaceEle
Definition: PipelineManager.hpp:35
MoFEM::FormsIntegrators::Assembly
Assembly methods.
Definition: FormsIntegrators.hpp:317
MoFEM::getFTensor0FromVec
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
MoFEM::OpCalculateScalarFieldValues
Get value at integration points for scalar field.
Definition: UserDataOperators.hpp:82
MoFEM::PipelineManager::createKSP
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
Definition: PipelineManager.cpp:49
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
MoFEM::FaceElementForcesAndSourcesCore
Face finite element.
Definition: FaceElementForcesAndSourcesCore.hpp:23
MoFEM::PipelineManager::setDomainRhsIntegrationRule
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
Definition: PipelineManager.hpp:530
AtomTest::AtomTest
AtomTest(MoFEM::Interface &m_field)
Definition: higher_derivatives.cpp:78
MoFEM::PipelineManager::getOpDomainLhsPipeline
boost::ptr_deque< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
Definition: PipelineManager.hpp:749
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
BiLinearForm
FIELD_NAME
constexpr char FIELD_NAME[]
Definition: higher_derivatives.cpp:16
FTensor::Index< 'i', 2 >
AtomTest::solveSystem
MoFEMErrorCode solveSystem()
[Push operators to pipeline]
Definition: child_and_parent.cpp:377
MoFEM::PipelineManager::setDomainLhsIntegrationRule
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Definition: PipelineManager.hpp:503
AtomTest::OpError::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Definition: higher_derivatives.cpp:122
diff_fun
auto diff_fun
Function derivative.
Definition: higher_derivatives.cpp:45
MoFEM::OpSetInvJacH1ForFace
Definition: UserDataOperators.hpp:2918
ElementsAndOps
Definition: child_and_parent.cpp:18
MoFEM::PipelineManager::getDomainLhsFE
boost::shared_ptr< FEMethod > & getDomainLhsFE()
Definition: PipelineManager.hpp:401
DomainEleOp
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
fun
auto fun
Function to approximate.
Definition: higher_derivatives.cpp:37
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::vectorDuplicate
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
Definition: PetscSmartObj.hpp:221
help
static char help[]
Definition: higher_derivatives.cpp:14
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
UserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Definition: HookeElement.hpp:75
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
eps
static const double eps
Definition: check_base_functions_derivatives_on_tet.cpp:11
OpDomainMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< BASE_DIM, FIELD_DIM > OpDomainMass
OPerator to integrate mass matrix for least square approximation.
Definition: higher_derivatives.cpp:67
CommonData
Definition: continuity_check_on_skeleton_with_simple_2d_for_h1.cpp:22
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
EntData
EntitiesFieldData::EntData EntData
Data on entities.
Definition: higher_derivatives.cpp:31
AtomTest::checkResults
MoFEMErrorCode checkResults()
[Check results]
Definition: dg_projection.cpp:213
MoFEM::createVectorMPI
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
Definition: PetscSmartObj.hpp:202
MoFEM::OpInvertMatrix
Definition: UserDataOperators.hpp:3249
MoFEM::OpCalculateHOJacForFace
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
Definition: HODataOperators.hpp:264
SPACE_DIM
constexpr int SPACE_DIM
Definition: higher_derivatives.cpp:19
AtomTest::readMesh
MoFEMErrorCode readMesh()
[Run programme]
Definition: child_and_parent.cpp:208
ReactionDiffusionEquation::D
const double D
diffusivity
Definition: reaction_diffusion.cpp:20
main
int main(int argc, char *argv[])
[Check results]
Definition: higher_derivatives.cpp:397
MOFEM_ATOM_TEST_INVALID
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
DomainEle
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Definition: child_and_parent.cpp:34
MoFEM::SmartPetscObj< Vec >
OpDomainSource
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
Definition: child_and_parent.cpp:55
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
OpDomainSource
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< BASE_DIM, FIELD_DIM > OpDomainSource
Operator to integrate the right hand side matrix for the problem.
Definition: higher_derivatives.cpp:74
convert.int
int
Definition: convert.py:64
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
AtomTest::setupProblem
MoFEMErrorCode setupProblem()
[Read mesh]
Definition: child_and_parent.cpp:268
AtomTest::runProblem
MoFEMErrorCode runProblem()
[Run programme]
Definition: child_and_parent.cpp:188
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
F
@ F
Definition: free_surface.cpp:394