v0.13.2
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Public Attributes | List of all members
ThermalElement::TimeSeriesMonitor Struct Reference

TS monitore it records temperature at time steps. More...

#include <users_modules/basic_finite_elements/src/ThermalElement.hpp>

Inheritance diagram for ThermalElement::TimeSeriesMonitor:
[legend]
Collaboration diagram for ThermalElement::TimeSeriesMonitor:
[legend]

Public Types

using Ele = ForcesAndSourcesCore
 
using VolEle = VolumeElementForcesAndSourcesCore
 
using VolOp = VolumeElementForcesAndSourcesCore::UserDataOperator
 
using SetPtsData = FieldEvaluatorInterface::SetPtsData
 

Public Member Functions

 TimeSeriesMonitor (MoFEM::Interface &m_field, const std::string series_name, const std::string temp_name, std::vector< std::array< double, 3 > > eval_points={})
 
MoFEMErrorCode postProcess ()
 

Public Attributes

MoFEM::InterfacemField
 
const std::string seriesName
 
const std::string tempName
 
BitRefLevel mask
 
boost::shared_ptr< SetPtsDatadataFieldEval
 
boost::shared_ptr< VectorDouble > tempPtr
 
std::vector< std::array< double, 3 > > evalPoints
 

Detailed Description

TS monitore it records temperature at time steps.

Definition at line 579 of file ThermalElement.hpp.

Member Typedef Documentation

◆ Ele

using ThermalElement::TimeSeriesMonitor::Ele = ForcesAndSourcesCore

Definition at line 586 of file ThermalElement.hpp.

◆ SetPtsData

using ThermalElement::TimeSeriesMonitor::SetPtsData = FieldEvaluatorInterface::SetPtsData

Definition at line 589 of file ThermalElement.hpp.

◆ VolEle

using ThermalElement::TimeSeriesMonitor::VolEle = VolumeElementForcesAndSourcesCore

Definition at line 587 of file ThermalElement.hpp.

◆ VolOp

using ThermalElement::TimeSeriesMonitor::VolOp = VolumeElementForcesAndSourcesCore::UserDataOperator

Definition at line 588 of file ThermalElement.hpp.

Constructor & Destructor Documentation

◆ TimeSeriesMonitor()

ThermalElement::TimeSeriesMonitor::TimeSeriesMonitor ( MoFEM::Interface m_field,
const std::string  series_name,
const std::string  temp_name,
std::vector< std::array< double, 3 > >  eval_points = {} 
)
inline

Definition at line 596 of file ThermalElement.hpp.

598 {})
599 : mField(m_field), seriesName(series_name), tempName(temp_name),
600 evalPoints(eval_points),
601 dataFieldEval(m_field.getInterface<FieldEvaluatorInterface>()
602 ->getData<VolEle>()) {
603 mask.set();
604
605 if (!evalPoints.empty()) {
606 ierr = m_field.getInterface<FieldEvaluatorInterface>()->buildTree3D(
607 dataFieldEval, "THERMAL_FE");
608 CHKERRABORT(PETSC_COMM_WORLD, ierr);
609
610 auto no_rule = [](int, int, int) { return -1; };
611
612 tempPtr = boost::make_shared<VectorDouble>();
613
614 boost::shared_ptr<Ele> vol_ele(dataFieldEval->feMethodPtr.lock());
615 vol_ele->getRuleHook = no_rule;
616
617 vol_ele->getOpPtrVector().push_back(
618 new OpCalculateScalarFieldValues("TEMP", tempPtr));
619 }
620 }
static PetscErrorCode ierr
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
boost::shared_ptr< SetPtsData > dataFieldEval
std::vector< std::array< double, 3 > > evalPoints
boost::shared_ptr< VectorDouble > tempPtr

Member Function Documentation

◆ postProcess()

MoFEMErrorCode ThermalElement::TimeSeriesMonitor::postProcess ( )

Definition at line 454 of file ThermalElement.cpp.

