v0.9.1
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, 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 613 of file MatrixManager.cpp.

614  : cOre(const_cast<MoFEM::Core &>(core)) {
615  PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
617  PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
619  PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
621  PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
623 }
PetscLogEvent MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays

◆ ~MatrixManager()

MoFEM::MatrixManager::~MatrixManager ( )

Destructor.

Definition at line 624 of file MatrixManager.cpp.

624 {}

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 741 of file MatrixManager.cpp.

742  {
743  MoFEM::CoreInterface &m_field = cOre;
745 
746  PetscLogEventBegin(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
747 
748  struct TestMatrixFillIn : public FEMethod {
749  CoreInterface *mFieldPtr;
750 
751  Mat A;
752 
753  int rowPrint, colPrint;
754 
755  TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
756  int col_print)
757  : mFieldPtr(m_field_ptr), A(a), rowPrint(row_print),
758  colPrint(col_print){};
759 
760  MoFEMErrorCode preProcess() {
763  }
764 
765  MoFEMErrorCode operator()() {
767 
768  if (refinedFiniteElementsPtr->find(
769  numeredEntFiniteElementPtr->getEnt()) ==
770  refinedFiniteElementsPtr->end()) {
771  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
772  "data inconsistency");
773  }
774 
775  for (FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
776  cit != colPtr->end(); cit++) {
777 
778  if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
779  refinedEntitiesPtr->end()) {
780  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
781  "data inconsistency");
782  }
783  if (!(*cit)->getActive()) {
784  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
785  "data inconsistency");
786  }
787 
788  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
789  Composite_Unique_mi_tag>::type::iterator ait;
790  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
791  boost::make_tuple((*cit)->getFieldEntityPtr()->getGlobalUniqueId(),
792  numeredEntFiniteElementPtr->getGlobalUniqueId()));
793  if (ait == adjacenciesPtr->end()) {
794  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
795  "adjacencies data inconsistency");
796  } else {
797  UId uid = ait->getEntUniqueId();
798  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
799  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
800  "data inconsistency");
801  }
802  if (dofsPtr->find((*cit)->getGlobalUniqueId()) == dofsPtr->end()) {
803  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
804  "data inconsistency");
805  }
806  }
807 
808  if ((*cit)->getEntType() != MBVERTEX) {
809 
810  FENumeredDofEntity_multiIndex::index<
811  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
812  hi_dit;
813  dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
814  .lower_bound(boost::make_tuple(
815  (*cit)->getName(), (*cit)->getEntType(),
816  (*cit)->sideNumberPtr->side_number));
817  hi_dit = colPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
818  .upper_bound(boost::make_tuple(
819  (*cit)->getName(), (*cit)->getEntType(),
820  (*cit)->sideNumberPtr->side_number));
821  int nb_dofs_on_ent = std::distance(dit, hi_dit);
822 
823  int max_order = (*cit)->getMaxOrder();
824  if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
825  nb_dofs_on_ent) {
826  std::cerr << "Warning: Number of Dofs in Col diffrent than number "
827  "of dofs for given entity order "
828  << (*cit)->getNbOfCoeffs() *
829  (*cit)->getOrderNbDofs(max_order)
830  << " " << nb_dofs_on_ent << std::endl;
831  }
832  }
833  }
834 
835  FENumeredDofEntity_multiIndex::iterator rit = rowPtr->begin();
836  for (; rit != rowPtr->end(); rit++) {
837 
838  if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
839  refinedEntitiesPtr->end()) {
840  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
841  "data inconsistency");
842  }
843  if (!(*rit)->getActive()) {
844  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
845  "data inconsistency");
846  }
847 
848  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
849  Composite_Unique_mi_tag>::type::iterator ait;
850  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
851  boost::make_tuple((*rit)->getFieldEntityPtr()->getGlobalUniqueId(),
852  numeredEntFiniteElementPtr->getGlobalUniqueId()));
853  if (ait == adjacenciesPtr->end()) {
854  std::ostringstream ss;
855  ss << *(*rit) << std::endl;
856  ss << *numeredEntFiniteElementPtr << std::endl;
857  ss << "dof: " << (*rit)->getBitRefLevel() << std::endl;
858  ss << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel()
859  << std::endl;
860  ss << "problem: " << problemPtr->getBitRefLevel() << std::endl;
861  ss << "problem mask: " << problemPtr->getMaskBitRefLevel()
862  << std::endl;
863  PetscPrintf(mFieldPtr->get_comm(), "%s", ss.str().c_str());
864  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
865  "adjacencies data inconsistency");
866  } else {
867  UId uid = ait->getEntUniqueId();
868  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
869  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
870  "data inconsistency");
871  }
872  if (dofsPtr->find((*rit)->getGlobalUniqueId()) == dofsPtr->end()) {
873  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
874  "data inconsistency");
875  }
876  }
877  int row = (*rit)->getPetscGlobalDofIdx();
878 
879  FENumeredDofEntity_multiIndex::iterator cit = colPtr->begin();
880  for (; cit != colPtr->end(); cit++) {
881 
882  int col = (*cit)->getPetscGlobalDofIdx();
883 
884  if (row == rowPrint && col == colPrint) {
885 
886  std::ostringstream ss;
887  ss << "fe:\n" << *numeredEntFiniteElementPtr << std::endl;
888  ss << "row:\n" << *(*rit) << std::endl;
889  ss << "col:\n" << *(*cit) << std::endl;
890 
891  ss << "fe:\n"
892  << numeredEntFiniteElementPtr->getBitRefLevel() << std::endl;
893  ss << "row:\n" << (*rit)->getBitRefLevel() << std::endl;
894  ss << "col:\n" << (*cit)->getBitRefLevel() << std::endl;
895 
896  std::cerr << ss.str() << std::endl;
897 
898  // PetscPrintf(mFieldPtr->get_comm(),"%s\n",ss.str().c_str());
899  }
900 
901  CHKERR MatSetValue(A, row, col, 1, INSERT_VALUES);
902  }
903 
904  if ((*rit)->getEntType() != MBVERTEX) {
905 
906  FENumeredDofEntity_multiIndex::index<
907  Composite_Name_Type_And_Side_Number_mi_tag>::type::iterator dit,
908  hi_dit;
909  dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
910  .lower_bound(boost::make_tuple(
911  (*rit)->getName(), (*rit)->getEntType(),
912  (*rit)->sideNumberPtr->side_number));
913  hi_dit = rowPtr->get<Composite_Name_Type_And_Side_Number_mi_tag>()
914  .upper_bound(boost::make_tuple(
915  (*rit)->getName(), (*rit)->getEntType(),
916  (*rit)->sideNumberPtr->side_number));
917  int nb_dofs_on_ent = std::distance(dit, hi_dit);
918 
919  int max_order = (*rit)->getMaxOrder();
920  if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
921  nb_dofs_on_ent) {
922  std::cerr << "Warning: Number of Dofs in Row diffrent than number "
923  "of dofs for given entity order "
924  << (*rit)->getNbOfCoeffs() *
925  (*rit)->getOrderNbDofs(max_order)
926  << " " << nb_dofs_on_ent << std::endl;
927  }
928  }
929  }
930 
932  }
933 
934  MoFEMErrorCode postProcess() {
936 
937  // cerr << mFieldPtr->get_comm_rank() << endl;
938  CHKERR MatAssemblyBegin(A, MAT_FLUSH_ASSEMBLY);
939  CHKERR MatAssemblyEnd(A, MAT_FLUSH_ASSEMBLY);
940 
942  }
943  };
944 
945  // create matrix
946  Mat A;
947  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, &A, verb);
948  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
949 
950  if (verb >= VERY_VERBOSE) {
951  MatView(A, PETSC_VIEWER_STDOUT_WORLD);
952  }
953 
954  if (verb >= NOISY) {
955  MatView(A, PETSC_VIEWER_DRAW_WORLD);
956  std::string wait;
957  std::cin >> wait;
958  }
959 
960  TestMatrixFillIn method(&m_field, A, row_print, col_print);
961 
962  // get problem
963  auto problems_ptr = m_field.get_problems();
964  auto &prb_set = problems_ptr->get<Problem_mi_tag>();
965  auto p_miit = prb_set.find(problem_name);
966  if (p_miit == prb_set.end())
967  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
968  "problem < %s > not found (top tip: check spelling)",
969  problem_name.c_str());
970 
971  if (verb >= VERBOSE)
972  PetscPrintf(m_field.get_comm(), "check problem < %s >\n",
973  problem_name.c_str());
974 
975  // loop all elements in problem and check if assemble is without error
976  auto fe_ptr = m_field.get_finite_elements();
977  for (auto &fe : *fe_ptr) {
978  if (verb >= VERBOSE)
979  PetscPrintf(m_field.get_comm(), "\tcheck element %s\n",
980  fe->getName().c_str());
981  CHKERR m_field.loop_finite_elements(problem_name, fe->getName(), method,
982  nullptr, MF_EXIST, verb);
983  }
984 
985  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
986  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
987  CHKERR MatDestroy(&A);
988 
989  PetscLogEventEnd(MOFEM_EVENT_checkMPIAIJWithArraysMatrixFillIn, 0, 0, 0, 0);
990 
992 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
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:514
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
#define CHKERR
Inline error check.
Definition: definitions.h:602
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
#define MatrixManagerFunctionBegin
Create adjacent matrices using different indices.
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 663 of file MatrixManager.cpp.

