v0.9.0
Files | Classes | 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.
 

Classes

struct  MixTransport::MixTransportElement::BlockData
 data for evaluation of het conductivity and heat capacity elements More...
 
struct  ThermalElement
 structure grouping operators and data used for thermal problemsIn order to assemble matrices and right hand vectors, the loops over elements, entities within the element and finally loop over integration points are executed. More...
 
struct  ThermalElement::TimeSeriesMonitor
 TS monitore it records temperature at time steps. More...
 
struct  ThermalElement::UpdateAndControl
 this calass is to control time steppingIt is used to save data for temperature rate vector to MoFEM field. More...
 
struct  ThermalElement::OpConvectionLhs
 
struct  ThermalElement::OpConvectionRhs
 operator to calculate convection therms on body surface and assemble to rhs of equations More...
 
struct  ThermalElement::OpRadiationRhs
 operator to calculate radiation therms on body surface and assemble to rhs of transient equations(Residual Vector) More...
 
struct  ThermalElement::OpRadiationLhs
 
struct  ThermalElement::OpHeatFlux
 operator for calculate heat flux and assemble to right hand side More...
 
struct  ThermalElement::OpHeatCapacityLhs
 operator to calculate left hand side of heat capacity terms More...
 
struct  ThermalElement::OpHeatCapacityRhs
 operator to calculate right hand side of heat capacity terms More...
 
struct  ThermalElement::OpThermalLhs
 
struct  ThermalElement::OpThermalRhs
 
struct  ThermalElement::OpGetTetRateAtGaussPts
 operator to calculate temperature rate at Gauss pts More...
 
struct  ThermalElement::OpGetTriTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetTetTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetFieldAtGaussPts< OP >
 operator to calculate temperature and rate of temperature at Gauss points More...
 
struct  ThermalElement::OpGetGradAtGaussPts
 operator to calculate temperature gradient at Gauss points More...
 
struct  ThermalElement::CommonData
 common data used by volume elements More...
 
struct  ThermalElement::RadiationData
 data for radiation More...
 
struct  ThermalElement::ConvectionData
 data for convection More...
 
struct  ThermalElement::FluxData
 data for calculation heat flux More...
 
struct  ThermalElement::BlockData
 data for calculation heat conductivity and heat capacity elements More...
 
struct  ThermalElement::MyTriFE
 define surface element More...
 
struct  ThermalElement::MyVolumeFE
 definition of volume element More...
 
struct  ThermalStressElement
 Implentation of thermal stress element. More...
 
struct  ThermalStressElement::OpThermalStressRhs
 
struct  ThermalStressElement::OpGetTemperatureAtGaussPts
 
struct  ThermalStressElement::CommonData
 
struct  ThermalStressElement::BlockData
 
struct  ThermalStressElement::MyVolumeFE
 

Functions

