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

Analytical Dirichlet boundary conditions. More...

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

Collaboration diagram for AnalyticalDirihletBC:
[legend]

Classes

struct  ApproxField
 finite element to appeximate analytical solution on surface More...
 
struct  DirichletBC
 

Public Member Functions

 AnalyticalDirihletBC (MoFEM::Interface &m_field, Range &bc_tris)
 
PetscErrorCode setApproxOps (MoFEM::Interface &m_field, string re_field_name, double(*fun)(double x, double y, double z, bool use_real), bool useReal, string nodals_positions="MESH_NODE_POSITIONS")
 
PetscErrorCode initializeProblem (MoFEM::Interface &m_field, string problem, string fe, string re_field_name, 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)
 
PetscErrorCode destroyProblem ()
 
 AnalyticalDirihletBC (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 ()
 

Public Attributes

ApproxField approxField
 
Range tRis
 
Mat A
 
Vec D
 
Vec F
 
KSP solver
 
KSP kspSolver
 

Detailed Description

Analytical Dirichlet boundary conditions.

Definition at line 27 of file AnalyticalDirichletObsolete.hpp.

Constructor & Destructor Documentation

◆ AnalyticalDirihletBC() [1/2]

AnalyticalDirihletBC::AnalyticalDirihletBC ( MoFEM::Interface m_field,
Range &  bc_tris 
)

Definition at line 369 of file AnalyticalDirichletObsolete.hpp.

◆ AnalyticalDirihletBC() [2/2]

AnalyticalDirihletBC::AnalyticalDirihletBC ( MoFEM::Interface m_field)

Definition at line 370 of file AnalyticalDirihlet_obsolete.hpp.

370 : approxField(m_field) {};

Member Function Documentation

◆ destroyProblem() [1/2]

PetscErrorCode AnalyticalDirihletBC::destroyProblem ( )

Definition at line 463 of file AnalyticalDirihlet_obsolete.hpp.

463  {
464  PetscFunctionBegin;
465  PetscErrorCode ierr;
466  ierr = KSPDestroy(&kspSolver); CHKERRQ(ierr);
467  ierr = MatDestroy(&A); CHKERRQ(ierr);
468  ierr = VecDestroy(&F); CHKERRQ(ierr);
469  ierr = VecDestroy(&D); CHKERRQ(ierr);
470  PetscFunctionReturn(0);
471  }
CHKERRQ(ierr)

◆ destroyProblem() [2/2]

PetscErrorCode AnalyticalDirihletBC::destroyProblem ( )

Definition at line 473 of file AnalyticalDirichletObsolete.hpp.

473  {
474  PetscFunctionBegin;
475  PetscErrorCode ierr;
476  ierr = KSPDestroy(&solver); CHKERRQ(ierr);
477  ierr = MatDestroy(&A); CHKERRQ(ierr);
478  ierr = VecDestroy(&F); CHKERRQ(ierr);
479  ierr = VecDestroy(&D); CHKERRQ(ierr);
480  PetscFunctionReturn(0);
481  }
CHKERRQ(ierr)

◆ initializeProblem() [1/2]

PetscErrorCode AnalyticalDirihletBC::initializeProblem ( MoFEM::Interface m_field,
string  problem,
string  fe,
string  re_field_name,
string  nodals_positions = "MESH_NODE_POSITIONS" 
)

Definition at line 389 of file AnalyticalDirichletObsolete.hpp.

392  {
393  PetscFunctionBegin;
394  PetscErrorCode ierr;
395  //Add triangle elements
396  ierr = m_field.add_finite_element(fe,MF_ZERO); CHKERRQ(ierr);
397  ierr = m_field.modify_finite_element_add_field_row(fe,re_field_name); CHKERRQ(ierr);
398  ierr = m_field.modify_finite_element_add_field_col(fe,re_field_name); CHKERRQ(ierr);
399  ierr = m_field.modify_finite_element_add_field_data(fe,re_field_name); CHKERRQ(ierr);
400 
401 
402  if(m_field.check_field(nodals_positions)) {
403  ierr = m_field.modify_finite_element_add_field_data(fe,nodals_positions); CHKERRQ(ierr);
404  }
405  ierr = m_field.modify_problem_add_finite_element(problem,fe); CHKERRQ(ierr);
407 
408  PetscFunctionReturn(0);
409  }
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 modify_problem_add_finite_element(const std::string &name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
virtual DEPRECATED MoFEMErrorCode add_ents_to_finite_element_by_TRIs(const Range &tris, const std::string &name)=0
add TRI entities from range to finite element database given by name
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

◆ initializeProblem() [2/2]

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

Definition at line 389 of file AnalyticalDirihlet_obsolete.hpp.

391  {
392  PetscFunctionBegin;
393  PetscErrorCode ierr;
394  ierr = m_field.add_finite_element(fe,MF_ZERO); CHKERRQ(ierr);
399  if(m_field.check_field(nodals_positions)) {
400  ierr = m_field.modify_finite_element_add_field_data(fe,nodals_positions); CHKERRQ(ierr);
401  }
403  PetscFunctionReturn(0);
404  }
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 DEPRECATED MoFEMErrorCode add_ents_to_finite_element_by_TRIs(const Range &tris, const std::string &name)=0
add TRI entities from range to finite element database given by name
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() [1/2]

PetscErrorCode AnalyticalDirihletBC::setApproxOps ( MoFEM::Interface m_field,
string  re_field_name,
double(*)(double x, double y, double z, bool use_real)  fun,
bool  useReal,
string  nodals_positions = "MESH_NODE_POSITIONS" 
)

Definition at line 371 of file AnalyticalDirichletObsolete.hpp.

376  {
377  PetscFunctionBegin;
378  if(m_field.check_field(nodals_positions)) {
379  approxField.getLoopFeApproxTri().get_op_to_do_Rhs().push_back(new ApproxField::OpHoCoordTri(nodals_positions,approxField.hoCoordsTri));
380  }
381  //loop over triangles
382  approxField.getLoopFeApproxTri().get_op_to_do_Rhs().push_back(new ApproxField::OpRhsTri(re_field_name,approxField.hoCoordsTri,fun,useReal));
383  approxField.getLoopFeApproxTri().get_op_to_do_Lhs().push_back(new ApproxField::OpLhsTri(re_field_name,approxField.hoCoordsTri));
384  //loop over tets
385  PetscFunctionReturn(0);
386  }
virtual bool check_field(const std::string &name) const =0
check if field is in database

◆ setApproxOps() [2/2]

template<typename FUNEVAL >
PetscErrorCode AnalyticalDirihletBC::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 373 of file AnalyticalDirihlet_obsolete.hpp.

377  {
378  PetscFunctionBegin;
379  if(approxField.getLoopFeApprox().get_op_to_do_Rhs().empty()) {
380  if(m_field.check_field(nodals_positions)) {
381  approxField.getLoopFeApprox().get_op_to_do_Rhs().push_back(new ApproxField::OpHoCoord(nodals_positions,approxField.hoCoords));
382  }
383  approxField.getLoopFeApprox().get_op_to_do_Lhs().push_back(new ApproxField::OpLhs(field_name,approxField.hoCoords));
384  }
385  approxField.getLoopFeApprox().get_op_to_do_Rhs().push_back(new ApproxField::OpRhs<FUNEVAL>(field_name,tris,approxField.hoCoords,funtcion_evaluator,field_number));
386  PetscFunctionReturn(0);
387  }
virtual bool check_field(const std::string &name) const =0
check if field is in database

◆ setProblem() [1/2]

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

Definition at line 409 of file AnalyticalDirihlet_obsolete.hpp.

410  {
411  PetscFunctionBegin;
412  PetscErrorCode ierr;
413 
414  ierr = m_field.VecCreateGhost(problem,ROW,&F); CHKERRQ(ierr);
415  ierr = m_field.VecCreateGhost(problem,COL,&D); CHKERRQ(ierr);
416  ierr = m_field.MatCreateMPIAIJWithArrays(problem,&A); CHKERRQ(ierr);
417 
418  ierr = KSPCreate(PETSC_COMM_WORLD,&kspSolver); CHKERRQ(ierr);
419  ierr = KSPSetOperators(kspSolver,A,A); CHKERRQ(ierr);
420  ierr = KSPSetFromOptions(kspSolver); CHKERRQ(ierr);
421 
422  PC pc;
423  ierr = KSPGetPC(kspSolver,&pc); CHKERRQ(ierr);
424  ierr = PCSetType(pc,PCLU); CHKERRQ(ierr);
425  ierr = PCFactorSetMatSolverPackage(pc,MATSOLVERMUMPS); CHKERRQ(ierr);
426  ierr = PCFactorSetUpMatSolverPackage(pc); CHKERRQ(ierr);
427 
428  PetscFunctionReturn(0);
429  }
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)

◆ setProblem() [2/2]

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

Definition at line 416 of file AnalyticalDirichletObsolete.hpp.

417  {
418  PetscFunctionBegin;
419  PetscErrorCode ierr;
420 
421  ierr = m_field.VecCreateGhost(problem,ROW,&F); CHKERRQ(ierr);
422  ierr = m_field.VecCreateGhost(problem,COL,&D); CHKERRQ(ierr);
423  ierr = m_field.MatCreateMPIAIJWithArrays(problem,&A); CHKERRQ(ierr);
424 
425  ierr = KSPCreate(PETSC_COMM_WORLD,&solver); CHKERRQ(ierr);
426  ierr = KSPSetOperators(solver,A,A); CHKERRQ(ierr);
427  ierr = KSPSetFromOptions(solver); CHKERRQ(ierr);
428 
429  PetscFunctionReturn(0);
430  }
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() [1/2]

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

Definition at line 431 of file AnalyticalDirihlet_obsolete.hpp.

432  {
433  PetscFunctionBegin;
434  PetscErrorCode ierr;
435 
436  ierr = VecZeroEntries(F); CHKERRQ(ierr);
437  ierr = MatZeroEntries(A); CHKERRQ(ierr);
438 
442  ierr = VecGhostUpdateBegin(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
443  ierr = VecGhostUpdateEnd(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
444  ierr = VecAssemblyBegin(F); CHKERRQ(ierr);
445  ierr = VecAssemblyEnd(F); CHKERRQ(ierr);
446  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
447  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
448 
449  ierr = KSPSolve(kspSolver,F,D); CHKERRQ(ierr);
450  ierr = VecGhostUpdateBegin(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
451  ierr = VecGhostUpdateEnd(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
452 
453  ierr = m_field.set_global_ghost_vector(problem,ROW,D,INSERT_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
454 
455  bc.tRis_ptr = &tris;
456  bc.map_zero_rows.clear();
457  bc.dofsIndices.clear();
458  bc.dofsValues.clear();
459 
460  PetscFunctionReturn(0);
461  }
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 ...
Vec snes_f
residual
CHKERRQ(ierr)
Mat snes_B
preconditioner of jacobian matrix
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)

◆ solveProblem() [2/2]

PetscErrorCode AnalyticalDirihletBC::solveProblem ( MoFEM::Interface m_field,
string  problem,
string  fe,
DirichletBC bc 
)

Definition at line 432 of file AnalyticalDirichletObsolete.hpp.

433  {
434  PetscFunctionBegin;
435  PetscErrorCode ierr;
436 
437  ierr = VecZeroEntries(F); CHKERRQ(ierr);
438  ierr = MatZeroEntries(A); CHKERRQ(ierr);
439 
443 
444  ierr = VecGhostUpdateBegin(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
445  ierr = VecGhostUpdateEnd(F,ADD_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
446  ierr = VecAssemblyBegin(F); CHKERRQ(ierr);
447  ierr = VecAssemblyEnd(F); CHKERRQ(ierr);
448  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
449  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
450 
451  //int ii1,jj1;
452  //ierr=MatGetSize(A,&ii1,&jj1);
453 
454  //std::string wait;
455  //ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);
456  //ierr = MatView(A,PETSC_VIEWER_DRAW_WORLD);
457  //std::cin >> wait;
458 
459  ierr = KSPSolve(solver,F,D); CHKERRQ(ierr);
460  ierr = VecGhostUpdateBegin(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
461  ierr = VecGhostUpdateEnd(D,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr);
462 
463  ierr = m_field.set_global_ghost_vector(problem,ROW,D,INSERT_VALUES,SCATTER_REVERSE); CHKERRQ(ierr);
464 
465  bc.tRis_ptr = &tRis;
466  bc.map_zero_rows.clear();
467  bc.dofsIndices.clear();
468  bc.dofsValues.clear();
469 
470  PetscFunctionReturn(0);
471  }
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 ...
Vec snes_f
residual
CHKERRQ(ierr)
Mat snes_B
preconditioner of jacobian matrix
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)

Member Data Documentation

◆ A

Mat AnalyticalDirihletBC::A

Definition at line 412 of file AnalyticalDirichletObsolete.hpp.

◆ approxField

ApproxField AnalyticalDirihletBC::approxField

Definition at line 367 of file AnalyticalDirichletObsolete.hpp.

◆ D

Vec AnalyticalDirihletBC::D

Definition at line 413 of file AnalyticalDirichletObsolete.hpp.

◆ F

Vec AnalyticalDirihletBC::F

Definition at line 413 of file AnalyticalDirichletObsolete.hpp.

◆ kspSolver

KSP AnalyticalDirihletBC::kspSolver

Definition at line 408 of file AnalyticalDirihlet_obsolete.hpp.

◆ solver

KSP AnalyticalDirihletBC::solver

Definition at line 414 of file AnalyticalDirichletObsolete.hpp.

◆ tRis

Range AnalyticalDirihletBC::tRis

Definition at line 368 of file AnalyticalDirichletObsolete.hpp.


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