v0.15.0
Loading...
Searching...
No Matches
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.
 

Classes

struct  ThermalElement
 structure grouping operators and data used for thermal problems More...
 
struct  ThermalElement::MyVolumeFE
 definition of volume element More...
 
struct  ThermalElement::MyTriFE
 define surface element More...
 
struct  ThermalElement::BlockData
 data for calculation heat conductivity and heat capacity elements More...
 
struct  ThermalElement::FluxData
 data for calculation heat flux More...
 
struct  ThermalElement::ConvectionData
 data for convection More...
 
struct  ThermalElement::RadiationData
 data for radiation More...
 
struct  ThermalElement::CommonData
 common data used by volume elements More...
 
struct  ThermalElement::OpGetGradAtGaussPts
 operator to calculate temperature gradient at Gauss points More...
 
struct  ThermalElement::OpGetFieldAtGaussPts< OP >
 operator to calculate temperature and rate of temperature at Gauss points More...
 
struct  ThermalElement::OpGetTetTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetTriTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetTetRateAtGaussPts
 operator to calculate temperature rate at Gauss pts More...
 
struct  ThermalElement::OpThermalRhs
 
struct  ThermalElement::OpThermalLhs
 
struct  ThermalElement::OpHeatCapacityRhs
 operator to calculate right hand side of heat capacity terms More...
 
struct  ThermalElement::OpHeatCapacityLhs
 operator to calculate left hand side of heat capacity terms More...
 
struct  ThermalElement::OpHeatFlux
 operator for calculate heat flux and assemble to right hand side More...
 
struct  ThermalElement::OpRadiationLhs
 
struct  ThermalElement::OpRadiationRhs
 operator to calculate radiation therms on body surface and assemble to rhs of transient equations(Residual Vector) More...
 
struct  ThermalElement::OpConvectionRhs
 operator to calculate convection therms on body surface and assemble to rhs of equations More...
 
struct  ThermalElement::OpConvectionLhs
 
struct  ThermalElement::UpdateAndControl
 this calass is to control time stepping More...
 
struct  ThermalElement::TimeSeriesMonitor
 TS monitore it records temperature at time steps. More...
 
struct  ThermalStressElement
 Implentation of thermal stress element. More...
 
struct  ThermalStressElement::MyVolumeFE
 
struct  ThermalStressElement::BlockData
 
struct  ThermalStressElement::CommonData
 
struct  ThermalStressElement::OpGetTemperatureAtGaussPts
 
struct  ThermalStressElement::OpThermalStressRhs
 

Functions

MoFEMErrorCode ThermalElement::addThermalElements (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add thermal element on tets
 
MoFEMErrorCode ThermalElement::addThermalFluxElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add heat flux element
 
MoFEMErrorCode ThermalElement::addThermalConvectionElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add convection element
 
MoFEMErrorCode ThermalElement::addThermalRadiationElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add Non-linear Radiation element
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementRhsOperators (string field_name, Vec &F)
 this function is used in case of stationary problem to set elements for rhs
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementLhsOperators (string field_name, Mat A)
 this function is used in case of stationary heat conductivity problem for lhs
 
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
 
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
 
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
 

Detailed Description

Function Documentation

◆ addThermalConvectionElement()

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

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

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 622 of file ThermalElement.cpp.

623 {
625
626 CHKERR mField.add_finite_element("THERMAL_CONVECTION_FE", MF_ZERO);
627 CHKERR mField.modify_finite_element_add_field_row("THERMAL_CONVECTION_FE",
628 field_name);
629 CHKERR mField.modify_finite_element_add_field_col("THERMAL_CONVECTION_FE",
630 field_name);
631 CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
632 field_name);
633 if (mField.check_field(mesh_nodals_positions)) {
634 CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
635 mesh_nodals_positions);
636 }
637
638 // this is alternative method for setting boundary conditions, to bypass bu
639 // in cubit file reader. not elegant, but good enough
641 if (std::regex_match(it->getName(), std::regex("(.*)CONVECTION(.*)"))) {
642 std::vector<double> data;
643 CHKERR it->getAttributes(data);
644 if (data.size() != 2) {
645 SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
646 }
647 setOfConvection[it->getMeshsetId()].cOnvection = data[0];
648 setOfConvection[it->getMeshsetId()].tEmperature = data[1];
649 CHKERR mField.get_moab().get_entities_by_type(
650 it->meshset, MBTRI, setOfConvection[it->getMeshsetId()].tRis, true);
652 setOfConvection[it->getMeshsetId()].tRis, MBTRI,
653 "THERMAL_CONVECTION_FE");
654 }
655 }
656
658}
@ MF_ZERO
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ BLOCKSET
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
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 modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
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
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
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
#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.
constexpr auto field_name
virtual moab::Interface & get_moab()=0
MoFEM::Interface & mField
std::map< int, ConvectionData > setOfConvection