MoFEMErrorCode ThermalElement::addThermalElements (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add thermal element on tetsIt get data from block set and define element in moab w More...
 
MoFEMErrorCode ThermalElement::addThermalFluxElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add heat flux elementIt get data from heat flux set and define element in moab. Alternatively uses block set with name HEAT_FLUX. More...
 
MoFEMErrorCode ThermalElement::addThermalConvectionElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add convection elementIt get data from convection set and define element in moab. Alternatively uses block set with name CONVECTION. More...
 
MoFEMErrorCode ThermalElement::addThermalRadiationElement (const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
 add Non-linear Radiation elementIt get data from Radiation set and define element in moab. Alternatively uses block set with name RADIATION. 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 elementIt get data from convection set and define element in moab. Alternatively uses block set with name CONVECTION.

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

Definition at line 571 of file ThermalElement.cpp.

572  {
574 
575  CHKERR mField.add_finite_element("THERMAL_CONVECTION_FE", MF_ZERO);
576  CHKERR mField.modify_finite_element_add_field_row("THERMAL_CONVECTION_FE",
577  field_name);
578  CHKERR mField.modify_finite_element_add_field_col("THERMAL_CONVECTION_FE",
579  field_name);
580  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
581  field_name);
582  if (mField.check_field(mesh_nodals_positions)) {
583  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
584  mesh_nodals_positions);
585  }
586 
587  // this is alternative method for setting boundary conditions, to bypass bu
588  // in cubit file reader. not elegant, but good enough
590  if (it->getName().compare(0, 10, "CONVECTION") == 0) {
591 
592  std::vector<double> data;
593  CHKERR it->getAttributes(data);
594  if (data.size() != 2) {
595  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
596  }
597  setOfConvection[it->getMeshsetId()].cOnvection = data[0];
598  setOfConvection[it->getMeshsetId()].tEmperature = data[1];
599  CHKERR mField.get_moab().get_entities_by_type(
600  it->meshset, MBTRI, setOfConvection[it->getMeshsetId()].tRis, true);
602  setOfConvection[it->getMeshsetId()].tRis, MBTRI,
603  "THERMAL_CONVECTION_FE");
604  }
605  }
606 
608 }
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 moab::Interface & get_moab()=0
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
std::map< int, ConvectionData > setOfConvection
MoFEM::Interface & mField
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#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.
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
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

◆ addThermalElements()

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

add thermal element on tetsIt 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 479 of file ThermalElement.cpp.

480  {
482 
483  CHKERR mField.add_finite_element("THERMAL_FE", MF_ZERO);
484  CHKERR mField.modify_finite_element_add_field_row("THERMAL_FE", field_name);
485  CHKERR mField.modify_finite_element_add_field_col("THERMAL_FE", field_name);
486  CHKERR mField.modify_finite_element_add_field_data("THERMAL_FE", field_name);
487  if (mField.check_field(mesh_nodals_positions)) {
489  mesh_nodals_positions);
490  }
491 
492  // takes skin of block of entities
493  // Skinner skin(&mField.get_moab());
494  // loop over all blocksets and get data which name is FluidPressure
496  mField, BLOCKSET | MAT_THERMALSET, it)) {
497 
498  Mat_Thermal temp_data;
499  ierr = it->getAttributeDataStructure(temp_data);
500 
501  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.resize(
502  3, 3); //(3X3) conductivity matrix
503  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.clear();
504  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(0, 0) =
505  temp_data.data.Conductivity;
506  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(1, 1) =
507  temp_data.data.Conductivity;
508  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(2, 2) =
509  temp_data.data.Conductivity;
510  // setOfBlocks[it->getMeshsetId()].cOnductivity =
511  // temp_data.data.Conductivity;
512 
513  setOfBlocks[it->getMeshsetId()].cApacity = temp_data.data.HeatCapacity;
514  CHKERR mField.get_moab().get_entities_by_type(
515  it->meshset, MBTET, setOfBlocks[it->getMeshsetId()].tEts, true);
517  setOfBlocks[it->getMeshsetId()].tEts, MBTET, "THERMAL_FE");
518  }
519 
521 }
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 moab::Interface & get_moab()=0
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEM::Interface & mField
Thermal material data structure.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
block name is "MAT_THERMAL"
Definition: definitions.h:224
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ addThermalFluxElement()

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

add heat flux elementIt get data from heat flux set and define element in moab. Alternatively uses block set with name HEAT_FLUX.

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

Definition at line 524 of file ThermalElement.cpp.

