v0.13.0
radiation.cpp
Go to the documentation of this file.
1 /**
2  * \file lesson6_radiation.cpp
3  * \example lesson6_radiation.cpp
4  *
5  * Using PipelineManager interface calculate the divergence of base functions,
6  * and integral of flux on the boundary. Since the h-div space is used, volume
7  * integral and boundary integral should give the same result.
8  */
9 
10 /* This file is part of MoFEM.
11  * MoFEM is free software: you can redistribute it and/or modify it under
12  * the terms of the GNU Lesser General Public License as published by the
13  * Free Software Foundation, either version 3 of the License, or (at your
14  * option) any later version.
15  *
16  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19  * License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23 
24 #include <MoFEM.hpp>
25 
26 using namespace MoFEM;
27 
28 static char help[] = "...\n\n";
29 
30 #include <BasicFiniteElements.hpp>
31 
42 
43 // Units
44 // Temperature: Kelvins
45 // Length: 1 km
46 // Time: 1 s
47 
48 constexpr double heat_conductivity = ((0.4+0.7)/2) * 1e3;
49 
50 constexpr double emissivity = 1;
51 constexpr double boltzmann_constant = 5.670374419e-2;
52 constexpr double Beta = emissivity * boltzmann_constant;
53 
54 constexpr double T_ambient = 2.7;
55 struct Example {
56 
57  Example(MoFEM::Interface &m_field) : mField(m_field) {}
58 
60 
61 private:
62  MoFEM::Interface &mField;
63 
64  static int integrationRule(int, int, int p_data) { return 2 * p_data; };
65 
70  MoFEMErrorCode kspSolve();
73 
74  boost::shared_ptr<VectorDouble> approxVals;
75  boost::shared_ptr<MatrixDouble> approxGradVals;
76 
77  struct OpRadiationLhs : public OpBase {
78 
79  private:
80  boost::shared_ptr<VectorDouble> approxVals;
81 
82  public:
83  OpRadiationLhs(boost::shared_ptr<VectorDouble> &approx_vals)
84  : OpBase("T", "T", OpBase::OPROWCOL), approxVals(approx_vals) {
85  this->sYmm = false;
86  }
87 
88  MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
89  };
90 
91  struct OpRadiationRhs : public OpBase {
92 
93  private:
94  boost::shared_ptr<VectorDouble> approxVals;
95 
96  public:
97  OpRadiationRhs(boost::shared_ptr<VectorDouble> &approx_vals)
98  : OpBase("T", "T", OpBase::OPROW), approxVals(approx_vals) {}
99 
100  MoFEMErrorCode iNtegrate(EntData &row_data);
101  };
102 
103  struct OpFluxRhs : public OpBase {
104 
105  private:
106  FTensor::Index<'i', 2> i; ///< summit Index
107 
108  public:
109  OpFluxRhs() : OpBase("T", "T", OpBase::OPROW) {}
110 
111  MoFEMErrorCode iNtegrate(EntData &row_data);
112  };
113 
115 
116  private:
117  boost::shared_ptr<VectorDouble> approxVals;
118  double &sumTemperature;
119  double &surfaceArea;
120 
121 
122  public:
124  boost::shared_ptr<VectorDouble> &approx_vals, double &sum_temp,
125  double &surf)
126  : EdgeEleOp("T", "T", OpBase::OPROW), approxVals(approx_vals),
127  sumTemperature(sum_temp), surfaceArea(surf) {}
128 
129  MoFEMErrorCode doWork(int side, EntityType type,
131 
132 
133  };
134 
135 
136 };
137 
140  CHKERR setupProblem();
141  CHKERR createCommonData();
142  CHKERR bC();
143  CHKERR OPs();
144  CHKERR kspSolve();
145  CHKERR postProcess();
146  CHKERR checkResults();
148 }
149 
150 //! [Set up problem]
153  Simple *simple = mField.getInterface<Simple>();
154  // Add field
155  CHKERR simple->addDomainField("T", H1, AINSWORTH_LEGENDRE_BASE, 1);
156  CHKERR simple->addBoundaryField("T", H1, AINSWORTH_LEGENDRE_BASE, 1);
157  constexpr int order = 3;
158  CHKERR simple->setFieldOrder("T", order);
159  CHKERR simple->setUp();
161 }
162 //! [Set up problem]
163 
164 //! [Create common data]
167  approxVals = boost::make_shared<VectorDouble>();
168  approxGradVals = boost::make_shared<MatrixDouble>();
170 }
171 //! [Create common data]
172 
173 //! [Boundary condition]
176  // Set initial values
177  auto set_initial_temperature = [&](VectorAdaptor &&field_data, double *xcoord,
178  double *ycoord, double *zcoord) {
180  field_data[0] = T_ambient;
182  };
183  FieldBlas *field_blas;
184  CHKERR mField.getInterface(field_blas);
185  CHKERR field_blas->setVertexDofs(set_initial_temperature, "T");
187 }
188 //! [Boundary condition]
189 
190 //! [Push operators to pipeline]
193  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
194  auto beta = [](const double r, const double, const double) {
195  return heat_conductivity * (2 * M_PI * r);
196  };
197 
198  auto det_ptr = boost::make_shared<VectorDouble>();
199  auto jac_ptr = boost::make_shared<MatrixDouble>();
200  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
201 
202  pipeline_mng->getOpDomainLhsPipeline().push_back(new OpSetHOWeightsOnFace());
203  pipeline_mng->getOpDomainLhsPipeline().push_back(
204  new OpCalculateHOJacForFace(jac_ptr));
205  pipeline_mng->getOpDomainLhsPipeline().push_back(
206  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
207  pipeline_mng->getOpDomainLhsPipeline().push_back(
208  new OpSetInvJacH1ForFace(inv_jac_ptr));
209  pipeline_mng->getOpDomainLhsPipeline().push_back(
210  new OpDomainGradGrad("T", "T", beta));
211  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integrationRule);
212 
213  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpSetHOWeightsOnFace());
214  pipeline_mng->getOpDomainRhsPipeline().push_back(
215  new OpCalculateHOJacForFace(jac_ptr));
216  pipeline_mng->getOpDomainRhsPipeline().push_back(
217  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
218  pipeline_mng->getOpDomainRhsPipeline().push_back(
219  new OpSetInvJacH1ForFace(inv_jac_ptr));
220  pipeline_mng->getOpDomainRhsPipeline().push_back(
221  new OpCalculateScalarFieldGradient<2>("T", approxGradVals));
222  pipeline_mng->getOpDomainRhsPipeline().push_back(
223  new OpDomainGradTimesVec("T", approxGradVals, beta));
224  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integrationRule);
225 
226  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
227  new OpCalculateScalarFieldValues("T", approxVals));
228  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
229  new OpRadiationRhs(approxVals));
230  pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpFluxRhs());
231  CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integrationRule);
232 
233  pipeline_mng->getOpBoundaryLhsPipeline().push_back(
234  new OpCalculateScalarFieldValues("T", approxVals));
235  pipeline_mng->getOpBoundaryLhsPipeline().push_back(
236  new OpRadiationLhs(approxVals));
237  CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(integrationRule);
239 }
240 //! [Push operators to pipeline]
241 
242 //! [Solve]
245  Simple *simple = mField.getInterface<Simple>();
246  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
247  auto ts = pipeline_mng->createTSIM();
248 
249  double ftime = 1;
250  CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
251  CHKERR TSSetFromOptions(ts);
252  CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
253 
254  auto T = smartCreateDMVector(simple->getDM());
255  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), T, INSERT_VALUES,
256  SCATTER_FORWARD);
257 
258  CHKERR TSSolve(ts, T);
259  CHKERR TSGetTime(ts, &ftime);
260 
261  PetscInt steps, snesfails, rejects, nonlinits, linits;
262  CHKERR TSGetTimeStepNumber(ts, &steps);
263  CHKERR TSGetSNESFailures(ts, &snesfails);
264  CHKERR TSGetStepRejections(ts, &rejects);
265  CHKERR TSGetSNESIterations(ts, &nonlinits);
266  CHKERR TSGetKSPIterations(ts, &linits);
267  MOFEM_LOG_C("EXAMPLE", Sev::inform,
268  "steps %d (%d rejected, %d SNES fails), ftime %g, nonlinits "
269  "%d, linits %d",
270  steps, rejects, snesfails, ftime, nonlinits, linits);
271 
273 }
274 //! [Solve]
275 
276 //! [Postprocess results]
279  PipelineManager *pipeline_mng = mField.getInterface<PipelineManager>();
280  pipeline_mng->getDomainLhsFE().reset();
281  pipeline_mng->getBoundaryLhsFE().reset();
282  pipeline_mng->getBoundaryRhsFE().reset();
283  auto post_proc_fe = boost::make_shared<PostProcFaceOnRefinedMesh>(mField);
284  post_proc_fe->generateReferenceElementMesh();
285  post_proc_fe->addFieldValuesPostProc("T");
286  pipeline_mng->getDomainRhsFE() = post_proc_fe;
287 
288  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
289  new OpCalculateScalarFieldValues("T", approxVals));
290 
291  double sum_temperature;
292  double surface_area;
293  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
294  new OpCalcSurfaceAverageTemperature(approxVals, sum_temperature,
295  surface_area));
296  auto calc_surfcae_area_op = pipeline_mng->getOpBoundaryRhsPipeline().back();
297 
298  sum_temperature = 0;
299  surface_area = 0;
300  CHKERR pipeline_mng->loopFiniteElements();
301  CHKERR post_proc_fe->writeFile("out_radiation.h5m");
302 
303  MOFEM_LOG_C("EXAMPLE", Sev::inform, "Surface area %3.4e [km]", surface_area);
304  MOFEM_LOG_C("EXAMPLE", Sev::inform,
305  "Average subsurface temperatute %3.4e [K]",
306  sum_temperature / surface_area);
307 
309 }
310 //! [Postprocess results]
311 
312 //! [Check results]
314 //! [Check results]
315 
316 int main(int argc, char *argv[]) {
317 
318  // Initialisation of MoFEM/PETSc and MOAB data structures
319  const char param_file[] = "param_file.petsc";
320  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
321 
322  // Add logging channel for example
323  auto core_log = logging::core::get();
324  core_log->add_sink(
326  LogManager::setLog("EXAMPLE");
327  MOFEM_LOG_TAG("EXAMPLE", "example");
328 
329  try {
330 
331  //! [Register MoFEM discrete manager in PETSc]
332  DMType dm_name = "DMMOFEM";
333  CHKERR DMRegister_MoFEM(dm_name);
334  //! [Register MoFEM discrete manager in PETSc
335 
336  //! [Create MoAB]
337  moab::Core mb_instance; ///< mesh database
338  moab::Interface &moab = mb_instance; ///< mesh database interface
339  //! [Create MoAB]
340 
341  //! [Create MoFEM]
342  MoFEM::Core core(moab); ///< finite element database
343  MoFEM::Interface &m_field = core; ///< finite element database insterface
344  //! [Create MoFEM]
345 
346  //! [Load mesh]
347  Simple *simple = m_field.getInterface<Simple>();
348  CHKERR simple->getOptions();
349  CHKERR simple->loadFile("");
350  //! [Load mesh]
351 
352  //! [Example]
353  Example ex(m_field);
354  CHKERR ex.runProblem();
355  //! [Example]
356  }
357  CATCH_ERRORS;
358 
360 }
361 
362 //! [Radiation Lhs]
364  EntData &col_data) {
366  // get element volume
367  const double vol = this->getMeasure();
368  // get integration weights
369  auto t_w = getFTensor0IntegrationWeight();
370  // get base function gradient on rows
371  auto t_row_base = row_data.getFTensor0N();
372  // gat temperature at integration points
373  auto t_val = getFTensor0FromVec(*(approxVals));
374  // get coordinate at integration points
375  auto t_coords = getFTensor1CoordsAtGaussPts();
376 
377  // loop over integration points
378  for (int gg = 0; gg != nbIntegrationPts; gg++) {
379 
380  // Cylinder radius
381  const double r_cylinder = t_coords(0);
382 
383  // take into account Jacobean
384  const double alpha = t_w * vol * Beta * (2 * M_PI * r_cylinder);
385  // loop over rows base functions
386  for (int rr = 0; rr != nbRows; ++rr) {
387  auto t_col_base = col_data.getFTensor0N(gg, 0);
388  // loop over columns
389  for (int cc = 0; cc != nbCols; cc++) {
390  if (std::abs(t_coords(0)) > std::numeric_limits<double>::epsilon()) {
391  locMat(rr, cc) += alpha * t_row_base * t_col_base * 4 * pow(t_val, 3);
392  }
393  ++t_col_base;
394  }
395  ++t_row_base;
396  }
397 
398  ++t_val;
399  ++t_coords;
400  ++t_w; // move to another integration weight
401  }
403 }
404 //! [Radiation Lhs]
405 
406 //! [Radiation Lhs]
409  // get element volume
410  const double vol = getMeasure();
411  // get integration weights
412  auto t_w = getFTensor0IntegrationWeight();
413  // get base function gradient on rows
414  auto t_row_base = row_data.getFTensor0N();
415  // gat temperature at integration points
416  auto t_val = getFTensor0FromVec(*(approxVals));
417  // get coordinate at integration points
418  auto t_coords = getFTensor1CoordsAtGaussPts();
419 
420  // loop over integration points
421  for (int gg = 0; gg != nbIntegrationPts; gg++) {
422 
423  // Cylinder radius
424  const double r_cylinder = t_coords(0);
425 
426  // take into account Jacobean
427  const double alpha = t_w * vol * Beta * (2 * M_PI * r_cylinder);
428  // loop over rows base functions
429  for (int rr = 0; rr != nbRows; ++rr) {
430  if (std::abs(t_coords(0)) > std::numeric_limits<double>::epsilon()) {
431  locF[rr] += alpha * t_row_base * (pow(t_val, 4) - pow(T_ambient, 4));
432  }
433  ++t_row_base;
434  }
435  ++t_coords;
436  ++t_val;
437  ++t_w; // move to another integration weight
438  }
439 
441 }
442 //! [Radiation Lhs]
443 
444 //! [Flux Rhs]
447  // get element volume
448  const double vol = getMeasure();
449  // get integration weights
450  auto t_w = getFTensor0IntegrationWeight();
451  // get base function gradient on rows
452  auto t_row_base = row_data.getFTensor0N();
453  // get coordinate at integration points
454  auto t_coords = getFTensor1CoordsAtGaussPts();
455  // // get time
456  const double time = getFEMethod()->ts_t;
457 
458  // Look to https://doi.org/10.1016/j.icarus.2014.12.028s
459  constexpr double flux_p = -0.03e6;
460  constexpr double flux_c = -0.23e6;
461 
462  // loop over integration points
463  for (int gg = 0; gg != nbIntegrationPts; gg++) {
464 
465  // Cylinder radius
466  const double r_cylinder = t_coords(0);
467 
468  const double r = std::sqrt(t_coords(i) * t_coords(i));
469  const double s = std::abs(t_coords(1)) / r;
470 
471  // take into account Jacobean
472  const double alpha = t_w * vol * (2 * M_PI * r_cylinder);
473  // loop over rows base functions
474  for (int rr = 0; rr != nbRows; ++rr) {
475  locF[rr] += alpha * t_row_base * (s * flux_p + flux_c) * time;
476  ++t_row_base;
477  }
478  ++t_coords;
479  ++t_w; // move to another integration weight
480  }
481 
483 }
484 //! [Flux Rhs]
485 
486 //! [Ave Temp]
488  int side, EntityType type, EntitiesFieldData::EntData &data) {
489 
491  if (type == MBVERTEX) {
492  // get element volume
493  const double vol = getMeasure();
494  // get integration weights
495  auto t_w = getFTensor0IntegrationWeight();
496  // gat temperature at integration points
497  auto t_val = getFTensor0FromVec(*(approxVals));
498  // get coordinate at integration points
499  auto t_coords = getFTensor1CoordsAtGaussPts();
500  // number of integration pts
501  size_t nb_integration_pts = getGaussPts().size2();
502 
503  // loop over integration points
504  for (auto gg = 0; gg != nb_integration_pts; ++gg) {
505 
506  // Cylinder radius
507  const double r_cylinder = t_coords(0);
508 
509  // take into account Jacobean
510  const double alpha = t_w * vol * (2 * M_PI * r_cylinder);
511 
512  sumTemperature += alpha * t_val;
513  surfaceArea += alpha;
514 
515  ++t_coords;
516  ++t_val;
517  ++t_w; // move to another integration weight
518  }
519  }
521 }
522 
523 //! [Ave Temp]
std::string param_file
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
constexpr double alpha
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:956
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:481
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::ptr_vector< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
SmartPetscObj< TS > createTSIM(SmartPetscObj< DM > dm=nullptr)
Create TS (time) implicit solver.
boost::ptr_vector< UserDataOperator > & getOpBoundaryLhsPipeline()
Get the Op Boundary Lhs Pipeline object.
boost::ptr_vector< UserDataOperator > & getOpBoundaryRhsPipeline()
Get the Op Boundary Rhs Pipeline object.
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
FTensor::Index< 'i', SPACE_DIM > i
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, 2 > OpDomainGradGrad
Definition: helmholtz.cpp:37
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
const double T
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:126
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
CoreTmp< 0 > Core
Definition: Core.hpp:1096
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:149
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
const double r
rate factor
int main(int argc, char *argv[])
[Check results]
Definition: radiation.cpp:316
EntitiesFieldData::EntData EntData
Definition: radiation.cpp:36
static char help[]
Definition: radiation.cpp:28
constexpr double Beta
Definition: radiation.cpp:52
constexpr double boltzmann_constant
Definition: radiation.cpp:51
constexpr double T_ambient
Definition: radiation.cpp:54
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, 2 > OpDomainGradTimesVec
Definition: radiation.cpp:40
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, 2 > OpDomainGradGrad
Definition: radiation.cpp:38
constexpr double emissivity
Definition: radiation.cpp:50
constexpr double heat_conductivity
Definition: radiation.cpp:48
OpCalcSurfaceAverageTemperature(boost::shared_ptr< VectorDouble > &approx_vals, double &sum_temp, double &surf)
Definition: radiation.cpp:123
boost::shared_ptr< VectorDouble > approxVals
Definition: radiation.cpp:117
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
[Flux Rhs]
Definition: radiation.cpp:487
FTensor::Index< 'i', 2 > i
summit Index
Definition: radiation.cpp:106
MoFEMErrorCode iNtegrate(EntData &row_data)
[Radiation Lhs]
Definition: radiation.cpp:445
boost::shared_ptr< VectorDouble > approxVals
Definition: radiation.cpp:80
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
[Radiation Lhs]
Definition: radiation.cpp:363
OpRadiationLhs(boost::shared_ptr< VectorDouble > &approx_vals)
Definition: radiation.cpp:83
OpRadiationRhs(boost::shared_ptr< VectorDouble > &approx_vals)
Definition: radiation.cpp:97
boost::shared_ptr< VectorDouble > approxVals
Definition: radiation.cpp:94
MoFEMErrorCode iNtegrate(EntData &row_data)
[Radiation Lhs]
Definition: radiation.cpp:407
[Example]
static int integrationRule(int, int, int p_data)
Definition: radiation.cpp:64
boost::shared_ptr< VectorDouble > approxVals
Definition: radiation.cpp:74
MoFEMErrorCode kspSolve()
[Push operators to pipeline]
Definition: radiation.cpp:243
MoFEMErrorCode checkResults()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
Definition: radiation.cpp:57
MoFEMErrorCode OPs()
MoFEMErrorCode runProblem()
boost::shared_ptr< MatrixDouble > approxGradVals
Definition: radiation.cpp:75
MoFEMErrorCode postProcess()
MoFEMErrorCode setupProblem()
MoFEMErrorCode bC()
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
Basic algebra on fields.
Definition: FieldBlas.hpp:31
MoFEMErrorCode setVertexDofs(VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
Set DOFs on vertices using user function.
Definition: FieldBlas.cpp:333
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get value at integration points for scalar field.
Modify integration weights on face to take in account higher-order geometry.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
MoFEM::FaceElementForcesAndSourcesCore FaceEle
boost::shared_ptr< FEMethod > & getBoundaryLhsFE()
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
VolumeElementForcesAndSourcesCoreBase::UserDataOperator UserDataOperator