16static char help[] = 
"...\n\n";
 
  111        boost::shared_ptr<VectorDouble> &approx_vals, 
double &sum_temp,
 
 
 
  140  constexpr int order = 3;
 
  150  approxVals = boost::make_shared<VectorDouble>();
 
  160  auto set_initial_temperature = [&](
VectorAdaptor &&field_data, 
double *xcoord,
 
  161                                     double *ycoord, 
double *zcoord) {
 
  181  auto det_ptr = boost::make_shared<VectorDouble>();
 
  182  auto jac_ptr = boost::make_shared<MatrixDouble>();
 
  183  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
 
  237  CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
 
  238  CHKERR TSSetFromOptions(ts);
 
  239  CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
 
  246  CHKERR TSGetTime(ts, &ftime);
 
  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);
 
  255              "steps %d (%d rejected, %d SNES fails), ftime %g, nonlinits " 
  257              steps, rejects, snesfails, ftime, nonlinits, linits);
 
 
  271      boost::make_shared<PostProcBrokenMeshInMoab<DomainEle>>(
mField);
 
  273  auto t_ptr = boost::make_shared<VectorDouble>();
 
  274  post_proc_fe->getOpPtrVector().push_back(
 
  279  post_proc_fe->getOpPtrVector().push_back(
 
  281      new OpPPMap(post_proc_fe->getPostProcMesh(),
 
  282                  post_proc_fe->getMapGaussPts(),
 
  297  double sum_temperature;
 
  300      new OpCalcSurfaceAverageTemperature(approxVals, sum_temperature,
 
  307  CHKERR post_proc_fe->writeFile(
"out_radiation.h5m");
 
  309  MOFEM_LOG_C(
"EXAMPLE", Sev::inform, 
"Surface area %3.4e [km]", surface_area);
 
  311              "Average subsurface temperatute %3.4e [K]",
 
  312              sum_temperature / surface_area);
 
  322int main(
int argc, 
char *argv[]) {
 
  325  const char param_file[] = 
"param_file.petsc";
 
  329  auto core_log = logging::core::get();
 
  338    DMType dm_name = 
"DMMOFEM";
 
  343    moab::Core mb_instance;              
 
  344    moab::Interface &moab = mb_instance; 
 
 
  373  const double vol = this->getMeasure();
 
  375  auto t_w = getFTensor0IntegrationWeight();
 
  381  auto t_coords = getFTensor1CoordsAtGaussPts();
 
  387    const double r_cylinder = t_coords(0);
 
  390    const double alpha = t_w * vol * 
Beta * (2 * M_PI * r_cylinder);
 
  392    for (
int rr = 0; rr != 
nbRows; ++rr) {
 
  395      for (
int cc = 0; cc != 
nbCols; cc++) {
 
  396        if (std::abs(t_coords(0)) > std::numeric_limits<double>::epsilon()) {
 
  397          locMat(rr, cc) += alpha * t_row_base * t_col_base * 4 *
 
  398                            std::pow(
static_cast<double>(t_val), 3);
 
 
  417  const double vol = getMeasure();
 
  419  auto t_w = getFTensor0IntegrationWeight();
 
  425  auto t_coords = getFTensor1CoordsAtGaussPts();
 
  428  for (
int gg = 0; gg != nbIntegrationPts; gg++) {
 
  431    const double r_cylinder = t_coords(0);
 
  434    const double alpha = t_w * vol * 
Beta * (2 * M_PI * r_cylinder);
 
  436    for (
int rr = 0; rr != nbRows; ++rr) {
 
  437      if (std::abs(t_coords(0)) > std::numeric_limits<double>::epsilon()) {
 
  438        locF[rr] += alpha * t_row_base *
 
  439                    (std::pow(
static_cast<double>(t_val), 4) -
 
 
  457  const double vol = getMeasure();
 
  459  auto t_w = getFTensor0IntegrationWeight();
 
  463  auto t_coords = getFTensor1CoordsAtGaussPts();
 
  465  const double time = getFEMethod()->ts_t;
 
  468  constexpr double flux_p = -0.03e6;
 
  469  constexpr double flux_c = -0.23e6;
 
  472  for (
int gg = 0; gg != nbIntegrationPts; gg++) {
 
  475    const double r_cylinder = t_coords(0);
 
  477    const double r = std::sqrt(t_coords(
i) * t_coords(
i));
 
  478    const double s = std::abs(t_coords(1)) / r;
 
  481    const double alpha = t_w * vol * (2 * M_PI * r_cylinder);
 
  483    for (
int rr = 0; rr != nbRows; ++rr) {
 
  484      locF[rr] += alpha * t_row_base * (s * flux_p + flux_c) * time;
 
 
  500  if (type == MBVERTEX) {
 
  502    const double vol = getMeasure();
 
  504    auto t_w = getFTensor0IntegrationWeight();
 
  508    auto t_coords = getFTensor1CoordsAtGaussPts();
 
  510    size_t nb_integration_pts = getGaussPts().size2();
 
  513    for (
auto gg = 0; gg != nb_integration_pts; ++gg) {
 
  516      const double r_cylinder = t_coords(0);
 
  519      const double alpha = t_w * vol * (2 * M_PI * r_cylinder);
 
  521      sumTemperature += alpha * t_val;
 
  522      surfaceArea += alpha;
 
 
#define MOFEM_LOG_C(channel, severity, format,...)
void simple(double P1[], double P2[], double P3[], double c[], const int N)
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
auto createDMVector(DM dm)
Get smart vector from DM.
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::ptr_deque< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
SmartPetscObj< TS > createTSIM(SmartPetscObj< DM > dm=nullptr)
Create TS (time) implicit solver.
boost::ptr_deque< UserDataOperator > & getOpBoundaryLhsPipeline()
Get the Op Boundary Lhs Pipeline object.
boost::ptr_deque< UserDataOperator > & getOpBoundaryRhsPipeline()
Get the Op Boundary Rhs Pipeline object.
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
FTensor::Index< 'i', SPACE_DIM > i
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, SPACE_DIM > OpDomainGradGrad
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, 1, SPACE_DIM > OpDomainGradTimesVec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorShallowArrayAdaptor< double > VectorAdaptor
implementation of Data Operators for Forces and Sources
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
constexpr double boltzmann_constant
constexpr double T_ambient
constexpr double emissivity
constexpr double heat_conductivity
OpCalcSurfaceAverageTemperature(boost::shared_ptr< VectorDouble > &approx_vals, double &sum_temp, double &surf)
boost::shared_ptr< VectorDouble > approxVals
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
[Flux Rhs]
FTensor::Index< 'i', 2 > i
summit Index
MoFEMErrorCode iNtegrate(EntData &row_data)
[Radiation Lhs]
boost::shared_ptr< VectorDouble > approxVals
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
[Radiation Lhs]
OpRadiationLhs(boost::shared_ptr< VectorDouble > &approx_vals)
OpRadiationRhs(boost::shared_ptr< VectorDouble > &approx_vals)
boost::shared_ptr< VectorDouble > approxVals
MoFEMErrorCode iNtegrate(EntData &row_data)
[Radiation Lhs]
static int integrationRule(int, int, int p_data)
boost::shared_ptr< VectorDouble > approxVals
MoFEMErrorCode kspSolve()
[Push operators to pipeline]
MoFEMErrorCode checkResults()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
MoFEMErrorCode runProblem()
boost::shared_ptr< MatrixDouble > approxGradVals
MoFEM::Interface & mField
Reference to MoFEM interface.
MoFEMErrorCode postProcess()
MoFEMErrorCode setupProblem()
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
default operator for EDGE element
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.
MoFEMErrorCode setVertexDofs(VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
Set DOFs on vertices using user function.
@ OPROW
operator doWork function is executed on FE rows
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
int nbRows
number of dofs on rows
int nbIntegrationPts
number of integration points
MatrixDouble locMat
local entity block matrix
int nbCols
number if dof on column
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Specialization for double precision scalar field values calculation.
Operator for inverting matrices at integration points.
Post post-proc data at points from hash maps.
Set inverse jacobian to base functions.
Modify integration weights on face to take into account higher-order geometry.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
Get domain right-hand side finite element.
boost::shared_ptr< FEMethod > & getDomainLhsFE()
Get domain left-hand side finite element.
MoFEM::FaceElementForcesAndSourcesCore FaceEle
boost::shared_ptr< FEMethod > & getBoundaryLhsFE()
Get boundary left-hand side finite element.
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
Set integration rule for domain right-hand side finite element.
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
Set integration rule for boundary left-hand side finite element.
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
Set integration rule for boundary right-hand side finite element.
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()
Get boundary right-hand side finite element.
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Set integration rule for domain left-hand side finite element.
Simple interface for fast problem set-up.
MoFEMErrorCode addDomainField(const std::string name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
MoFEMErrorCode getOptions()
get options
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.