664  {
665  MoFEM::CoreInterface &m_field = cOre;
666  CreateRowComressedADJMatrix *core_ptr =
667  static_cast<CreateRowComressedADJMatrix *>(&cOre);
669  PetscLogEventBegin(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
670 
671  auto problems_ptr = m_field.get_problems();
672  auto &prb = problems_ptr->get<Problem_mi_tag>();
673  auto p_miit = prb.find(name);
674  if (p_miit == prb.end()) {
675  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
676  "problem < %s > is not found (top tip: check spelling)",
677  name.c_str());
678  }
679 
680  std::vector<int> i_vec, j_vec;
681  CHKERR core_ptr->createMatArrays<Idx_mi_tag>(p_miit, MATMPIADJ, i_vec, j_vec,
682  true, verb);
683  int *_i, *_j;
684  CHKERR PetscMalloc(i_vec.size() * sizeof(int), &_i);
685  CHKERR PetscMalloc(j_vec.size() * sizeof(int), &_j);
686  copy(i_vec.begin(), i_vec.end(), _i);
687  copy(j_vec.begin(), j_vec.end(), _j);
688 
689  int nb_col_dofs = p_miit->getNbDofsCol();
690  CHKERR MatCreateMPIAdj(m_field.get_comm(), i_vec.size() - 1, nb_col_dofs, _i,
691  _j, PETSC_NULL, Adj);
692  CHKERR MatSetOption(*Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
693 
694  PetscLogEventEnd(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
696 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
#define CHKERR
Inline error check.
Definition: definitions.h:602
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
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:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ 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 627 of file MatrixManager.cpp.

628  {
629  MoFEM::CoreInterface &m_field = cOre;
630  CreateRowComressedADJMatrix *core_ptr =
631  static_cast<CreateRowComressedADJMatrix *>(&cOre);
633  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
634 
635  auto problems_ptr = m_field.get_problems();
636  auto &prb = problems_ptr->get<Problem_mi_tag>();
637  auto p_miit = prb.find(name);
638  if (p_miit == prb.end()) {
639  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
640  "problem < %s > is not found (top tip: check spelling)",
641  name.c_str());
642  }
643 
644  std::vector<int> i_vec, j_vec;
645  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
646  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
647 
648  int nb_row_dofs = p_miit->getNbDofsRow();
649  int nb_col_dofs = p_miit->getNbDofsCol();
650  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
651  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
652 
653  CHKERR ::MatCreateMPIAIJWithArrays(
654  m_field.get_comm(), nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
655  nb_col_dofs, &*i_vec.begin(), &*j_vec.begin(), PETSC_NULL, Aij);
656 
657  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
659 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
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:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ 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 699 of file MatrixManager.cpp.

700  {
701  MoFEM::CoreInterface &m_field = cOre;
702  CreateRowComressedADJMatrix *core_ptr =
703  static_cast<CreateRowComressedADJMatrix *>(&cOre);
705  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
706 
707  auto problems_ptr = m_field.get_problems();
708  auto &prb = problems_ptr->get<Problem_mi_tag>();
709  auto p_miit = prb.find(name);
710  if (p_miit == prb.end()) {
711  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
712  "problem < %s > is not found (top tip: check spelling)",
713  name.c_str());
714  }
715 
716  std::vector<int> i_vec, j_vec;
717  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(p_miit, MATAIJ, i_vec,
718  j_vec, false, verb);
719 
720  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
721  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
722 
723  double *_a;
724  CHKERR PetscMalloc(j_vec.size() * sizeof(double), &_a);
725 
726  Mat tmpMat;
727  CHKERR ::MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, nb_local_dofs_row,
728  nb_local_dofs_col, &*i_vec.begin(),
729  &*j_vec.begin(), _a, &tmpMat);
730  CHKERR MatDuplicate(tmpMat, MAT_SHARE_NONZERO_PATTERN, Aij);
731  CHKERR MatDestroy(&tmpMat);
732 
733  CHKERR PetscFree(_a);
734 
735  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
737 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
InterfaceThis interface is used by user to:
Definition: Interface.hpp:42
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
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 601 of file MatrixManager.cpp.

602  {
604  *iface = NULL;
605  if (uuid == IDD_MOFEMMatrixManager) {
606  *iface = const_cast<MatrixManager *>(this);
608  }
609  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
611 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
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: