v0.8.23
Public Member Functions | Public Attributes | List of all members
DirichletFixFieldAtEntitiesBc Struct Reference

Fix dofs on entities. More...

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

Inheritance diagram for DirichletFixFieldAtEntitiesBc:
[legend]
Collaboration diagram for DirichletFixFieldAtEntitiesBc:
[legend]

Public Member Functions

 DirichletFixFieldAtEntitiesBc (MoFEM::Interface &m_field, const std::string &field_name, Mat aij, Vec x, Vec f, Range &ents)
 
 DirichletFixFieldAtEntitiesBc (MoFEM::Interface &m_field, const std::string &field_name, Range &ents)
 
MoFEMErrorCode iNitalize ()
 
MoFEMErrorCode preProcess ()
 function is run at the beginning of loop More...
 
MoFEMErrorCode postProcess ()
 function is run at the end of loop More...
 
- Public Member Functions inherited from DirichletDisplacementBc
 DirichletDisplacementBc (MoFEM::Interface &m_field, const std::string &field_name, Mat Aij, Vec X, Vec F)
 
 DirichletDisplacementBc (MoFEM::Interface &m_field, const std::string &field_name)
 
- Public Member Functions inherited from MoFEM::FEMethod
MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 FEMethod ()
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name, const EntityType type, const int side_number) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name, const EntityType type, const int side_number) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const EntityHandle ent) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const EntityHandle ent) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name, const EntityHandle ent) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name, const EntityHandle ent) const
 
- Public Member Functions inherited from MoFEM::BasicMethod
 BasicMethod ()
 
virtual ~BasicMethod ()
 
int getNinTheLoop () const
 get number of evaluated element in the loop More...
 
int getLoopSize () const
 get loop size More...
 
MoFEMErrorCode copyBasicMethod (const BasicMethod &basic)
 Copy data from other base method to this base method. More...
 
virtual MoFEMErrorCode operator() ()
 function is run for every finite element More...
 
- Public Member Functions inherited from MoFEM::KspMethod
 KspMethod ()
 
virtual ~KspMethod ()
 
MoFEMErrorCode setKspCtx (const KSPContext &ctx)
 set operator type More...
 
MoFEMErrorCode setKsp (KSP ksp)
 set solver More...
 
MoFEMErrorCode copyKsp (const KspMethod &ksp)
 copy data form another method More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 
- Public Member Functions inherited from MoFEM::SnesMethod
 SnesMethod ()
 
virtual ~SnesMethod ()
 
MoFEMErrorCode setSnesCtx (const SNESContext &ctx)
 Set SNES context. More...
 
MoFEMErrorCode setSnes (SNES snes)
 Set SNES instance. More...
 
MoFEMErrorCode copySnes (const SnesMethod &snes)
 Copy snes data. More...
 
- Public Member Functions inherited from MoFEM::TSMethod
 TSMethod ()
 
virtual ~TSMethod ()
 
MoFEMErrorCode setTsCtx (const TSContext &ctx)
 Set Ts context. More...
 
MoFEMErrorCode copyTs (const TSMethod &ts)
 Copy TS solver data. More...
 
MoFEMErrorCode setTs (TS _ts)
 Set TS solver. More...
 

Public Attributes

Range eNts
 
std::vector< std::string > fieldNames
 
- Public Attributes inherited from DirichletDisplacementBc
MoFEM::InterfacemField
 
const std::string fieldName
 field name to set Dirichlet BC More...
 
double dIag
 diagonal value set on zeroed column and rows More...
 
std::map< DofIdx, FieldData > mapZeroRows
 
std::vector< intdofsIndices
 
std::vector< doubledofsValues
 
std::vector< doubledofsXValues
 
boost::ptr_vector< MethodForForceScalingmethodsOp
 
- Public Attributes inherited from MoFEM::FEMethod
std::string feName
 Name of finite element. More...
 
