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

Analytical Dirichlet boundary conditions. More...

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

Collaboration diagram for AnalyticalDirichletBC:
[legend]

Classes

struct  ApproxField
 finite element to approximate analytical solution on surface More...
 
struct  DirichletBC
 Structure used to enforce analytical boundary conditions. More...
 

Public Member Functions

 AnalyticalDirichletBC (MoFEM::Interface &m_field)
 
template<typename FUNEVAL >
MoFEMErrorCode setApproxOps (MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< FUNEVAL > function_evaluator, const int field_number=0, const string nodals_positions="MESH_NODE_POSITIONS")
 Set operators used to calculate the rhs vector and the lhs matrix. More...
 
MoFEMErrorCode setFiniteElement (MoFEM::Interface &m_field, string fe, string field, Range &tris, string nodals_positions="MESH_NODE_POSITIONS")
 set finite element More...
 
MoFEMErrorCode setUpProblem (MoFEM::Interface &m_field, string problem)
 set problem solver and create matrices and vectors More...
 
MoFEMErrorCode solveProblem (MoFEM::Interface &m_field, string problem, string fe, DirichletBC &bc, Range &tris)
 solve boundary problem More...
 
MoFEMErrorCode solveProblem (MoFEM::Interface &m_field, string problem, string fe, DirichletBC &bc)
 solve boundary problem More...
 
MoFEMErrorCode destroyProblem ()
 Destroy problem. More...
 

Public Attributes

ApproxField approxField
 
Mat A
 
Vec D
 
Vec F
 
KSP kspSolver
 

Detailed Description

Analytical Dirichlet boundary conditions.

Definition at line 30 of file AnalyticalDirichlet.hpp.

Constructor & Destructor Documentation

◆ AnalyticalDirichletBC()

AnalyticalDirichletBC::AnalyticalDirichletBC ( MoFEM::Interface m_field)

Definition at line 227 of file AnalyticalDirichlet.cpp.

228  : approxField(m_field){};

Member Function Documentation

◆ destroyProblem()

MoFEMErrorCode AnalyticalDirichletBC::destroyProblem ( )

Destroy problem.

Destroy matrices and vectors used to solve boundary problem, i.e. finding values of DOFs which approximate analytical solution.

Returns
error code

Definition at line 307 of file AnalyticalDirichlet.cpp.

307  {
309  CHKERR KSPDestroy(&kspSolver);
310  CHKERR MatDestroy(&A);
311  CHKERR VecDestroy(&F);
312  CHKERR VecDestroy(&D);
314 }
#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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ setApproxOps()

template<typename FUNEVAL >
MoFEMErrorCode AnalyticalDirichletBC::setApproxOps ( MoFEM::Interface m_field,
const std::string  field_name,
boost::shared_ptr< FUNEVAL >  function_evaluator,
const int  field_number = 0,
const string  nodals_positions = "MESH_NODE_POSITIONS" 
)

Set operators used to calculate the rhs vector and the lhs matrix.

To enforce analytical function on boundary, first function has to be approximated by finite element base functions. This is done by solving system of linear equations, Following function set finite element operators to calculate the left hand side matrix and the left hand side matrix.

Parameters
m_fieldinterface
field_namefield name
function_evaluatoranalytical function to evaluate
field_numberfield index
nodals_positionsname of the field for ho-geometry description
Returns
error code

Definition at line 193 of file AnalyticalDirichlet.hpp.

196  {
198  if (approxField.getLoopFeApprox().getOpPtrVector().empty()) {
199  if (m_field.check_field(nodals_positions)) {
201  new ApproxField::OpHoCoord(nodals_positions, approxField.hoCoords));
202  }
204  new ApproxField::OpLhs(field_name, approxField.hoCoords));
205  }
207  new ApproxField::OpRhs<FUNEVAL>(field_name, approxField.hoCoords,
208  function_evaluator, field_number));
210  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:500
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:507
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
virtual bool check_field(const std::string &name) const =0
check if field is in database

