v0.9.0
Public Member Functions | Public Attributes | Private Attributes | List of all members
MoFEM::MatrixManager Struct Reference

Matrix manager is used to build and partition problems. More...

#include <src/interfaces/MatrixManager.hpp>

Inheritance diagram for MoFEM::MatrixManager:
[legend]
Collaboration diagram for MoFEM::MatrixManager:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 MatrixManager (const MoFEM::Core &core)
 
 ~MatrixManager ()
 Destructor. More...
 
template<class Tag >
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, Mat *Aij, int verb=QUIET)
 Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows. More...
 
template<class Tag >
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
 Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows. More...
 
template<class Tag >
MoFEMErrorCode createMPIAdjWithArrays (const std::string name, Mat *Adj, int verb=QUIET)
 Creates a sparse matrix representing an adjacency list. More...
 
template<class Tag >
MoFEMErrorCode createSeqAIJWithArrays (const std::string name, Mat *Aij, int verb=QUIET)
 Create sequencial matrix. More...
 
template<class Tag >
MoFEMErrorCode createSeqAIJWithArrays (const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
 Create sequencial matrix. More...
 
template<class Tag >
MoFEMErrorCode checkMPIAIJWithArraysMatrixFillIn (const std::string problem_name, int row_print, int col_print, int verb=QUIET)
 check if matrix fill in correspond to finite element indices More...
 
template<>
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAdjWithArrays (const std::string name, Mat *Adj, int verb)
 
template<>
MoFEMErrorCode createSeqAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode checkMPIAIJWithArraysMatrixFillIn (const std::string problem_name, int row_print, int col_print, int verb)
 
template<>
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAdjWithArrays (const std::string name, Mat *Adj, int verb)
 
template<>
MoFEMErrorCode createSeqAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode checkMPIAIJWithArraysMatrixFillIn (const std::string problem_name, int row_print, int col_print, int verb)
 
- 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 ()=default
 
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 Attributes

MoFEM::CorecOre
 

Private Attributes

PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
 
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
 
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
 
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
 

Additional Inherited Members

- 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

Matrix manager is used to build and partition problems.

Examples
build_problems.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, hcurl_check_approx_in_2d.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, and remove_entities_from_problem.cpp.

Definition at line 34 of file MatrixManager.hpp.

Constructor & Destructor Documentation

◆ MatrixManager()

MoFEM::MatrixManager::MatrixManager ( const MoFEM::Core core)

Definition at line 625 of file MatrixManager.cpp.

626  : cOre(const_cast<MoFEM::Core &>(core)) {
627  PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
629  PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
631  PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
633  PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
635 }
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays

◆ ~MatrixManager()

MoFEM::MatrixManager::~MatrixManager ( )

Destructor.

Definition at line 636 of file MatrixManager.cpp.

636 {}

Member Function Documentation

◆ checkMPIAIJWithArraysMatrixFillIn() [1/3]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::checkMPIAIJWithArraysMatrixFillIn ( const std::string  problem_name,
int  row_print,
int  col_print,
int  verb = QUIET 
)

check if matrix fill in correspond to finite element indices

This is used to check consistency of code. If problem is notices with additional non-zero elements in matrix, this function can help detect problem. Should be used as a part of atom tests

Template Parameters
Tag
Parameters
problem_name
rowprint info at particular row
colprint info at particular col
Returns
MoFEMErrorCode

Definition at line 200 of file MatrixManager.hpp.

202  {
203  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
204  return 0;
205  }

◆ checkMPIAIJWithArraysMatrixFillIn() [2/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::checkMPIAIJWithArraysMatrixFillIn ( const std::string  problem_name,
int  row_print,
int  col_print,
int  verb 
)

◆ checkMPIAIJWithArraysMatrixFillIn() [3/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::checkMPIAIJWithArraysMatrixFillIn ( const std::string  problem_name,
int  row_print,
int  col_print,
int  verb 
)

Definition at line 756 of file MatrixManager.cpp.

757  {
758  MoFEM::CoreInterface &m_field = cOre;
760  PetscLogEventBegin(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
761 
762  struct TestMatrixFillIn : public FEMethod {
763  CoreInterface *mFieldPtr;
764 
765  Mat A;
766 
767  int rowPrint, colPrint;
768 
769  TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
770  int col_print)
771  : mFieldPtr(m_field_ptr), A(a), rowPrint(row_print),
772  colPrint(col_print){};
773 
774  MoFEMErrorCode preProcess() {
777  }
778 
779  MoFEMErrorCode operator()() {
781 
782  if (refinedFiniteElementsPtr->find(
783  numeredEntFiniteElementPtr->getEnt()) ==
784  refinedFiniteElementsPtr->end()) {
785  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
786  "data inconsistency");
787  }
788 
789  for (FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
790  cit != colPtr->end(); cit++) {
791 
792  if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
793  refinedEntitiesPtr->end()) {
794  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
795  "data inconsistency");
796  }
797  if (!(*cit)->getActive()) {
798  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
799  "data inconsistency");
800  }
801 
802  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
803  Composite_Unique_mi_tag>::type::iterator ait;
804  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
805  boost::make_tuple((*cit)->getFieldEntityPtr()->getGlobalUniqueId(),
806  numeredEntFiniteElementPtr->getGlobalUniqueId()));
807  if (ait == adjacenciesPtr->end()) {
808  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
809  "adjacencies data inconsistency");
810  } else {
811  UId uid = ait->getEntUniqueId();
812  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
813  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
814  "data inconsistency");
815  }
816  if (dofsPtr->find((*cit)->getGlobalUniqueId()) == dofsPtr->end()) {
817  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
818  "data inconsistency");
819  }
820  }
821 
822  if ((*cit)->getEntType() != MBVERTEX) {
823 
824  FENumeredDofEntity_multiIndex::index<
825  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
826  hi_dit;
827  dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
828  .lower_bound(boost::make_tuple(
829  (*cit)->getName(), (*cit)->getEntType(),
830  (*cit)->sideNumberPtr->side_number));
831  hi_dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
832  .upper_bound(boost::make_tuple(
833  (*cit)->getName(), (*cit)->getEntType(),
834  (*cit)->sideNumberPtr->side_number));
835  int nb_dofs_on_ent = std::distance(dit, hi_dit);
836 
837  int max_order = (*cit)->getMaxOrder();
838  if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
839  nb_dofs_on_ent) {
840  std::cerr << "Warning: Number of Dofs in Col diffrent than number "
841  "of dofs for given entity order "
842  << (*cit)->getNbOfCoeffs() *
843  (*cit)->getOrderNbDofs(max_order)
844  << " " << nb_dofs_on_ent << std::endl;
845  }
846  }
847  }
848 
849  FENumeredDofEntity_multiIndex::iterator rit = rowPtr->begin();
850  for (; rit != rowPtr->end(); rit++) {
851 
852  if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
853  refinedEntitiesPtr->end()) {
854  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
855  "data inconsistency");
856  }
857  if (!(*rit)->getActive()) {
858  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
859  "data inconsistency");
860  }
861 
862  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
863  Composite_Unique_mi_tag>::type::iterator ait;
864  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
865  boost::make_tuple((*rit)->getFieldEntityPtr()->getGlobalUniqueId(),
866  numeredEntFiniteElementPtr->getGlobalUniqueId()));
867  if (ait == adjacenciesPtr->end()) {
868  std::ostringstream ss;
869  ss << *(*rit) << std::endl;
870  ss << *numeredEntFiniteElementPtr << std::endl;
871  ss << "dof: " << (*rit)->getBitRefLevel() << std::endl;
872  ss << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel()
873  << std::endl;
874  ss << "problem: " << problemPtr->getBitRefLevel() << std::endl;
875  ss << "problem mask: " << problemPtr->getMaskBitRefLevel()
876  << std::endl;
877  PetscPrintf(mFieldPtr->get_comm(), "%s", ss.str().c_str());
878  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
879  "adjacencies data inconsistency");
880  } else {
881  UId uid = ait->getEntUniqueId();
882  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
883  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
884  "data inconsistency");
885  }
886  if (dofsPtr->find((*rit)->getGlobalUniqueId()) == dofsPtr->end()) {
887  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
888  "data inconsistency");
889  }
890  }
891  int row = (*rit)->getPetscGlobalDofIdx();
892 
893  FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
894  for (; cit != colPtr->end(); cit++) {
895 
896  int col = (*cit)->getPetscGlobalDofIdx();
897 
898  if (row == rowPrint && col == colPrint) {
899 
900  std::ostringstream ss;
901  ss << "fe:\n" << *numeredEntFiniteElementPtr << std::endl;
902  ss << "row:\n" << *(*rit) << std::endl;
903  ss << "col:\n" << *(*cit) << std::endl;
904 
905  ss << "fe:\n"
906  << numeredEntFiniteElementPtr->getBitRefLevel() << std::endl;
907  ss << "row:\n" << (*rit)->getBitRefLevel() << std::endl;
908  ss << "col:\n" << (*cit)->getBitRefLevel() << std::endl;
909 
910  std::cerr << ss.str() << std::endl;
911 
912  // PetscPrintf(mFieldPtr->get_comm(),"%s\n",ss.str().c_str());
913  }
914 
915  CHKERR MatSetValue(A, row, col, 1, INSERT_VALUES);
916  }
917 
918  if ((*rit)->getEntType() != MBVERTEX) {
919 
920  FENumeredDofEntity_multiIndex::index<
921  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
922  hi_dit;
923  dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
924  .lower_bound(boost::make_tuple(
925  (*rit)->getName(), (*rit)->getEntType(),
926  (*rit)->sideNumberPtr->side_number));
927  hi_dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
928  .upper_bound(boost::make_tuple(
929  (*rit)->getName(), (*rit)->getEntType(),
930  (*rit)->sideNumberPtr->side_number));
931  int nb_dofs_on_ent = std::distance(dit, hi_dit);
932 
933  int max_order = (*rit)->getMaxOrder();
934  if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
935  nb_dofs_on_ent) {
936  std::cerr << "Warning: Number of Dofs in Row diffrent than number "
937  "of dofs for given entity order "
938  << (*rit)->getNbOfCoeffs() *
939  (*rit)->getOrderNbDofs(max_order)
940  << " " << nb_dofs_on_ent << std::endl;
941  }
942  }
943  }
944 
946  }
947 
948  MoFEMErrorCode postProcess() {
950 
951  // cerr << mFieldPtr->get_comm_rank() << endl;
952  CHKERR MatAssemblyBegin(A, MAT_FLUSH_ASSEMBLY);
953  CHKERR MatAssemblyEnd(A, MAT_FLUSH_ASSEMBLY);
954 
956  }
957  };
958 
959  // create matrix
960  Mat A;
961  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, &A, verb);
962  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
963 
964  if (verb >= VERY_VERBOSE) {
965  MatView(A, PETSC_VIEWER_STDOUT_WORLD);
966  }
967 
968  if (verb >= NOISY) {
969  MatView(A, PETSC_VIEWER_DRAW_WORLD);
970  std::string wait;
971  std::cin >> wait;
972  }
973 
974  TestMatrixFillIn method(&m_field, A, row_print, col_print);
975 
976  // get problem
977  const Problem_multiIndex *problems_ptr;
978  CHKERR m_field.get_problems(&problems_ptr);
979  auto &prb_set = problems_ptr->get<Problem_mi_tag>();
980  auto p_miit = prb_set.find(problem_name);
981  if (p_miit == prb_set.end())
982  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
983  "problem < %s > not found (top tip: check spelling)",
984  problem_name.c_str());
985 
986  if (verb >= VERBOSE)
987  PetscPrintf(m_field.get_comm(), "check problem < %s >\n",
988  problem_name.c_str());
989 
990  // loop all elements in problem and check if assemble is without error
991  const FiniteElement_multiIndex *fe_ptr;
992  CHKERR m_field.get_finite_elements(&fe_ptr);
993  for (auto &fe : *fe_ptr) {
994  if (verb >= VERBOSE)
995  PetscPrintf(m_field.get_comm(), "\tcheck element %s\n",
996  fe->getName().c_str());
997  CHKERR m_field.loop_finite_elements(problem_name, fe->getName(), method,
998  nullptr, MF_EXIST, verb);
999  }
1000 
1001  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
1002  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
1003  CHKERR MatDestroy(&A);
1004 
1005  PetscLogEventEnd(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
1006 
1008 }
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
virtual MoFEMErrorCode loop_finite_elements(const std::string &problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
virtual MoFEMErrorCode get_finite_elements(const FiniteElement_multiIndex **fe_ptr) const =0
Get finite elements multi-index.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
virtual MoFEMErrorCode get_problems(const Problem_multiIndex **problems_ptr) const =0
Get pointer to problems multi-index.
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
virtual MPI_Comm & get_comm() const =0
uint128_t UId
Unique Id.
Definition: Types.hpp:41

◆ createMPIAdjWithArrays() [1/3]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::createMPIAdjWithArrays ( const std::string  name,
Mat *  Adj,
int  verb = QUIET 
)

Creates a sparse matrix representing an adjacency list.

The matrix does not have numerical values associated with it, but is intended for ordering (to reduce bandwidth etc) and partitioning.

See PETSc for details

Note
This matrix object does not support most matrix operations, include MatSetValues(). You must NOT free the ii, values and jj arrays yourself. PETSc will free them when the matrix is destroyed; you must allocate them with PetscMalloc(). If you call from Fortran you need not create the arrays with PetscMalloc(). Should not include the matrix diagonals.
Template Parameters
Tag
Parameters
name
Adj
verb
Returns
MoFEMErrorCode

Definition at line 117 of file MatrixManager.hpp.

118  {
119  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
120  return MOFEM_NOT_IMPLEMENTED;
121  }

◆ createMPIAdjWithArrays() [2/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createMPIAdjWithArrays ( const std::string  name,
Mat *  Adj,
int  verb 
)

◆ createMPIAdjWithArrays() [3/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createMPIAdjWithArrays ( const std::string  name,
Mat *  Adj,
int  verb 
)

Definition at line 676 of file MatrixManager.cpp.

677  {
678  MoFEM::CoreInterface &m_field = cOre;
679  CreateRowComressedADJMatrix *core_ptr =
680  static_cast<CreateRowComressedADJMatrix *>(&cOre);
682  PetscLogEventBegin(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
683 
684  const Problem_multiIndex *problems_ptr;
685  CHKERR m_field.get_problems(&problems_ptr);
686  auto &prb = problems_ptr->get<Problem_mi_tag>();
687  auto p_miit = prb.find(name);
688  if (p_miit == prb.end()) {
689  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
690  "problem < %s > is not found (top tip: check spelling)",
691  name.c_str());
692  }
693 
694  std::vector<int> i_vec, j_vec;
695  CHKERR core_ptr->createMatArrays<Idx_mi_tag>(p_miit, MATMPIADJ, i_vec, j_vec,
696  true, verb);
697  int *_i, *_j;
698  CHKERR PetscMalloc(i_vec.size() * sizeof(int), &_i);
699  CHKERR PetscMalloc(j_vec.size() * sizeof(int), &_j);
700  copy(i_vec.begin(), i_vec.end(), _i);
701  copy(j_vec.begin(), j_vec.end(), _j);
702 
703  int nb_col_dofs = p_miit->getNbDofsCol();
704  CHKERR MatCreateMPIAdj(m_field.get_comm(), i_vec.size() - 1, nb_col_dofs, _i,
705  _j, PETSC_NULL, Adj);
706  CHKERR MatSetOption(*Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
707 
708  PetscLogEventEnd(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
710 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
virtual MoFEMErrorCode get_problems(const Problem_multiIndex **problems_ptr) const =0
Get pointer to problems multi-index.
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
virtual MPI_Comm & get_comm() const =0

◆ createMPIAIJWithArrays() [1/4]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::createMPIAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb = QUIET 
)

Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij
verb
Returns
MoFEMErrorCode

Definition at line 62 of file MatrixManager.hpp.

63  {
64  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
65  return 0;
66  }

◆ createMPIAIJWithArrays() [2/4]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::createMPIAIJWithArrays ( const std::string  name,
SmartPetscObj< Mat > &  aij_ptr,
int  verb = QUIET 
)

Creates a MPI AIJ matrix using arrays that contain in standard CSR format the local rows.

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij(SmartPetscObj)
verb
Returns
MoFEMErrorCode

Definition at line 83 of file MatrixManager.hpp.

85  {
87  Mat aij;
88  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
89  aij_ptr.reset(aij, false);
91  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ createMPIAIJWithArrays() [3/4]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createMPIAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb 
)

◆ createMPIAIJWithArrays() [4/4]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createMPIAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb 
)

Definition at line 639 of file MatrixManager.cpp.

640  {
641  MoFEM::CoreInterface &m_field = cOre;
642  CreateRowComressedADJMatrix *core_ptr =
643  static_cast<CreateRowComressedADJMatrix *>(&cOre);
645  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
646 
647  const Problem_multiIndex *problems_ptr;
648  CHKERR m_field.get_problems(&problems_ptr);
649  auto &prb = problems_ptr->get<Problem_mi_tag>();
650  auto p_miit = prb.find(name);
651  if (p_miit == prb.end()) {
652  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
653  "problem < %s > is not found (top tip: check spelling)",
654  name.c_str());
655  }
656 
657  std::vector<int> i_vec, j_vec;
658  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
659  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
660 
661  int nb_row_dofs = p_miit->getNbDofsRow();
662  int nb_col_dofs = p_miit->getNbDofsCol();
663  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
664  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
665 
666  CHKERR ::MatCreateMPIAIJWithArrays(
667  m_field.get_comm(), nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
668  nb_col_dofs, &*i_vec.begin(), &*j_vec.begin(), PETSC_NULL, Aij);
669 
670  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
672 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
virtual MoFEMErrorCode get_problems(const Problem_multiIndex **problems_ptr) const =0
Get pointer to problems multi-index.
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
virtual MPI_Comm & get_comm() const =0

◆ createSeqAIJWithArrays() [1/4]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::createSeqAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb = QUIET 
)

Create sequencial matrix.

Creates a sparse matrix in AIJ (compressed row) format (the default parallel PETSc format). For good matrix assembly performance the user should preallocate the matrix storage by setting the parameter nz (or the array nnz). By setting these parameters accurately, performance during matrix assembly can be increased by more than a factor of 50.

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij
i
j
v
verb
Returns
MoFEMErrorCode

Definition at line 146 of file MatrixManager.hpp.

147  {
148  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
149  return 0;
150  }

◆ createSeqAIJWithArrays() [2/4]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::createSeqAIJWithArrays ( const std::string  name,
SmartPetscObj< Mat > &  aij_ptr,
int  verb = QUIET 
)

Create sequencial matrix.

Creates a sparse matrix in AIJ (compressed row) format (the default parallel PETSc format). For good matrix assembly performance the user should preallocate the matrix storage by setting the parameter nz (or the array nnz). By setting these parameters accurately, performance during matrix assembly can be increased by more than a factor of 50.

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij(SmartPetscObj)
i
j
v
verb
Returns
MoFEMErrorCode

Definition at line 175 of file MatrixManager.hpp.

177  {
179  Mat aij;
180  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
181  aij_ptr.reset(aij, false);
183  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ createSeqAIJWithArrays() [3/4]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createSeqAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb 
)

◆ createSeqAIJWithArrays() [4/4]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createSeqAIJWithArrays ( const std::string  name,
Mat *  Aij,
int  verb 
)

Definition at line 713 of file MatrixManager.cpp.

714  {
715  MoFEM::CoreInterface &m_field = cOre;
716  CreateRowComressedADJMatrix *core_ptr =
717  static_cast<CreateRowComressedADJMatrix *>(&cOre);
719  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
720 
721  const Problem_multiIndex *problems_ptr;
722  CHKERR m_field.get_problems(&problems_ptr);
723  auto &prb = problems_ptr->get<Problem_mi_tag>();
724  auto p_miit = prb.find(name);
725  if (p_miit == prb.end()) {
726  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
727  "problem < %s > is not found (top tip: check spelling)",
728  name.c_str());
729  }
730 
731  std::vector<int> i_vec, j_vec;
732  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(p_miit, MATAIJ, i_vec,
733  j_vec, false, verb);
734 
735  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
736  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
737 
738  double *_a;
739  CHKERR PetscMalloc(j_vec.size() * sizeof(double), &_a);
740 
741  Mat tmpMat;
742  CHKERR ::MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, nb_local_dofs_row,
743  nb_local_dofs_col, &*i_vec.begin(),
744  &*j_vec.begin(), _a, &tmpMat);
745  CHKERR MatDuplicate(tmpMat, MAT_SHARE_NONZERO_PATTERN, Aij);
746  CHKERR MatDestroy(&tmpMat);
747 
748  CHKERR PetscFree(_a);
749 
750  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
752 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
virtual MoFEMErrorCode get_problems(const Problem_multiIndex **problems_ptr) const =0
Get pointer to problems multi-index.
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
virtual MPI_Comm & get_comm() const =0

◆ query_interface()

MoFEMErrorCode MoFEM::MatrixManager::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 613 of file MatrixManager.cpp.

614  {
616  *iface = NULL;
617  if (uuid == IDD_MOFEMMatrixManager) {
618  *iface = const_cast<MatrixManager *>(this);
620  }
621  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
623 }
#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
static const MOFEMuuid IDD_MOFEMMatrixManager

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::MatrixManager::cOre

Definition at line 39 of file MatrixManager.hpp.

◆ MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
private

Definition at line 211 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAdjWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
private

Definition at line 209 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
private

Definition at line 208 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
private

Definition at line 210 of file MatrixManager.hpp.


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