v0.14.0
Files | Functions
Thermal element
Collaboration diagram for Thermal element:

Files

file  ThermalElement.cpp
 
file  ThermalElement.hpp
 Operators and data structures for thermal analysis.
 
file  ThermalStressElement.hpp
 Implementation of thermal stresses element.
 
file  thermal_steady.cpp
 Example of steady thermal analysis.
 
file  thermal_unsteady.cpp
 Example of thermal unsteady analyze.
 

Functions

MoFEMErrorCode ThermalElement::addThermalElements (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add thermal element on tets More...
 
MoFEMErrorCode ThermalElement::addThermalFluxElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add heat flux element More...
 
MoFEMErrorCode ThermalElement::addThermalConvectionElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add convection element More...
 
MoFEMErrorCode ThermalElement::addThermalRadiationElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add Non-linear Radiation element More...
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementRhsOperators (string field_name, Vec &F)
 this function is used in case of stationary problem to set elements for rhs More...
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementLhsOperators (string field_name, Mat A)
 this function is used in case of stationary heat conductivity problem for lhs More...
 
MoFEMErrorCode ThermalElement::setThermalFluxFiniteElementRhsOperators (string field_name, Vec &F, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 this function is used in case of stationary problem for heat flux terms More...
 
MoFEMErrorCode ThermalElement::setTimeSteppingProblem (string field_name, string rate_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 set up operators for unsteady heat flux; convection; radiation problem More...
 
MoFEMErrorCode ThermalElement::setTimeSteppingProblem (TsCtx &ts_ctx, string field_name, string rate_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 set up operators for unsteady heat flux; convection; radiation problem More...
 

Detailed Description

Function Documentation

◆ addThermalConvectionElement()

MoFEMErrorCode ThermalElement::addThermalConvectionElement ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

add convection element

It get data from convection set and define element in moab. Alternatively uses block set with name including substring CONVECTION.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 621 of file ThermalElement.cpp.

622  {
624 
625  CHKERR mField.add_finite_element("THERMAL_CONVECTION_FE", MF_ZERO);
626  CHKERR mField.modify_finite_element_add_field_row("THERMAL_CONVECTION_FE",
627  field_name);
628  CHKERR mField.modify_finite_element_add_field_col("THERMAL_CONVECTION_FE",
629  field_name);
630  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
631  field_name);
632  if (mField.check_field(mesh_nodals_positions)) {
633  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
634  mesh_nodals_positions);
635  }
636 
637  // this is alternative method for setting boundary conditions, to bypass bu
638  // in cubit file reader. not elegant, but good enough
640  if (std::regex_match(it->getName(), std::regex("(.*)CONVECTION(.*)"))) {
641  std::vector<double> data;
642  CHKERR it->getAttributes(data);
643  if (data.size() != 2) {
644  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
645  }
646  setOfConvection[it->getMeshsetId()].cOnvection = data[0];
647  setOfConvection[it->getMeshsetId()].tEmperature = data[1];
648  CHKERR mField.get_moab().get_entities_by_type(
649  it->meshset, MBTRI, setOfConvection[it->getMeshsetId()].tRis, true);
651  setOfConvection[it->getMeshsetId()].tRis, MBTRI,
652  "THERMAL_CONVECTION_FE");
653  }
654  }
655 
657 }

◆ addThermalElements()

MoFEMErrorCode ThermalElement::addThermalElements ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

add thermal element on tets

It get data from block set and define element in moab w

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 527 of file ThermalElement.cpp.

528  {
530 
531  CHKERR mField.add_finite_element("THERMAL_FE", MF_ZERO);
535  if (mField.check_field(mesh_nodals_positions)) {
537  mesh_nodals_positions);
538  }
539 
540  // takes skin of block of entities
541  // Skinner skin(&mField.get_moab());
542  // loop over all blocksets and get data which name is FluidPressure
544  mField, BLOCKSET | MAT_THERMALSET, it)) {
545 
546  Mat_Thermal temp_data;
547  ierr = it->getAttributeDataStructure(temp_data);
548 
549  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.resize(
550  3, 3); //(3X3) conductivity matrix
551  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.clear();
552  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(0, 0) =
553  temp_data.data.Conductivity;
554  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(1, 1) =
555  temp_data.data.Conductivity;
556  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(2, 2) =
557  temp_data.data.Conductivity;
558  // setOfBlocks[it->getMeshsetId()].cOnductivity =
559  // temp_data.data.Conductivity;
560  setOfBlocks[it->getMeshsetId()].cApacity = temp_data.data.HeatCapacity;
561  if (temp_data.data.User2 != 0) {
562  setOfBlocks[it->getMeshsetId()].initTemp = temp_data.data.User2;
563  }
564  CHKERR mField.get_moab().get_entities_by_type(
565  it->meshset, MBTET, setOfBlocks[it->getMeshsetId()].tEts, true);
567  setOfBlocks[it->getMeshsetId()].tEts, MBTET, "THERMAL_FE");
568  }
569 
571 }