◆ setFiniteElement()

MoFEMErrorCode AnalyticalDirichletBC::setFiniteElement ( MoFEM::Interface m_field,
string  fe,
string  field,
Range &  tris,
string  nodals_positions = "MESH_NODE_POSITIONS" 
)

set finite element

Deprecated:
no need to use function with argument of triangle range
Parameters
m_fieldmofem interface
fefinite element name
fieldfield name
trisfaces where analytical boundary is given
nodals_positionsfield having higher order geometry description
Returns
error code

Definition at line 231 of file AnalyticalDirichlet.cpp.

233  {
235 
236  CHKERR m_field.add_finite_element(fe, MF_ZERO);
237  CHKERR m_field.modify_finite_element_add_field_row(fe, field);
238  CHKERR m_field.modify_finite_element_add_field_col(fe, field);
239  CHKERR m_field.modify_finite_element_add_field_data(fe, field);
240  if (m_field.check_field(nodals_positions)) {
241  CHKERR m_field.modify_finite_element_add_field_data(fe, nodals_positions);
242  }
243  CHKERR m_field.add_ents_to_finite_element_by_type(tris, MBTRI, fe);
245 }
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
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
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:595
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:406

◆ setUpProblem()

MoFEMErrorCode AnalyticalDirichletBC::setUpProblem ( MoFEM::Interface m_field,
string  problem 
)

set problem solver and create matrices and vectors

Parameters
m_fieldmofem interface
problemproblem name
Returns
error code

Definition at line 247 of file AnalyticalDirichlet.cpp.

248  {
250  CHKERR m_field.getInterface<VecManager>()->vecCreateGhost(problem, ROW, &F);
251  CHKERR m_field.getInterface<VecManager>()->vecCreateGhost(problem, COL, &D);
253  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem, &A);
254 
255  CHKERR KSPCreate(PETSC_COMM_WORLD, &kspSolver);
256  CHKERR KSPSetOperators(kspSolver, A, A);
257  CHKERR KSPSetFromOptions(kspSolver);
259 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:36
#define CHKERR
Inline error check.
Definition: definitions.h:595
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
Matrix manager is used to build and partition problems.

◆ solveProblem() [1/2]

MoFEMErrorCode AnalyticalDirichletBC::solveProblem ( MoFEM::Interface m_field,
string  problem,
string  fe,
DirichletBC bc,
Range &  tris 
)

solve boundary problem

Deprecated:
use setUpProblem instead

This functions solve for DOFs on boundary where analytical solution is given. i.e. finding values of DOFs which approximate analytical solution.

Parameters
m_fieldmofem interface
problemproblem name
fefinite element name
bcDriblet boundary structure used to apply boundary conditions
tristriangles on boundary
Returns
error code

◆ solveProblem() [2/2]

MoFEMErrorCode AnalyticalDirichletBC::solveProblem ( MoFEM::Interface m_field,
string  problem,
string  fe,
DirichletBC bc 
)

solve boundary problem

This functions solve for DOFs on boundary where analytical solution is given.

Parameters
m_fieldmofem interface
problemproblem name
fefinite element name
bcDriblet boundary structure used to apply boundary conditions
Returns
[description]

Member Data Documentation

◆ A

Mat AnalyticalDirichletBC::A
Deprecated:
use setFiniteElement instead

Definition at line 255 of file AnalyticalDirichlet.hpp.

◆ approxField

ApproxField AnalyticalDirichletBC::approxField

Definition at line 172 of file AnalyticalDirichlet.hpp.

◆ D

Vec AnalyticalDirichletBC::D

Definition at line 256 of file AnalyticalDirichlet.hpp.

◆ F

Vec AnalyticalDirichletBC::F

Definition at line 256 of file AnalyticalDirichlet.hpp.

◆ kspSolver

KSP AnalyticalDirichletBC::kspSolver

Definition at line 257 of file AnalyticalDirichlet.hpp.


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