v0.9.0
Classes | Public Member Functions | Public Attributes | List of all members
AnalyticalDirichletHelmholtzBC Struct Reference

Analytical Dirichlet boundary conditions. More...

#include <users_modules/helmholtz/src/AnalyticalDirichletHelmholtz.hpp>

Collaboration diagram for AnalyticalDirichletHelmholtzBC:
[legend]

Classes

struct  ApproxField
 finite element to approximate analytical solution on surface More...
 
struct  DirichletBC
 
struct  DirichletBC_LhsOnly
 
struct  DirichletBC_RhsOnly
 

Public Member Functions

 AnalyticalDirichletHelmholtzBC (MoFEM::Interface &m_field)
 
template<typename FUNEVAL >
PetscErrorCode setApproxOps (MoFEM::Interface &m_field, string field_name, Range &tris, boost::shared_ptr< FUNEVAL > funtcion_evaluator, int field_number=0, string nodals_positions="MESH_NODE_POSITIONS")
 
PetscErrorCode initializeProblem (MoFEM::Interface &m_field, string fe, string field, Range &tris, string nodals_positions="MESH_NODE_POSITIONS")
 
PetscErrorCode setProblem (MoFEM::Interface &m_field, string problem)
 
PetscErrorCode solveProblem (MoFEM::Interface &m_field, string problem, string fe, DirichletBC &bc, Range &tris)
 
PetscErrorCode destroyProblem ()
 
PetscErrorCode sortBc (DirichletBC &bc, Range &tris)
 

Public Attributes

ApproxField approxField
 
Mat A
 
Vec D
 
Vec F
 
KSP kspSolver
 

Detailed Description

Analytical Dirichlet boundary conditions.

Definition at line 30 of file AnalyticalDirichletHelmholtz.hpp.

Constructor & Destructor Documentation

◆ AnalyticalDirichletHelmholtzBC()

AnalyticalDirichletHelmholtzBC::AnalyticalDirichletHelmholtzBC ( MoFEM::Interface m_field)

Definition at line 275 of file AnalyticalDirichletHelmholtz.cpp.

Member Function Documentation

◆ destroyProblem()

PetscErrorCode AnalyticalDirichletHelmholtzBC::destroyProblem ( )

Definition at line 449 of file AnalyticalDirichletHelmholtz.cpp.

449  {
450  PetscFunctionBegin;
451  PetscErrorCode ierr;
452  ierr = KSPDestroy(&kspSolver); CHKERRQ(ierr);
453  ierr = MatDestroy(&A); CHKERRQ(ierr);
454  ierr = VecDestroy(&F); CHKERRQ(ierr);
455  ierr = VecDestroy(&D); CHKERRQ(ierr);
456  PetscFunctionReturn(0);
457 
458  }
CHKERRQ(ierr)

◆ initializeProblem()

PetscErrorCode AnalyticalDirichletHelmholtzBC::initializeProblem ( MoFEM::Interface m_field,
string  fe,
string  field,
Range &  tris,
string  nodals_positions = "MESH_NODE_POSITIONS" 
)

Definition at line 277 of file AnalyticalDirichletHelmholtz.cpp.

279  {
280  PetscFunctionBegin;
281  PetscErrorCode ierr;
282  ierr = m_field.add_finite_element(fe,MF_ZERO); CHKERRQ(ierr);
286  if(m_field.check_field(nodals_positions)) {
287  ierr = m_field.modify_finite_element_add_field_data(fe,nodals_positions); CHKERRQ(ierr);
288  }
289  ierr = m_field.add_ents_to_finite_element_by_type(tris,MBTRI,fe); CHKERRQ(ierr);
290  PetscFunctionReturn(0);
291  }
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
CHKERRQ(ierr)
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_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

◆ setApproxOps()

template<typename FUNEVAL >
PetscErrorCode AnalyticalDirichletHelmholtzBC::setApproxOps ( MoFEM::Interface m_field,
string  field_name,
Range &  tris,
boost::shared_ptr< FUNEVAL >  funtcion_evaluator,
int  field_number = 0,
string  nodals_positions = "MESH_NODE_POSITIONS" 
)

Definition at line 242 of file AnalyticalDirichletHelmholtz.hpp.

246  {
247  PetscFunctionBegin;
248  if(approxField.getLoopFeApprox().getOpPtrVector().empty()) {
249  if(m_field.check_field(nodals_positions)) {
250  approxField.getLoopFeApprox().getOpPtrVector().push_back(new ApproxField::OpHoCoord(nodals_positions,approxField.hoCoords));
251  }
252  approxField.getLoopFeApprox().getOpPtrVector().push_back(new ApproxField::OpLhs(field_name,approxField.hoCoords));
253  }
254  approxField.getLoopFeApprox().getOpPtrVector().push_back(new ApproxField::OpRhs<FUNEVAL>(field_name,tris,approxField.hoCoords,funtcion_evaluator,field_number));
255  PetscFunctionReturn(0);
256  }
virtual bool check_field(const std::string &name) const =0
check if field is in database

◆ setProblem()