525  {
527 
528  CHKERR mField.add_finite_element("THERMAL_FLUX_FE", MF_ZERO);
530  field_name);
532  field_name);
534  field_name);
535  if (mField.check_field(mesh_nodals_positions)) {
537  mesh_nodals_positions);
538  }
539 
541  it)) {
542  CHKERR it->getBcDataStructure(setOfFluxes[it->getMeshsetId()].dAta);
543  CHKERR mField.get_moab().get_entities_by_type(
544  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
546  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
547  }
548 
549  // this is alternative method for setting boundary conditions, to bypass bu
550  // in cubit file reader. not elegant, but good enough
552  if (it->getName().compare(0, 9, "HEAT_FLUX") == 0) {
553  std::vector<double> data;
554  CHKERR it->getAttributes(data);
555  if (data.size() != 1) {
556  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
557  }
558  strcpy(setOfFluxes[it->getMeshsetId()].dAta.data.name, "HeatFlu");
559  setOfFluxes[it->getMeshsetId()].dAta.data.flag1 = 1;
560  setOfFluxes[it->getMeshsetId()].dAta.data.value1 = data[0];
561  CHKERR mField.get_moab().get_entities_by_type(
562  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
564  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
565  }
566  }
567 
569 }
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
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 moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEM::Interface & mField
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#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.
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ addThermalRadiationElement()

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

add Non-linear Radiation elementIt get data from Radiation set and define element in moab. Alternatively uses block set with name RADIATION.

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

Definition at line 610 of file ThermalElement.cpp.

611  {
613 
614  CHKERR mField.add_finite_element("THERMAL_RADIATION_FE", MF_ZERO);
615  CHKERR mField.modify_finite_element_add_field_row("THERMAL_RADIATION_FE",
616  field_name);
617  CHKERR mField.modify_finite_element_add_field_col("THERMAL_RADIATION_FE",
618  field_name);
619  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
620  field_name);
621  if (mField.check_field(mesh_nodals_positions)) {
622  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
623  mesh_nodals_positions);
624  }
625 
626  // this is alternative method for setting boundary conditions, to bypass bu
627  // in cubit file reader. not elegant, but good enough
629  if (it->getName().compare(0, 9, "RADIATION") == 0) {
630  std::vector<double> data;
631  ierr = it->getAttributes(data);
632  if (data.size() != 3) {
633  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
634  }
635  setOfRadiation[it->getMeshsetId()].sIgma = data[0];
636  setOfRadiation[it->getMeshsetId()].eMissivity = data[1];
637  setOfRadiation[it->getMeshsetId()].aMbienttEmp = data[2];
638  CHKERR mField.get_moab().get_entities_by_type(
639  it->meshset, MBTRI, setOfRadiation[it->getMeshsetId()].tRis, true);
641  setOfRadiation[it->getMeshsetId()].tRis, MBTRI,
642  "THERMAL_RADIATION_FE");
643  }
644  }
645 
647 }
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 moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEM::Interface & mField
std::map< int, RadiationData > setOfRadiation
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#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.
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

664  {
666  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
667  for (; sit != setOfBlocks.end(); sit++) {
668  // add finite elemen
669  feLhs.getOpPtrVector().push_back(
670  new OpThermalLhs(field_name, A, sit->second, commonData));
671  }
673 }
CommonData commonData
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MyVolumeFE feLhs

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

650  {
652  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
653  feRhs.getOpPtrVector().push_back(
654  new OpGetGradAtGaussPts(field_name, commonData));
655  for (; sit != setOfBlocks.end(); sit++) {
656  // add finite element
657  feRhs.getOpPtrVector().push_back(
658  new OpThermalRhs(field_name, F, sit->second, commonData));
659  }
661 }
CommonData commonData
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MyVolumeFE feRhs
cauclate right hand side for tetrahedral elements
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

676  {
678  bool hoGeometry = false;
679  if (mField.check_field(mesh_nodals_positions)) {
680  hoGeometry = true;
681  }
682  std::map<int, FluxData>::iterator sit = setOfFluxes.begin();
683  for (; sit != setOfFluxes.end(); sit++) {
684  // add finite element
685  feFlux.getOpPtrVector().push_back(
686  new OpHeatFlux(field_name, F, sit->second, hoGeometry));
687  }
689 }
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEM::Interface & mField
virtual bool check_field(const std::string &name) const =0
check if field is in database
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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

◆ 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