v0.12.1
photon_diffusion.cpp
Go to the documentation of this file.
1 /**
2  * \file photon_diffusion.cpp
3  * \example photon_diffusion.cpp
4  *
5  **/
6 
7 /* This file is part of MoFEM.
8  * MoFEM is free software: you can redistribute it and/or modify it under
9  * the terms of the GNU Lesser General Public License as published by the
10  * Free Software Foundation, either version 3 of the License, or (at your
11  * option) any later version.
12  *
13  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
14  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  * License for more details.
17  *
18  * You should have received a copy of the GNU Lesser General Public
19  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
20 
21 #include <stdlib.h>
22 #include <cmath>
23 #include <BasicFiniteElements.hpp>
24 
25 using namespace MoFEM;
26 
27 static char help[] = "...\n\n";
28 
29 template <int DIM> struct ElementsAndOps {};
30 
31 //! [Define dimension]
32 constexpr int SPACE_DIM = 3; //< Space dimension of problem, mesh
33 //! [Define dimension]
34 
40 
51 
58 
59 double n = 1.44; ///< refractive index of diffusive medium
60 double c = 30.; ///< speed of light (cm/ns)
61 double v = c / n; ///< phase velocity of light in medium (cm/ns)
62 double mu_a = 0.09; ///< absorption coefficient (cm^-1)
63 double mu_sp = 16.5; ///< scattering coefficient (cm^-1)
64 double flux = 1e3; ///< impulse magnitude
65 double duration = 0.05; ///< impulse duration (ns)
66 
67 PetscBool from_initial = PETSC_TRUE;
68 PetscBool output_volume = PETSC_FALSE;
69 
70 int order = 3;
72 
73 double A = 3.0;
74 
75 double h = 0.5 / A; ///< convective heat coefficient
76 double D = 1. / (3. * (mu_a + mu_sp));
77 double inv_v = 1. / v;
78 
79 /**
80  * @brief Monitor solution
81  *
82  * This functions is called by TS solver at the end of each step. It is used
83  * to output results to the hard drive.
84  */
85 struct Monitor : public FEMethod {
86 
87  Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc,
88  boost::shared_ptr<PostProcFaceOnRefinedMesh> skin_post_proc)
89  : dM(dm), postProc(post_proc), skinPostProc(skin_post_proc){};
90 
91  MoFEMErrorCode preProcess() { return 0; }
92  MoFEMErrorCode operator()() { return 0; }
93 
96  if (ts_step % save_every_nth_step == 0) {
97  if (output_volume) {
98  CHKERR DMoFEMLoopFiniteElements(dM, "dFE", postProc);
99  CHKERR postProc->writeFile(
100  "out_volume_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
101  }
102  if (skinPostProc) {
103  CHKERR DMoFEMLoopFiniteElements(dM, "CAMERA_FE", skinPostProc);
104  CHKERR skinPostProc->writeFile(
105  "out_camera_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
106  }
107  }
109  }
110 
111 private:
113  boost::shared_ptr<PostProcEle> postProc;
114  boost::shared_ptr<PostProcFaceOnRefinedMesh> skinPostProc;
115 };
116 
117 struct PhotonDiffusion {
118 public:
120 
121  // Declaration of the main function to run analysis
123 
124 private:
125  // Declaration of other main functions called in runProgram()
134 
135  // Main interfaces
136  MoFEM::Interface &mField;
137 
138  // Object to mark boundary entities for the assembling of domain elements
139  boost::shared_ptr<std::vector<unsigned char>> boundaryMarker;
140 
141  boost::shared_ptr<FEMethod> domianLhsFEPtr;
142  boost::shared_ptr<FEMethod> boundaryLhsFEPtr;
143  boost::shared_ptr<FEMethod> boundaryRhsFEPtr;
144 };
145 
146 PhotonDiffusion::PhotonDiffusion(MoFEM::Interface &m_field) : mField(m_field) {}
147 
150 
151  auto *simple = mField.getInterface<Simple>();
153  CHKERR simple->getOptions();
154  CHKERR simple->loadFile();
155 
157 }
158 
161 
162  auto *simple = mField.getInterface<Simple>();
163  CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 1);
164  CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, 1);
165 
166  CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-flux", &flux, PETSC_NULL);
167  CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-duration", &duration,
168  PETSC_NULL);
169  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-from_initial", &from_initial,
170  PETSC_NULL);
171  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-output_volume", &output_volume,
172  PETSC_NULL);
173 
174  MOFEM_LOG("PHOTON", Sev::inform) << "Refractive index: " << n;
175  MOFEM_LOG("PHOTON", Sev::inform) << "Speed of light (cm/ns): " << c;
176  MOFEM_LOG("PHOTON", Sev::inform) << "Phase velocity in medium (cm/ns): " << v;
177  MOFEM_LOG("PHOTON", Sev::inform)
178  << "Absorption coefficient (cm^-1): " << mu_a;
179  MOFEM_LOG("PHOTON", Sev::inform)
180  << "Scattering coefficient (cm^-1): " << mu_sp;
181  MOFEM_LOG("PHOTON", Sev::inform) << "Impulse magnitude: " << flux;
182  MOFEM_LOG("PHOTON", Sev::inform) << "Impulse duration (ns): " << duration;
183 
184  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
185  CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-save_step", &save_every_nth_step,
186  PETSC_NULL);
187 
188  MOFEM_LOG("PHOTON", Sev::inform) << "Approximation order: " << order;
189  MOFEM_LOG("PHOTON", Sev::inform) << "Save step: " << save_every_nth_step;
190 
191  CHKERR simple->setFieldOrder("U", order);
192 
193  auto set_camera_skin_fe = [&]() {
195  auto meshset_mng = mField.getInterface<MeshsetsManager>();
196 
197  Range camera_surface;
198  const std::string block_name = "CAM";
199  bool add_fe = false;
200 
202  if (bit->getName().compare(0, block_name.size(), block_name) == 0) {
203  MOFEM_LOG("PHOTON", Sev::inform) << "Found CAM block";
204  CHKERR mField.get_moab().get_entities_by_dimension(
205  bit->getMeshset(), 2, camera_surface, true);
206  add_fe = true;
207  }
208  }
209 
210  MOFEM_LOG("PHOTON", Sev::noisy) << "CAM block entities:\n"
211  << camera_surface;
212 
213  if (add_fe) {
214  CHKERR mField.add_finite_element("CAMERA_FE");
217  "CAMERA_FE");
218  }
220  };
221 
222  auto my_simple_set_up = [&]() {
224  CHKERR simple->defineFiniteElements();
225  CHKERR simple->defineProblem(PETSC_TRUE);
226  CHKERR simple->buildFields();
227  CHKERR simple->buildFiniteElements();
228 
229  if (mField.check_finite_element("CAMERA_FE")) {
230  CHKERR mField.build_finite_elements("CAMERA_FE");
231  CHKERR DMMoFEMAddElement(simple->getDM(), "CAMERA_FE");
232  }
233 
234  CHKERR simple->buildProblem();
236  };
237 
238  CHKERR set_camera_skin_fe();
239  CHKERR my_simple_set_up();
240 
242 }
243 
246 
247  auto integration_rule = [](int o_row, int o_col, int approx_order) {
248  return 2 * approx_order;
249  };
250 
251  auto *pipeline_mng = mField.getInterface<PipelineManager>();
252  CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
253  CHKERR pipeline_mng->setDomainLhsIntegrationRule(integration_rule);
254  CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(integration_rule);
255  CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
256 
258 }
259 
262 
264 }
265 
268  auto bc_mng = mField.getInterface<BcManager>();
269  auto *simple = mField.getInterface<Simple>();
270  CHKERR bc_mng->pushMarkDOFsOnEntities(simple->getProblemName(), "MIX", "U", 0,
271  0, false);
272  CHKERR bc_mng->pushMarkDOFsOnEntities(simple->getProblemName(), "SPOT", "U",
273  0, 0, false);
275 }
276 
279 
280  auto bc_mng = mField.getInterface<BcManager>();
281  auto &bc_map = bc_mng->getBcMapByBlockName();
282 
283  auto add_domain_base_ops = [&](auto &pipeline) {
284  auto jac_ptr = boost::make_shared<MatrixDouble>();
285  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
286  auto det_ptr = boost::make_shared<VectorDouble>();
287  pipeline.push_back(new OpCalculateHOJacVolume(jac_ptr));
288  pipeline.push_back(new OpInvertMatrix<3>(jac_ptr, det_ptr, inv_jac_ptr));
289  pipeline.push_back(new OpSetHOInvJacToScalarBases(H1, inv_jac_ptr));
290  pipeline.push_back(new OpSetHOWeights(det_ptr));
291  };
292 
293  auto add_domain_lhs_ops = [&](auto &pipeline) {
294  pipeline.push_back(new OpDomainGradGrad(
295  "U", "U", [](double, double, double) -> double { return D; }));
296  auto get_mass_coefficient = [&](const double, const double, const double) {
297  return inv_v * domianLhsFEPtr->ts_a + mu_a;
298  };
299  pipeline.push_back(new OpDomainMass("U", "U", get_mass_coefficient));
300  };
301 
302  auto add_domain_rhs_ops = [&](auto &pipeline) {
303  auto grad_u_at_gauss_pts = boost::make_shared<MatrixDouble>();
304  auto u_at_gauss_pts = boost::make_shared<VectorDouble>();
305  auto dot_u_at_gauss_pts = boost::make_shared<VectorDouble>();
306  pipeline.push_back(new OpCalculateScalarFieldGradient<SPACE_DIM>(
307  "U", grad_u_at_gauss_pts));
308  pipeline.push_back(new OpCalculateScalarFieldValues("U", u_at_gauss_pts));
309  pipeline.push_back(
310  new OpCalculateScalarFieldValuesDot("U", dot_u_at_gauss_pts));
311  pipeline.push_back(new OpDomainGradTimesVec(
312  "U", grad_u_at_gauss_pts,
313  [](double, double, double) -> double { return D; }));
314  pipeline.push_back(new OpDomainTimesScalarField(
315  "U", dot_u_at_gauss_pts,
316  [](const double, const double, const double) { return inv_v; }));
317  pipeline.push_back(new OpDomainTimesScalarField(
318  "U", u_at_gauss_pts,
319  [](const double, const double, const double) { return mu_a; }));
320  auto source_term = [&](const double, const double, const double) {
321  return 0;
322  };
323  pipeline.push_back(new OpDomainSource("U", source_term));
324  };
325 
326  auto add_boundary_base_ops = [&](auto &pipeline) {
327  pipeline.push_back(new OpSetHOWeightsOnFace());
328  };
329 
330  auto add_lhs_base_ops = [&](auto &pipeline) {
331  for (auto b : bc_map) {
332  if (std::regex_match(b.first, std::regex("(.*)EXT(.*)"))) {
333  pipeline.push_back(new OpBoundaryMass(
334  "U", "U",
335 
336  [](const double, const double, const double) { return h; },
337 
338  b.second->getBcEntsPtr()));
339  }
340  }
341  };
342 
343  auto add_rhs_base_ops = [&](auto &pipeline) {
344  auto u_at_gauss_pts = boost::make_shared<VectorDouble>();
345  pipeline.push_back(new OpCalculateScalarFieldValues("U", u_at_gauss_pts));
346  for (auto b : bc_map) {
347  if (std::regex_match(b.first, std::regex("(.*)EXT(.*)"))) {
348  pipeline.push_back(new OpBoundaryTimeScalarField(
349  "U", u_at_gauss_pts,
350 
351  [](const double, const double, const double) { return h; },
352 
353  b.second->getBcEntsPtr()));
354  }
355  }
356  };
357 
358  auto pipeline_mng = mField.getInterface<PipelineManager>();
359  add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
360  add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
361  add_domain_lhs_ops(pipeline_mng->getOpDomainLhsPipeline());
362  add_domain_rhs_ops(pipeline_mng->getOpDomainRhsPipeline());
363 
364  add_boundary_base_ops(pipeline_mng->getOpBoundaryLhsPipeline());
365  add_boundary_base_ops(pipeline_mng->getOpBoundaryRhsPipeline());
366  add_lhs_base_ops(pipeline_mng->getOpBoundaryLhsPipeline());
367  add_rhs_base_ops(pipeline_mng->getOpBoundaryRhsPipeline());
368 
369  domianLhsFEPtr = pipeline_mng->getDomainLhsFE();
370  boundaryLhsFEPtr = pipeline_mng->getBoundaryLhsFE();
371  boundaryRhsFEPtr = pipeline_mng->getBoundaryRhsFE();
372 
374 }
375 
378 
379  auto *simple = mField.getInterface<Simple>();
380  auto *pipeline_mng = mField.getInterface<PipelineManager>();
381 
382  auto create_post_process_element = [&]() {
383  auto post_froc_fe = boost::make_shared<PostProcEle>(mField);
384  post_froc_fe->generateReferenceElementMesh();
385  post_froc_fe->addFieldValuesPostProc("U");
386  post_froc_fe->addFieldValuesGradientPostProc("U");
387  return post_froc_fe;
388  };
389 
390  auto create_post_process_camera_element = [&]() {
391  if (mField.check_finite_element("CAMERA_FE")) {
392  auto post_proc_skin =
393  boost::make_shared<PostProcFaceOnRefinedMesh>(mField);
394  post_proc_skin->generateReferenceElementMesh();
395  CHKERR post_proc_skin->addFieldValuesPostProc("U");
396  CHKERR post_proc_skin->addFieldValuesGradientPostProcOnSkin(
397  "U", simple->getDomainFEName());
398  return post_proc_skin;
399  } else {
400  return boost::shared_ptr<PostProcFaceOnRefinedMesh>();
401  }
402  };
403 
404  auto set_time_monitor = [&](auto dm, auto solver) {
406  boost::shared_ptr<Monitor> monitor_ptr(
407  new Monitor(dm, create_post_process_element(),
408  create_post_process_camera_element()));
409  boost::shared_ptr<ForcesAndSourcesCore> null;
410  CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
411  monitor_ptr, null, null);
413  };
414 
415  auto dm = simple->getDM();
416  auto D = smartCreateDMVector(dm);
417 
418  if (from_initial) {
419 
420  MOFEM_LOG("PHOTON", Sev::inform)
421  << "reading vector in binary from vector.dat ...";
422  PetscViewer viewer;
423  PetscViewerBinaryOpen(PETSC_COMM_WORLD, "initial_vector.dat",
424  FILE_MODE_READ, &viewer);
425  VecLoad(D, viewer);
426 
427  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
428  }
429 
430  auto solver = pipeline_mng->createTS();
431 
432  CHKERR TSSetSolution(solver, D);
433  CHKERR set_time_monitor(dm, solver);
434  CHKERR TSSetSolution(solver, D);
435  CHKERR TSSetFromOptions(solver);
436  CHKERR TSSetUp(solver);
437  CHKERR TSSolve(solver, NULL);
438 
439  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
440  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
441  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
442 
444 }
445 
448 
449  // Processes to set output results are integrated in solveSystem()
450 
452 }
453 
456 
457  CHKERR readMesh();
465 
467 }
468 
469 int main(int argc, char *argv[]) {
470 
471  // Initialisation of MoFEM/PETSc and MOAB data structures
472  const char param_file[] = "param_file.petsc";
473  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
474 
475  // Add logging channel for example
476  auto core_log = logging::core::get();
477  core_log->add_sink(
478  LogManager::createSink(LogManager::getStrmWorld(), "PHOTON"));
479  LogManager::setLog("PHOTON");
480  MOFEM_LOG_TAG("PHOTON", "photon_diffusion")
481 
482  // Error handling
483  try {
484  // Register MoFEM discrete manager in PETSc
485  DMType dm_name = "DMMOFEM";
486  CHKERR DMRegister_MoFEM(dm_name);
487 
488  // Create MOAB instance
489  moab::Core mb_instance; // mesh database
490  moab::Interface &moab = mb_instance; // mesh database interface
491 
492  // Create MoFEM instance
493  MoFEM::Core core(moab); // finite element database
494  MoFEM::Interface &m_field = core; // finite element interface
495 
496  // Run the main analysis
497  PhotonDiffusion heat_problem(m_field);
498  CHKERR heat_problem.runProgram();
499  }
500  CATCH_ERRORS;
501 
502  // Finish work: cleaning memory, getting statistics, etc.
504 
505  return 0;
506 }
std::string param_file
ForcesAndSourcesCore::UserDataOperator UserDataOperator
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, FIELD_DIM > OpDomainMass
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ 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
@ BLOCKSET
Definition: definitions.h:161
#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
auto integration_rule
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:956
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:458
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:478
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:541
virtual MoFEMErrorCode add_ents_to_finite_element_by_dim(const EntityHandle entities, const int dim, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:309
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:340
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, 2 > OpDomainGradGrad
Definition: helmholtz.cpp:37
FormsIntegrators< EdgeEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpBoundarySource
Definition: helmholtz.cpp:43
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1 > OpDomainTimesScalarField
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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: DMMMoFEM.cpp:991
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
CoreTmp< 0 > Core
Definition: Core.hpp:1095
FaceElementForcesAndSourcesCore BoundaryEle
OpCalculateScalarFieldValuesFromPetscVecImpl< PetscData::CTX_SET_X_T > OpCalculateScalarFieldValuesDot
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1949
int main(int argc, char *argv[])
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1 > OpBoundaryTimeScalarField
double mu_sp
scattering coefficient (cm^-1)
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 1 > OpDomainMass
double A
static char help[]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1 > OpDomainTimesScalarField
double c
speed of light (cm/ns)
double v
phase velocity of light in medium (cm/ns)
double flux
impulse magnitude
double inv_v
constexpr int SPACE_DIM
[Define dimension]
double h
convective heat coefficient
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, SPACE_DIM > OpDomainGradGrad
double duration
impulse duration (ns)
int save_every_nth_step
double n
refractive index of diffusive medium
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 1 > OpBoundaryMass
int order
PetscBool from_initial
double D
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpDomainSource
PetscBool output_volume
double mu_a
absorption coefficient (cm^-1)
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpBoundaryMass
[Only used with Hencky/nonlinear material]
Definition: plastic.cpp:93
static constexpr int approx_order
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Simple interface for fast problem set-up.
Definition: BcManager.hpp:27
BcMapByBlockName & getBcMapByBlockName()
Get the bc map.
Definition: BcManager.hpp:131
virtual bool check_finite_element(const std::string &name) const =0
Check if finite element is in database.
virtual moab::Interface & get_moab()=0
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.
structure for User Loop Methods on finite elements
Interface for managing meshsets containing materials and boundary conditions.
Calculate jacobian on Hex or other volume which is not simplex.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get value at integration points for scalar field.
Set inverse jacobian to base functions.
Set inverse jacobian to base functions.
Modify integration weights on face to take in account higher-order geometry.
PipelineManager interface.
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
[Push operators to pipeline]
Monitor(SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc, boost::shared_ptr< PostProcFaceOnRefinedMesh > skin_post_proc)
boost::shared_ptr< PostProcFaceOnRefinedMesh > skinPostProc
MoFEMErrorCode postProcess()
function is run at the end of loop
MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode assembleSystem()
boost::shared_ptr< FEMethod > boundaryRhsFEPtr
MoFEMErrorCode solveSystem()
MoFEM::Interface & mField
MoFEMErrorCode readMesh()
MoFEMErrorCode outputResults()
MoFEMErrorCode initialCondition()
boost::shared_ptr< FEMethod > domianLhsFEPtr
PhotonDiffusion(MoFEM::Interface &m_field)
MoFEMErrorCode runProgram()
boost::shared_ptr< FEMethod > boundaryLhsFEPtr
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode setIntegrationRules()
MoFEMErrorCode setupProblem()
Post processing.