PetscErrorCode AnalyticalDirichletHelmholtzBC::setProblem ( MoFEM::Interface m_field,
string  problem 
)

Definition at line 391 of file AnalyticalDirichletHelmholtz.cpp.

393  {
394  PetscFunctionBegin;
395  PetscErrorCode ierr;
396 
397  ierr = m_field.VecCreateGhost(problem,ROW,&F); CHKERRQ(ierr);
398  ierr = m_field.VecCreateGhost(problem,COL,&D); CHKERRQ(ierr);
399  ierr = m_field.MatCreateMPIAIJWithArrays(problem,&A); CHKERRQ(ierr);
400 
401  ierr = KSPCreate(PETSC_COMM_WORLD,&kspSolver); CHKERRQ(ierr);
402  ierr = KSPSetOperators(kspSolver,A,A); CHKERRQ(ierr);
403  ierr = KSPSetFromOptions(kspSolver); CHKERRQ(ierr);
404 
405  //PC pc;
406  //ierr = KSPGetPC(kspSolver,&pc); CHKERRQ(ierr);
407  //ierr = PCSetType(pc,PCLU); CHKERRQ(ierr);
408  //ierr = PCFactorSetMatSolverPackage(pc,MATSOLVERMUMPS); CHKERRQ(ierr);
409  //ierr = PCFactorSetUpMatSolverPackage(pc); CHKERRQ(ierr);
410 
411  PetscFunctionReturn(0);
412  }
CHKERRQ(ierr)
DEPRECATED MoFEMErrorCode VecCreateGhost(const std::string &name, RowColData rc, Vec *V) const
create ghost vector for problem (collective)collective - need to be run on all processors in communic...
DEPRECATED MoFEMErrorCode MatCreateMPIAIJWithArrays(const std::string &name, Mat *Aij, int verb=DEFAULT_VERBOSITY)

◆ solveProblem()

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

Definition at line 414 of file AnalyticalDirichletHelmholtz.cpp.

416  {
417  PetscFunctionBegin;
418  PetscErrorCode ierr;
419 
420  ierr = VecZeroEntries(F); CHKERRQ(ierr);
421  ierr = MatZeroEntries(A); CHKERRQ(ierr);
422  ierr = VecGhostUpdateBegin(F,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
423  ierr = VecGhostUpdateEnd(F,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
424 
425  approxField.getLoopFeApprox().snes_B = A;
426  approxField.getLoopFeApprox().snes_f = F;
428  ierr = VecGhostUpdateBegin(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
429  ierr = VecGhostUpdateEnd(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
430  ierr = VecAssemblyBegin(F); CHKERRQ(ierr);
431  ierr = VecAssemblyEnd(F); CHKERRQ(ierr);
432  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
433  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
434 
435  ierr = KSPSolve(kspSolver,F,D); CHKERRQ(ierr);
436  ierr = VecGhostUpdateBegin(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
437  ierr = VecGhostUpdateEnd(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
438 
439  ierr = m_field.set_global_ghost_vector(problem,ROW,D,INSERT_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
440 
441  bc.tRis_ptr = &tris;
442  bc.mapZeroRows.clear();
443  bc.dofsIndices.clear();
444  bc.dofsValues.clear();
445 
446  PetscFunctionReturn(0);
447  }
DEPRECATED MoFEMErrorCode set_global_ghost_vector(const Problem *problem_ptr, RowColData rc, Vec V, InsertMode mode, ScatterMode scatter_mode) const
set values of vector from/to mesh database (collective)collective - need tu be run on all processors ...
CHKERRQ(ierr)
DEPRECATED MoFEMErrorCode loop_finite_elements(const Problem *problem_ptr, const std::string &fe_name, FEMethod &method, int lower_rank, int upper_rank, MoFEMTypes bh, int verb=DEFAULT_VERBOSITY)

◆ sortBc()

PetscErrorCode AnalyticalDirichletHelmholtzBC::sortBc ( DirichletBC bc,
Range &  tris 
)

Definition at line 460 of file AnalyticalDirichletHelmholtz.cpp.

462  {
463  PetscFunctionBegin;
464  PetscErrorCode ierr;
465  bc.tRis_ptr = &tris;
466  bc.mapZeroRows.clear();
467  bc.dofsIndices.clear();
468  bc.dofsValues.clear();
469  PetscFunctionReturn(0);
470 
471  }

Member Data Documentation

◆ A

Mat AnalyticalDirichletHelmholtzBC::A

Definition at line 263 of file AnalyticalDirichletHelmholtz.hpp.

◆ approxField

ApproxField AnalyticalDirichletHelmholtzBC::approxField

Definition at line 238 of file AnalyticalDirichletHelmholtz.hpp.

◆ D

Vec AnalyticalDirichletHelmholtzBC::D

Definition at line 264 of file AnalyticalDirichletHelmholtz.hpp.

◆ F

Vec AnalyticalDirichletHelmholtzBC::F

Definition at line 264 of file AnalyticalDirichletHelmholtz.hpp.

◆ kspSolver

KSP AnalyticalDirichletHelmholtzBC::kspSolver

Definition at line 265 of file AnalyticalDirichletHelmholtz.hpp.


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