◆ addThermalFluxElement()

MoFEMErrorCode ThermalElement::addThermalFluxElement ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

add heat flux element

It get data from heat flux set and define element in moab. Alternatively uses block set with name including substring HEAT_FLUX.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 574 of file ThermalElement.cpp.

575  {
577 
578  CHKERR mField.add_finite_element("THERMAL_FLUX_FE", MF_ZERO);
580  field_name);
582  field_name);
584  field_name);
585  if (mField.check_field(mesh_nodals_positions)) {
587  mesh_nodals_positions);
588  }
589 
591  it)) {
592  CHKERR it->getBcDataStructure(setOfFluxes[it->getMeshsetId()].dAta);
593  CHKERR mField.get_moab().get_entities_by_type(
594  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
596  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
597  }
598 
599  // this is alternative method for setting boundary conditions, to bypass bu
600  // in cubit file reader. not elegant, but good enough
602  if (std::regex_match(it->getName(), std::regex("(.*)HEAT_FLUX(.*)"))) {
603  std::vector<double> data;
604  CHKERR it->getAttributes(data);
605  if (data.size() != 1) {
606  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
607  }
608  strcpy(setOfFluxes[it->getMeshsetId()].dAta.data.name, "HeatFlu");
609  setOfFluxes[it->getMeshsetId()].dAta.data.flag1 = 1;
610  setOfFluxes[it->getMeshsetId()].dAta.data.value1 = data[0];
611  CHKERR mField.get_moab().get_entities_by_type(
612  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
614  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
615  }
616  }
617 
619 }

◆ addThermalRadiationElement()

MoFEMErrorCode ThermalElement::addThermalRadiationElement ( const std::string  field_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

add Non-linear Radiation element

It get data from Radiation set and define element in moab. Alternatively uses block set with name including substring RADIATION.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 659 of file ThermalElement.cpp.

660  {
662 
663  CHKERR mField.add_finite_element("THERMAL_RADIATION_FE", MF_ZERO);
664  CHKERR mField.modify_finite_element_add_field_row("THERMAL_RADIATION_FE",
665  field_name);
666  CHKERR mField.modify_finite_element_add_field_col("THERMAL_RADIATION_FE",
667  field_name);
668  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
669  field_name);
670  if (mField.check_field(mesh_nodals_positions)) {
671  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
672  mesh_nodals_positions);
673  }
674 
675  // this is alternative method for setting boundary conditions, to bypass bu
676  // in cubit file reader. not elegant, but good enough
678  if (std::regex_match(it->getName(), std::regex("(.*)RADIATION(.*)"))) {
679  std::vector<double> data;
680  ierr = it->getAttributes(data);
681  if (data.size() != 3) {
682  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
683  }
684  setOfRadiation[it->getMeshsetId()].sIgma = data[0];
685  setOfRadiation[it->getMeshsetId()].eMissivity = data[1];
686  setOfRadiation[it->getMeshsetId()].aMbienttEmp = data[2];
687  CHKERR mField.get_moab().get_entities_by_type(
688  it->meshset, MBTRI, setOfRadiation[it->getMeshsetId()].tRis, true);
690  setOfRadiation[it->getMeshsetId()].tRis, MBTRI,
691  "THERMAL_RADIATION_FE");
692  }
693  }
694 
696 }

◆ setThermalFiniteElementLhsOperators()

MoFEMErrorCode ThermalElement::setThermalFiniteElementLhsOperators ( string  field_name,
Mat  A 
)

this function is used in case of stationary heat conductivity problem for lhs