◆ addThermalElements()

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

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

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 528 of file ThermalElement.cpp.

529 {
531
536 if (mField.check_field(mesh_nodals_positions)) {
538 mesh_nodals_positions);
539 }
540
541 // takes skin of block of entities
542 // Skinner skin(&mField.get_moab());
543 // loop over all blocksets and get data which name is FluidPressure
546
547 Mat_Thermal temp_data;
548 ierr = it->getAttributeDataStructure(temp_data);
549
550 setOfBlocks[it->getMeshsetId()].cOnductivity_mat.resize(
551 3, 3); //(3X3) conductivity matrix
552 setOfBlocks[it->getMeshsetId()].cOnductivity_mat.clear();
553 setOfBlocks[it->getMeshsetId()].cOnductivity_mat(0, 0) =
554 temp_data.data.Conductivity;
555 setOfBlocks[it->getMeshsetId()].cOnductivity_mat(1, 1) =
556 temp_data.data.Conductivity;
557 setOfBlocks[it->getMeshsetId()].cOnductivity_mat(2, 2) =
558 temp_data.data.Conductivity;
559 // setOfBlocks[it->getMeshsetId()].cOnductivity =
560 // temp_data.data.Conductivity;
561 setOfBlocks[it->getMeshsetId()].cApacity = temp_data.data.HeatCapacity;
562 if (temp_data.data.User2 != 0) {
563 setOfBlocks[it->getMeshsetId()].initTemp = temp_data.data.User2;
564 }
565 CHKERR mField.get_moab().get_entities_by_type(
566 it->meshset, MBTET, setOfBlocks[it->getMeshsetId()].tEts, true);
568 setOfBlocks[it->getMeshsetId()].tEts, MBTET, "THERMAL_FE");
569 }
570
572}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MAT_THERMALSET
block name is "MAT_THERMAL"
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Thermal material data structure.
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData

◆ addThermalFluxElement()

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

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

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 575 of file ThermalElement.cpp.

576 {
578
579 CHKERR mField.add_finite_element("THERMAL_FLUX_FE", MF_ZERO);
581 field_name);
583 field_name);
585 field_name);
586 if (mField.check_field(mesh_nodals_positions)) {
588 mesh_nodals_positions);
589 }
590
592 it)) {
593 CHKERR it->getBcDataStructure(setOfFluxes[it->getMeshsetId()].dAta);
594 CHKERR mField.get_moab().get_entities_by_type(
595 it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
597 setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
598 }
599
600 // this is alternative method for setting boundary conditions, to bypass bu
601 // in cubit file reader. not elegant, but good enough
603 if (std::regex_match(it->getName(), std::regex("(.*)HEAT_FLUX(.*)"))) {
604 std::vector<double> data;
605 CHKERR it->getAttributes(data);
606 if (data.size() != 1) {
607 SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
608 }
609 strcpy(setOfFluxes[it->getMeshsetId()].dAta.data.name, "HeatFlu");
610 setOfFluxes[it->getMeshsetId()].dAta.data.flag1 = 1;
611 setOfFluxes[it->getMeshsetId()].dAta.data.value1 = data[0];
612 CHKERR mField.get_moab().get_entities_by_type(
613 it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
615 setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
616 }
617 }
618
620}
@ HEATFLUXSET
@ SIDESET
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData

◆ addThermalRadiationElement()

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

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

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 660 of file ThermalElement.cpp.

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

◆ setThermalFiniteElementLhsOperators()

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

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

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

Definition at line 714 of file ThermalElement.cpp.

714 {
716 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
717 for (; sit != setOfBlocks.end(); sit++) {
718 // add finite elemen
719 feLhs.getOpPtrVector().push_back(
720 new OpThermalLhs(field_name, A, sit->second, commonData));
721 }
723}
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
CommonData commonData

◆ setThermalFiniteElementRhsOperators()

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

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

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

Definition at line 700 of file ThermalElement.cpp.

700 {
702 std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
703 feRhs.getOpPtrVector().push_back(
704 new OpGetGradAtGaussPts(field_name, commonData));
705 for (; sit != setOfBlocks.end(); sit++) {
706 // add finite element
707 feRhs.getOpPtrVector().push_back(
708 new OpThermalRhs(field_name, F, sit->second, commonData));
709 }
711}
@ F
MyVolumeFE feRhs
cauclate right hand side for tetrahedral elements

◆ setThermalFluxFiniteElementRhsOperators()

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

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

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

Definition at line 725 of file ThermalElement.cpp.

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

◆ setTimeSteppingProblem() [1/2]

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

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

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

Definition at line 775 of file ThermalElement.cpp.

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

◆ setTimeSteppingProblem() [2/2]

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

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

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

Definition at line 860 of file ThermalElement.cpp.

862 {
864
865 CHKERR setTimeSteppingProblem(field_name, rate_name, mesh_nodals_positions);
866
867 // rhs
868 TsCtx::FEMethodsSequence &loops_to_do_Rhs =
870 loops_to_do_Rhs.push_back(TsCtx::PairNameFEMethodPtr("THERMAL_FE", &feRhs));
871 loops_to_do_Rhs.push_back(
872 TsCtx::PairNameFEMethodPtr("THERMAL_FLUX_FE", &feFlux));
873 if (mField.check_finite_element("THERMAL_CONVECTION_FE"))
874 loops_to_do_Rhs.push_back(
875 TsCtx::PairNameFEMethodPtr("THERMAL_CONVECTION_FE", &feConvectionRhs));
876 if (mField.check_finite_element("THERMAL_RADIATION_FE"))
877 loops_to_do_Rhs.push_back(
878 TsCtx::PairNameFEMethodPtr("THERMAL_RADIATION_FE", &feRadiationRhs));
879
880 // lhs
881 TsCtx::FEMethodsSequence &loops_to_do_Mat =
883 loops_to_do_Mat.push_back(TsCtx::PairNameFEMethodPtr("THERMAL_FE", &feLhs));
884 if (mField.check_finite_element("THERMAL_CONVECTION_FE"))
885 loops_to_do_Mat.push_back(
886 TsCtx::PairNameFEMethodPtr("THERMAL_CONVECTION_FE", &feConvectionLhs));
887 if (mField.check_finite_element("THERMAL_RADIATION_FE"))
888 loops_to_do_Mat.push_back(
889 TsCtx::PairNameFEMethodPtr("THERMAL_RADIATION_FE", &feRadiationLhs));
890
892}
MoFEMErrorCode 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
MoFEM::TsCtx * ts_ctx
virtual bool check_finite_element(const std::string &name) const =0
Check if finite element is in database.
MoFEM::FEMethodsSequence FEMethodsSequence
Definition TsCtx.hpp:26
FEMethodsSequence & getLoopsIFunction()
Get the loops to do IFunction object.
Definition TsCtx.hpp:63
FEMethodsSequence & getLoopsIJacobian()
Get the loops to do IJacobian object.
Definition TsCtx.hpp:83