boost::shared_ptr< const NumeredEntFiniteElementnumeredEntFiniteElementPtr
 
boost::shared_ptr< const FENumeredDofEntity_multiIndexrowPtr
 Pointer to finite element rows dofs view. More...
 
boost::shared_ptr< const FENumeredDofEntity_multiIndexcolPtr
 Pointer to finite element columns dofs view. More...
 
boost::shared_ptr< const FEDofEntity_multiIndexdataPtr
 Pointer to finite element data dofs. More...
 
boost::shared_ptr< const FieldEntity_vector_viewrowFieldEntsPtr
 Pointer to finite element field entities row view. More...
 
boost::shared_ptr< const FieldEntity_vector_viewcolFieldEntsPtr
 Pointer to finite element field entities column view. More...
 
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_viewdataFieldEntsPtr
 Pointer to finite element field entities data view. More...
 
- Public Attributes inherited from MoFEM::BasicMethod
int nInTheLoop
 number currently of processed method More...
 
int loopSize
 local number oe methods to process More...
 
int rAnk
 processor rank More...
 
int sIze
 number of processors in communicator More...
 
const RefEntity_multiIndexrefinedEntitiesPtr
 container of mofem dof entities More...
 
const RefElement_multiIndexrefinedFiniteElementsPtr
 container of mofem finite element entities More...
 
const ProblemproblemPtr
 raw pointer to problem More...
 
const Field_multiIndexfieldsPtr
 raw pointer to fields container More...
 
const FieldEntity_multiIndexentitiesPtr
 raw pointer to container of field entities More...
 
const DofEntity_multiIndexdofsPtr
 raw pointer container of dofs More...
 
const FiniteElement_multiIndexfiniteElementsPtr
 raw pointer to container finite elements More...
 
const EntFiniteElement_multiIndexfiniteElementsEntitiesPtr
 
const FieldEntityEntFiniteElementAdjacencyMap_multiIndexadjacenciesPtr
 
boost::function< MoFEMErrorCode()> preProcessHook
 Hook function for pre-processing. More...
 
boost::function< MoFEMErrorCode()> postProcessHook
 Hook function for post-processing. More...
 
boost::function< MoFEMErrorCode()> operatorHook
 Hook function for operator. More...
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 
- Public Attributes inherited from MoFEM::KspMethod
KSPContext ksp_ctx
 Context. More...
 
KSP ksp
 KSP solver. More...
 
Vec ksp_f
 the right hand side vector More...
 
Mat ksp_A
 matrix More...
 
Mat ksp_B
 preconditioner matrix More...
 
- Public Attributes inherited from MoFEM::SnesMethod
SNESContext snes_ctx
 
SNES snes
 
Vec snes_x
 
Vec snes_f
 
Mat snes_A
 
Mat snes_B
 
- Public Attributes inherited from MoFEM::TSMethod
TSContext ts_ctx
 
TS ts
 
Vec ts_u
 
Vec ts_u_t
 
Vec ts_F
 
Mat ts_A
 
Mat ts_B
 
PetscInt ts_step
 
PetscReal ts_a
 
PetscReal ts_t
 

Additional Inherited Members

- Public Types inherited from MoFEM::KspMethod
enum  KSPContext { CTX_SETFUNCTION, CTX_OPERATORS, CTX_KSPNONE }
 pass information about context of KSP/DM for with finite element is computed More...
 
- Public Types inherited from MoFEM::SnesMethod
enum  SNESContext { CTX_SNESSETFUNCTION, CTX_SNESSETJACOBIAN, CTX_SNESNONE }
 
- Public Types inherited from MoFEM::TSMethod
enum  TSContext {
  CTX_TSSETRHSFUNCTION, CTX_TSSETRHSJACOBIAN, CTX_TSSETIFUNCTION, CTX_TSSETIJACOBIAN,
  CTX_TSTSMONITORSET, CTX_TSNONE
}
 
- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Fix dofs on entities.

Examples
mesh_smoothing.cpp, and minimal_surface_area.cpp.

Definition at line 134 of file DirichletBC.hpp.

Constructor & Destructor Documentation

◆ DirichletFixFieldAtEntitiesBc() [1/2]

DirichletFixFieldAtEntitiesBc::DirichletFixFieldAtEntitiesBc ( MoFEM::Interface m_field,
const std::string &  field_name,
Mat  aij,
Vec  x,
Vec  f,
Range &  ents 
)

Definition at line 138 of file DirichletBC.hpp.

143  :
144  DirichletDisplacementBc(m_field,field_name,aij,x,f),
145  eNts(ents) {
146  fieldNames.push_back(fieldName);
147  }
const std::string fieldName
field name to set Dirichlet BC
Definition: DirichletBC.hpp:46
std::vector< std::string > fieldNames
DirichletDisplacementBc(MoFEM::Interface &m_field, const std::string &field_name, Mat Aij, Vec X, Vec F)
Definition: DirichletBC.cpp:50

◆ DirichletFixFieldAtEntitiesBc() [2/2]

DirichletFixFieldAtEntitiesBc::DirichletFixFieldAtEntitiesBc ( MoFEM::Interface m_field,
const std::string &  field_name,
Range &  ents 
)

Definition at line 149 of file DirichletBC.hpp.

151  :
152  DirichletDisplacementBc(m_field,field_name),eNts(ents) {
153  fieldNames.push_back(fieldName);
154  }
const std::string fieldName
field name to set Dirichlet BC
Definition: DirichletBC.hpp:46
std::vector< std::string > fieldNames
DirichletDisplacementBc(MoFEM::Interface &m_field, const std::string &field_name, Mat Aij, Vec X, Vec F)
Definition: DirichletBC.cpp:50

Member Function Documentation

◆ iNitalize()

MoFEMErrorCode DirichletFixFieldAtEntitiesBc::iNitalize ( )
virtual

Reimplemented from DirichletDisplacementBc.

Definition at line 516 of file DirichletBC.cpp.

516  {
518  if (mapZeroRows.empty()) {
519 
520  for (const auto &field_name : fieldNames) {
521 
522  auto for_each_dof = [&](auto &dof) {
524  mapZeroRows[dof->getPetscGlobalDofIdx()] = 0;
526  };
527 
529  for_each_dof);
530  }
531 
532  dofsIndices.resize(mapZeroRows.size());
533  dofsValues.resize(mapZeroRows.size());
534  int ii = 0;
535  for (auto mit = mapZeroRows.begin(); mit != mapZeroRows.end();
536  ++mit, ++ii) {
537  dofsIndices[ii] = mit->first;
538  dofsValues[ii] = mit->second;
539  }
540 
541  }
543 }
std::vector< std::string > fieldNames
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:500
std::map< DofIdx, FieldData > mapZeroRows
Definition: DirichletBC.hpp:57
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:507
const Problem * problemPtr
raw pointer to problem
std::vector< int > dofsIndices
Definition: DirichletBC.hpp:58
std::vector< double > dofsValues
Definition: DirichletBC.hpp:59
#define CHKERR
Inline error check.
Definition: definitions.h:595
static MoFEMErrorCode set_numered_dofs_on_ents(const Problem *problem_ptr, const string &field_name, Range &ents, boost::function< MoFEMErrorCode(const boost::shared_ptr< MoFEM::NumeredDofEntity > &dof)> for_each_dof)
Definition: DirichletBC.cpp:23
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ postProcess()

MoFEMErrorCode DirichletFixFieldAtEntitiesBc::postProcess ( )
virtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from DirichletDisplacementBc.

Definition at line 568 of file DirichletBC.cpp.