454 {
456
457 CHKERR mField.getInterface<VecManager>()->setGlobalGhostVector(
458 problemPtr, ROW, ts_u, INSERT_VALUES, SCATTER_REVERSE);
459
460 BitRefLevel proble_bit_level = problemPtr->getBitRefLevel();
461
462 SeriesRecorder *recorder_ptr = NULL;
463 CHKERR mField.getInterface(recorder_ptr);
464 CHKERR recorder_ptr->record_begin(seriesName);
465 CHKERR recorder_ptr->record_field(seriesName, tempName, proble_bit_level,
466 mask);
467 CHKERR recorder_ptr->record_end(seriesName, ts_t);
468
469 auto post_proc_at_points = [&](std::array<double, 3> point, int num) {
471
472 dataFieldEval->setEvalPoints(point.data(), point.size() / 3);
473
474 struct OpPrint : public VolOp {
475
476 std::array<double, 3> pointCoords;
477 int pointNum;
478 boost::shared_ptr<VectorDouble> tempPtr;
479
480 OpPrint(boost::shared_ptr<VectorDouble> temp_ptr,
481 std::array<double, 3> &point_coords, int point_num)
482 : VolOp("TEMP", VolOp::OPROW), tempPtr(temp_ptr),
483 pointCoords(point_coords), pointNum(point_num) {}
484
485 MoFEMErrorCode doWork(int side, EntityType type,
486 DataForcesAndSourcesCore::EntData &data) {
488 if (type == MBVERTEX) {
489 if (getGaussPts().size2()) {
490
491 auto t_p = getFTensor0FromVec(*tempPtr);
492
493 MOFEM_LOG("THERMALSYNC", Sev::inform)
494 << "Pnt: " << std::to_string(pointNum)
495 << " Crd: " << getVectorAdaptor(pointCoords.data(), 3)
496 << " Tmp: " << t_p;
497 }
498 }
500 }
501 };
502
503 if (auto fe_ptr = dataFieldEval->feMethodPtr.lock()) {
504 fe_ptr->getOpPtrVector().push_back(new OpPrint(tempPtr, point, num));
505 CHKERR mField.getInterface<FieldEvaluatorInterface>()->evalFEAtThePoint3D(
506 point.data(), 1e-12, "DMTHERMAL", "THERMAL_FE", dataFieldEval,
508 QUIET);
509 fe_ptr->getOpPtrVector().pop_back();
510 }
511
513 };
514
515 if (!evalPoints.empty()) {
516 int num = 0;
517 for (auto p : evalPoints)
518 CHKERR post_proc_at_points(p, num++);
519 }
520
522
524}
static Index< 'p', 3 > p
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:345
@ QUIET
Definition: definitions.h:208
@ ROW
Definition: definitions.h:123
@ MF_EXIST
Definition: definitions.h:100
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#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
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
virtual MoFEMErrorCode record_field(const std::string &serie_name, const std::string &field_name, const BitRefLevel &bit, const BitRefLevel &mask)
virtual MoFEMErrorCode record_begin(const std::string &serie_name)
virtual MoFEMErrorCode record_end(const std::string &serie_name, double time=0)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition: Templates.hpp:31
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
Field evaluator interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:23

Member Data Documentation

◆ dataFieldEval

boost::shared_ptr<SetPtsData> ThermalElement::TimeSeriesMonitor::dataFieldEval

Definition at line 591 of file ThermalElement.hpp.

◆ evalPoints

std::vector<std::array<double, 3> > ThermalElement::TimeSeriesMonitor::evalPoints

Definition at line 594 of file ThermalElement.hpp.

◆ mask

BitRefLevel ThermalElement::TimeSeriesMonitor::mask

Definition at line 584 of file ThermalElement.hpp.

◆ mField

MoFEM::Interface& ThermalElement::TimeSeriesMonitor::mField

Definition at line 581 of file ThermalElement.hpp.

◆ seriesName

const std::string ThermalElement::TimeSeriesMonitor::seriesName

Definition at line 582 of file ThermalElement.hpp.

◆ tempName

const std::string ThermalElement::TimeSeriesMonitor::tempName

Definition at line 583 of file ThermalElement.hpp.

◆ tempPtr

boost::shared_ptr<VectorDouble> ThermalElement::TimeSeriesMonitor::tempPtr

Definition at line 592 of file ThermalElement.hpp.


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