Definition at line 713 of file ThermalElement.cpp.

713  {
715  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
716  for (; sit != setOfBlocks.end(); sit++) {
717  // add finite elemen
718  feLhs.getOpPtrVector().push_back(
719  new OpThermalLhs(field_name, A, sit->second, commonData));
720  }
722 }

◆ setThermalFiniteElementRhsOperators()

MoFEMErrorCode ThermalElement::setThermalFiniteElementRhsOperators ( string  field_name,
Vec &  F 
)

this function is used in case of stationary problem to set elements for rhs

Definition at line 699 of file ThermalElement.cpp.

699  {
701  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
702  feRhs.getOpPtrVector().push_back(
703  new OpGetGradAtGaussPts(field_name, commonData));
704  for (; sit != setOfBlocks.end(); sit++) {
705  // add finite element
706  feRhs.getOpPtrVector().push_back(
707  new OpThermalRhs(field_name, F, sit->second, commonData));
708  }
710 }

◆ setThermalFluxFiniteElementRhsOperators()

MoFEMErrorCode ThermalElement::setThermalFluxFiniteElementRhsOperators ( string  field_name,
Vec &  F,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

this function is used in case of stationary problem for heat flux terms

Definition at line 724 of file ThermalElement.cpp.

725  {
727  bool hoGeometry = false;
728  if (mField.check_field(mesh_nodals_positions)) {
729  hoGeometry = true;
730  }
731  std::map<int, FluxData>::iterator sit = setOfFluxes.begin();
732  for (; sit != setOfFluxes.end(); sit++) {
733  // add finite element
734  feFlux.getOpPtrVector().push_back(
735  new OpHeatFlux(field_name, F, sit->second, hoGeometry));
736  }
738 }

◆ setTimeSteppingProblem() [1/2]

MoFEMErrorCode ThermalElement::setTimeSteppingProblem ( string  field_name,
string  rate_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

set up operators for unsteady heat flux; convection; radiation problem

Definition at line 774 of file ThermalElement.cpp.

776  {
778 
779  bool hoGeometry = false;
780  if (mField.check_field(mesh_nodals_positions)) {
781  hoGeometry = true;
782  }
783 
784  {
785  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
786  for (; sit != setOfBlocks.end(); sit++) {
787  // add finite element
788  // those methods are to calculate matrices on Lhs
789  // feLhs.getOpPtrVector().push_back(new
790  // OpGetTetTemperatureAtGaussPts(field_name,commonData));
791  feLhs.getOpPtrVector().push_back(
792  new OpThermalLhs(field_name, sit->second, commonData));
793  feLhs.getOpPtrVector().push_back(
794  new OpHeatCapacityLhs(field_name, sit->second, commonData));
795  // those methods are to calculate vectors on Rhs
796  feRhs.getOpPtrVector().push_back(
797  new OpGetTetTemperatureAtGaussPts(field_name, commonData));
798  feRhs.getOpPtrVector().push_back(
799  new OpGetTetRateAtGaussPts(rate_name, commonData));
800  feRhs.getOpPtrVector().push_back(
801  new OpGetGradAtGaussPts(field_name, commonData));
802  feRhs.getOpPtrVector().push_back(
803  new OpThermalRhs(field_name, sit->second, commonData));
804  feRhs.getOpPtrVector().push_back(
805  new OpHeatCapacityRhs(field_name, sit->second, commonData));
806  }
807  }
808 
809  // Flux
810  {
811  std::map<int, FluxData>::iterator sit = setOfFluxes.begin();
812  for (; sit != setOfFluxes.end(); sit++) {
813  feFlux.getOpPtrVector().push_back(
814  new OpHeatFlux(field_name, sit->second, hoGeometry));
815  }
816  }
817 
818  // Convection
819  {
820  std::map<int, ConvectionData>::iterator sit = setOfConvection.begin();
821  for (; sit != setOfConvection.end(); sit++) {
822  feConvectionRhs.getOpPtrVector().push_back(
823  new OpGetTriTemperatureAtGaussPts(field_name, commonData));
824  feConvectionRhs.getOpPtrVector().push_back(
825  new OpConvectionRhs(field_name, sit->second, commonData, hoGeometry));
826  }
827  }
828  {
829  std::map<int, ConvectionData>::iterator sit = setOfConvection.begin();
830  for (; sit != setOfConvection.end(); sit++) {
831  feConvectionLhs.getOpPtrVector().push_back(
832  new OpConvectionLhs(field_name, sit->second, hoGeometry));
833  }
834  }
835 
836  // Radiation
837  {
838  std::map<int, RadiationData>::iterator sit = setOfRadiation.begin();
839  for (; sit != setOfRadiation.end(); sit++) {
840  feRadiationRhs.getOpPtrVector().push_back(
841  new OpGetTriTemperatureAtGaussPts(field_name, commonData));
842  feRadiationRhs.getOpPtrVector().push_back(
843  new OpRadiationRhs(field_name, sit->second, commonData, hoGeometry));
844  }
845  }
846  {
847  std::map<int, RadiationData>::iterator sit = setOfRadiation.begin();
848  for (; sit != setOfRadiation.end(); sit++) {
849  feRadiationLhs.getOpPtrVector().push_back(
850  new OpGetTriTemperatureAtGaussPts(field_name, commonData));
851  feRadiationLhs.getOpPtrVector().push_back(
852  new OpRadiationLhs(field_name, sit->second, commonData, hoGeometry));
853  }
854  }
855 
857 }

◆ setTimeSteppingProblem() [2/2]

MoFEMErrorCode ThermalElement::setTimeSteppingProblem ( TsCtx &  ts_ctx,
string  field_name,
string  rate_name,
const std::string  mesh_nodals_positions = "MESH_NODE_POSITIONS" 
)

set up operators for unsteady heat flux; convection; radiation problem

MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
SIDESET
@ SIDESET
Definition: definitions.h:160
ThermalElement::setOfBlocks
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
Definition: ThermalElement.hpp:100
ThermalElement::setOfConvection
std::map< int, ConvectionData > setOfConvection
Definition: ThermalElement.hpp:124
MoFEM::CoreInterface::modify_finite_element_add_field_row
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
_IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
Definition: MeshsetsManager.hpp:49
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
ThermalElement::feRhs
MyVolumeFE feRhs
cauclate right hand side for tetrahedral elements
Definition: ThermalElement.hpp:51
MoFEM::Mat_Thermal::data
_data_ data
Definition: MaterialBlocks.hpp:217
ThermalElement::feRadiationLhs
MyTriFE feRadiationLhs
Definition: ThermalElement.hpp:77
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::add_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
ThermalElement::commonData
CommonData commonData
Definition: ThermalElement.hpp:152
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM::CoreInterface::modify_finite_element_add_field_col
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MAT_THERMALSET
@ MAT_THERMALSET
block name is "MAT_THERMAL"
Definition: definitions.h:174
ThermalElement::mField
MoFEM::Interface & mField
Definition: ThermalElement.hpp:83
MoFEM::CoreInterface::check_field
virtual bool check_field(const std::string &name) const =0
check if field is in database
ThermalElement::feFlux
MyTriFE feFlux
Definition: ThermalElement.hpp:66
ThermalElement::feConvectionLhs
MyTriFE feConvectionLhs
Definition: ThermalElement.hpp:70
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
ThermalElement::setOfFluxes
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
Definition: ThermalElement.hpp:111
ThermalElement::feLhs
MyVolumeFE feLhs
Definition: ThermalElement.hpp:53
ThermalElement::setOfRadiation
std::map< int, RadiationData > setOfRadiation
Definition: ThermalElement.hpp:139
MF_ZERO
@ MF_ZERO
Definition: definitions.h:111
_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_
#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.
Definition: MeshsetsManager.hpp:71
ThermalElement::feRadiationRhs
MyTriFE feRadiationRhs
Definition: ThermalElement.hpp:76
MoFEM::Mat_Thermal
Thermal material data structure.
Definition: MaterialBlocks.hpp:201
BLOCKSET
@ BLOCKSET
Definition: definitions.h:161
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
MoFEM::ForcesAndSourcesCore::getOpPtrVector
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
Definition: ForcesAndSourcesCore.hpp:83
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
HEATFLUXSET
@ HEATFLUXSET
Definition: definitions.h:169
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
ThermalElement::feConvectionRhs
MyTriFE feConvectionRhs
Definition: ThermalElement.hpp:69
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
F
@ F
Definition: free_surface.cpp:394