568  {
570  if (snes_ctx == CTX_SNESNONE && ts_ctx == CTX_TSNONE) {
571  if (snes_B) {
572  CHKERR MatAssemblyBegin(snes_B, MAT_FINAL_ASSEMBLY);
573  CHKERR MatAssemblyEnd(snes_B, MAT_FINAL_ASSEMBLY);
574  CHKERR MatZeroRowsColumns(snes_B, dofsIndices.size(),
575  dofsIndices.empty() ? PETSC_NULL
576  : &*dofsIndices.begin(),
577  dIag, PETSC_NULL, PETSC_NULL);
578  }
579  if (snes_f) {
580  CHKERR VecAssemblyBegin(snes_f);
581  CHKERR VecAssemblyEnd(snes_f);
582  int ii = 0;
583  for (std::vector<int>::iterator vit = dofsIndices.begin();
584  vit != dofsIndices.end(); vit++, ii++) {
585  CHKERR VecSetValue(snes_f, *vit, dofsValues[ii], INSERT_VALUES);
586  }
587  CHKERR VecAssemblyBegin(snes_f);
588  CHKERR VecAssemblyEnd(snes_f);
589  }
590  }
591 
592  switch (snes_ctx) {
593  case CTX_SNESNONE: {
594  } break;
595  case CTX_SNESSETFUNCTION: {
596  CHKERR VecAssemblyBegin(snes_f);
597  CHKERR VecAssemblyEnd(snes_f);
598  int ii = 0;
599  for (std::vector<int>::iterator vit = dofsIndices.begin();
600  vit != dofsIndices.end(); vit++, ii++) {
601  CHKERR VecSetValue(snes_f, *vit, dofsValues[ii], INSERT_VALUES);
602  }
603  CHKERR VecAssemblyBegin(snes_f);
604  CHKERR VecAssemblyEnd(snes_f);
605  } break;
606  case CTX_SNESSETJACOBIAN: {
607  CHKERR MatAssemblyBegin(snes_B, MAT_FINAL_ASSEMBLY);
608  CHKERR MatAssemblyEnd(snes_B, MAT_FINAL_ASSEMBLY);
609  CHKERR MatZeroRowsColumns(snes_B, dofsIndices.size(),
610  dofsIndices.empty() ? PETSC_NULL
611  : &*dofsIndices.begin(),
612  dIag, PETSC_NULL, PETSC_NULL);
613  } break;
614  default:
615  SETERRQ(PETSC_COMM_SELF, 1, "unknown snes stage");
616  }
617 
619 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
std::vector< int > dofsIndices
Definition: DirichletBC.hpp:58
double dIag
diagonal value set on zeroed column and rows
Definition: DirichletBC.hpp:47
std::vector< double > dofsValues
Definition: DirichletBC.hpp:59
#define CHKERR
Inline error check.
Definition: definitions.h:595
SNESContext snes_ctx
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
TSContext ts_ctx

◆ preProcess()

MoFEMErrorCode DirichletFixFieldAtEntitiesBc::preProcess ( )
virtual

function is run at the beginning of loop

It is used to zeroing matrices and vectors, calculation of shape functions on reference element, preprocessing boundary conditions, etc.

Reimplemented from DirichletDisplacementBc.

Definition at line 545 of file DirichletBC.cpp.

545  {
547 
548  switch (ts_ctx) {
549  case CTX_TSSETIFUNCTION: {
551  snes_x = ts_u;
552  snes_f = ts_F;
553  break;
554  }
555  case CTX_TSSETIJACOBIAN: {
557  snes_B = ts_B;
558  break;
559  }
560  default:
561  break;
562  }
563 
564  CHKERR iNitalize();
566 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
#define CHKERR
Inline error check.
Definition: definitions.h:595
SNESContext snes_ctx
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
TSContext ts_ctx

Member Data Documentation

◆ eNts

Range DirichletFixFieldAtEntitiesBc::eNts

Definition at line 136 of file DirichletBC.hpp.

◆ fieldNames

std::vector<std::string> DirichletFixFieldAtEntitiesBc::fieldNames

Definition at line 137 of file DirichletBC.hpp.


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