v0.14.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 (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 MatrixManager (const MoFEM::Core &core)
 
virtual ~MatrixManager ()=default
 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 createMPIAIJCUSPARSEWithArrays (const std::string name, Mat *Aij, int verb=QUIET)
 
template<class Tag >
MoFEMErrorCode createMPIAIJCUSPARSEWithArrays (const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
 
template<class Tag >
MoFEMErrorCode createSeqAIJCUSPARSEWithArrays (const std::string name, Mat *Aij, int verb=QUIET)
 
template<class Tag >
MoFEMErrorCode createSeqAIJCUSPARSEWithArrays (const std::string name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
 
template<class Tag >
MoFEMErrorCode createMPIAIJ (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 createMPIAIJ (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...
 
MoFEMErrorCode checkMatrixFillIn (const std::string problem_name, int row_print, int col_print, Mat A, int verb=QUIET)
 check if matrix fill in correspond to finite element indices 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<class Tag >
MoFEMErrorCode checkMPIAIJMatrixFillIn (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<class Tag >
MoFEMErrorCode createHybridL2MPIAIJ (const std::string problem_name, SmartPetscObj< Mat > &aij_ptr, int verb=QUIET)
 Create a Hybrid MPIAij object. More...
 
template<>
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAIJCUSPARSEWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createSeqAIJCUSPARSEWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAIJ (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 checkMPIAIJMatrixFillIn (const std::string problem_name, int row_print, int col_print, int verb)
 
template<>
MoFEMErrorCode createHybridL2MPIAIJ (const std::string problem_name, SmartPetscObj< Mat > &aij_ptr, int verb)
 
template<>
MoFEMErrorCode createMPIAIJWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAIJCUSPARSEWithArrays (const std::string name, Mat *Aij, int verb)
 
template<>
MoFEMErrorCode createMPIAIJ (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 createSeqAIJCUSPARSEWithArrays (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 checkMPIAIJMatrixFillIn (const std::string problem_name, int row_print, int col_print, int verb)
 
template<>
MoFEMErrorCode createHybridL2MPIAIJ (const std::string problem_name, SmartPetscObj< Mat > &aij_ptr, int verb)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference 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
 

Public Attributes

MoFEM::CorecOre
 

Private Attributes

PetscLogEvent MOFEM_EVENT_createMPIAIJ
 
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
 
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
 
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
 
PetscLogEvent MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
 
PetscLogEvent MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
 
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Detailed Description

Matrix manager is used to build and partition problems.

Examples
build_large_problem.cpp, build_problems.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, nonlinear_dynamics.cpp, prism_polynomial_approximation.cpp, quad_polynomial_approximation.cpp, remove_entities_from_problem.cpp, and remove_entities_from_problem_not_partitioned.cpp.

Definition at line 21 of file MatrixManager.hpp.

Constructor & Destructor Documentation

◆ MatrixManager()

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

Definition at line 606 of file MatrixManager.cpp.

607  : cOre(const_cast<MoFEM::Core &>(core)) {
608  PetscLogEventRegister("MatMngCrtMPIAIJ", 0, &MOFEM_EVENT_createMPIAIJ);
609  PetscLogEventRegister("MatMngCrtMPIAIJWthArr", 0,
611  PetscLogEventRegister("MatMngCrtMPIAdjWithArr", 0,
613  PetscLogEventRegister("MatMngCrtMPIAIJCUSPARSEWthArr", 0,
615  PetscLogEventRegister("MatMngCrtSeqAIJCUSPARSEWthArrs", 0,
617  PetscLogEventRegister("MatMngCrtSeqAIJWthArrs", 0,
619  PetscLogEventRegister("MatMngCrtCheckMPIAIJWthArrFillIn", 0,
621 }

◆ ~MatrixManager()

virtual MoFEM::MatrixManager::~MatrixManager ( )
virtualdefault

Destructor.

Member Function Documentation

◆ checkMatrixFillIn()

MoFEMErrorCode MoFEM::MatrixManager::checkMatrixFillIn ( const std::string  problem_name,
int  row_print,
int  col_print,
Mat  A,
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

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

Definition at line 899 of file MatrixManager.cpp.

901  {
902  MoFEM::CoreInterface &m_field = cOre;
904 
905  PetscLogEventBegin(MOFEM_EVENT_checkMatrixFillIn, 0, 0, 0, 0);
906 
907  struct TestMatrixFillIn : public FEMethod {
908  CoreInterface *mFieldPtr;
909 
910  Mat A;
911 
912  int rowPrint, colPrint;
913 
914  TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
915  int col_print)
916  : mFieldPtr(m_field_ptr), A(a), rowPrint(row_print),
917  colPrint(col_print){};
918 
919  MoFEMErrorCode preProcess() {
922  }
923 
924  MoFEMErrorCode operator()() {
926 
927  if (refinedFiniteElementsPtr->find(
928  numeredEntFiniteElementPtr->getEnt()) ==
929  refinedFiniteElementsPtr->end()) {
930  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
931  "data inconsistency");
932  }
933 
934  auto row_dofs = getRowDofsPtr();
935  auto col_dofs = getColDofsPtr();
936 
937  for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
938 
939  if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
940  refinedEntitiesPtr->end()) {
941  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
942  "data inconsistency");
943  }
944  if (!(*cit)->getActive()) {
945  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
946  "data inconsistency");
947  }
948 
949  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
950  Composite_Unique_mi_tag>::type::iterator ait;
951  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
952  boost::make_tuple((*cit)->getFieldEntityPtr()->getLocalUniqueId(),
953  numeredEntFiniteElementPtr->getLocalUniqueId()));
954  if (ait == adjacenciesPtr->end()) {
955  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
956  "adjacencies data inconsistency");
957  } else {
958  UId uid = ait->getEntUniqueId();
959  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
960  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
961  "data inconsistency");
962  }
963  if (dofsPtr->find((*cit)->getLocalUniqueId()) == dofsPtr->end()) {
964  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
965  "data inconsistency");
966  }
967  }
968 
969  if ((*cit)->getEntType() != MBVERTEX) {
970 
971  auto range =
972  col_dofs->get<Ent_mi_tag>().equal_range((*cit)->getEnt());
973  int nb_dofs_on_ent = std::distance(range.first, range.second);
974 
975  int max_order = (*cit)->getMaxOrder();
976  if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
977  nb_dofs_on_ent) {
978 
979  /* It could be that you have
980  removed DOFs from problem, and for example if this was vector filed
981  with components {Ux,Uy,Uz}, you removed on Uz element.*/
982 
983  MOFEM_LOG("SELF", Sev::warning)
984  << "Warning: Number of Dofs in Col different than number "
985  "of dofs for given entity order "
986  << (*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order)
987  << " " << nb_dofs_on_ent;
988  }
989  }
990  }
991 
992  for (auto rit = row_dofs->begin(); rit != row_dofs->end(); rit++) {
993 
994  if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
995  refinedEntitiesPtr->end()) {
996  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
997  "data inconsistency");
998  }
999  if (!(*rit)->getActive()) {
1000  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
1001  "data inconsistency");
1002  }
1003 
1004  FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
1005  Composite_Unique_mi_tag>::type::iterator ait;
1006  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
1007  boost::make_tuple((*rit)->getFieldEntityPtr()->getLocalUniqueId(),
1008  numeredEntFiniteElementPtr->getLocalUniqueId()));
1009  if (ait == adjacenciesPtr->end()) {
1011  MOFEM_LOG("SELF", Sev::error) << *(*rit);
1012  MOFEM_LOG("SELF", Sev::error) << *(*rit);
1013  MOFEM_LOG("SELF", Sev::error) << *numeredEntFiniteElementPtr;
1014  MOFEM_LOG("SELF", Sev::error) << "dof: " << (*rit)->getBitRefLevel();
1015  MOFEM_LOG("SELF", Sev::error)
1016  << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel();
1017  MOFEM_LOG("SELF", Sev::error)
1018  << "problem: " << problemPtr->getBitRefLevel();
1019  MOFEM_LOG("SELF", Sev::error)
1020  << "problem mask: " << problemPtr->getBitRefLevelMask();
1021  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
1022  "adjacencies data inconsistency");
1023  } else {
1024  UId uid = ait->getEntUniqueId();
1025  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
1026  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
1027  "data inconsistency");
1028  }
1029  if (dofsPtr->find((*rit)->getLocalUniqueId()) == dofsPtr->end()) {
1030  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
1031  "data inconsistency");
1032  }
1033  }
1034  int row = (*rit)->getPetscGlobalDofIdx();
1035 
1036  auto col_dofs = getColDofsPtr();
1037  for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
1038 
1039  int col = (*cit)->getPetscGlobalDofIdx();
1040 
1041  if (row == rowPrint && col == colPrint) {
1042  MOFEM_LOG("SELF", Sev::noisy) << "fe:\n"
1043  << *numeredEntFiniteElementPtr;
1044  MOFEM_LOG("SELF", Sev::noisy) << "row:\n" << *(*rit);
1045  MOFEM_LOG("SELF", Sev::noisy) << "col:\n" << *(*cit);
1046  MOFEM_LOG("SELF", Sev::noisy)
1047  << "fe:\n"
1048  << numeredEntFiniteElementPtr->getBitRefLevel();
1049  MOFEM_LOG("SELF", Sev::noisy) << "row:\n"
1050  << (*rit)->getBitRefLevel();
1051  MOFEM_LOG("SELF", Sev::noisy) << "col:\n"
1052  << (*cit)->getBitRefLevel();
1053  }
1054 
1055  CHKERR MatSetValue(A, row, col, 1, INSERT_VALUES);
1056  }
1057 
1058  if ((*rit)->getEntType() != MBVERTEX) {
1059 
1060  auto range =
1061  row_dofs->get<Ent_mi_tag>().equal_range((*rit)->getEnt());
1062  int nb_dofs_on_ent = std::distance(range.first, range.second);
1063 
1064  int max_order = (*rit)->getMaxOrder();
1065  if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
1066  nb_dofs_on_ent) {
1067 
1068  /* It could be that you have removed DOFs from problem, and for
1069  * example if this was vector filed with components {Ux,Uy,Uz}, you
1070  * removed on Uz element. */
1071 
1072  MOFEM_LOG("SELF", Sev::warning)
1073  << "Warning: Number of Dofs in Row different than number "
1074  "of dofs for given entity order "
1075  << (*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order)
1076  << " " << nb_dofs_on_ent;
1077  }
1078  }
1079  }
1080 
1082  }
1083 
1084  MoFEMErrorCode postProcess() {
1086 
1087  // cerr << mFieldPtr->get_comm_rank() << endl;
1088  CHKERR MatAssemblyBegin(A, MAT_FLUSH_ASSEMBLY);
1089  CHKERR MatAssemblyEnd(A, MAT_FLUSH_ASSEMBLY);
1090 
1092  }
1093  };
1094 
1095  // create matrix
1096  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1097 
1098  if (verb >= VERY_VERBOSE) {
1099  MatView(A, PETSC_VIEWER_STDOUT_WORLD);
1100  }
1101 
1102  if (verb >= NOISY) {
1103  MatView(A, PETSC_VIEWER_DRAW_WORLD);
1104  std::string wait;
1105  std::cin >> wait;
1106  }
1107 
1108  TestMatrixFillIn method(&m_field, A, row_print, col_print);
1109 
1110  // get problem
1111  auto problems_ptr = m_field.get_problems();
1112  auto &prb_set = problems_ptr->get<Problem_mi_tag>();
1113  auto p_miit = prb_set.find(problem_name);
1114  if (p_miit == prb_set.end())
1115  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
1116  "problem < %s > not found (top tip: check spelling)",
1117  problem_name.c_str());
1118  MOFEM_LOG_C("WORLD", Sev::inform, "check problem < %s >",
1119  problem_name.c_str());
1120 
1121  // loop all elements in problem and check if assemble is without error
1122  auto fe_ptr = m_field.get_finite_elements();
1123  for (auto &fe : *fe_ptr) {
1124  MOFEM_LOG_C("WORLD", Sev::verbose, "\tcheck element %s",
1125  fe->getName().c_str());
1126  CHKERR m_field.loop_finite_elements(problem_name, fe->getName(), method,
1127  nullptr, MF_EXIST,
1128  CacheTupleSharedPtr(), verb);
1129  }
1130 
1131  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
1132  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
1133 
1134  PetscLogEventEnd(MOFEM_EVENT_checkMatrixFillIn, 0, 0, 0, 0);
1135 
1137 }

◆ checkMPIAIJMatrixFillIn() [1/3]

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

◆ checkMPIAIJMatrixFillIn() [2/3]

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

Definition at line 1153 of file MatrixManager.cpp.

1154  {
1156  // create matrix
1157  SmartPetscObj<Mat> A;
1158  CHKERR createMPIAIJ<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1159  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1160  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1162 }

◆ checkMPIAIJMatrixFillIn() [3/3]

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

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 286 of file MatrixManager.hpp.

288  {
289  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
290  return 0;
291  }

◆ checkMPIAIJWithArraysMatrixFillIn() [1/3]

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

◆ checkMPIAIJWithArraysMatrixFillIn() [2/3]

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

Definition at line 1141 of file MatrixManager.cpp.

1142  {
1144  // create matrix
1145  SmartPetscObj<Mat> A;
1146  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1147  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1148  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1150 }

◆ checkMPIAIJWithArraysMatrixFillIn() [3/3]

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

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 265 of file MatrixManager.hpp.

267  {
268  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
269  return 0;
270  }

◆ createHybridL2MPIAIJ() [1/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createHybridL2MPIAIJ ( const std::string  problem_name,
SmartPetscObj< Mat > &  aij_ptr,
int  verb 
)

◆ createHybridL2MPIAIJ() [2/3]

template<>
MoFEMErrorCode MoFEM::MatrixManager::createHybridL2MPIAIJ ( const std::string  problem_name,
SmartPetscObj< Mat > &  aij_ptr,
int  verb 
)

Definition at line 1165 of file MatrixManager.cpp.

1166  {
1167  MoFEM::CoreInterface &m_field = cOre;
1169 
1170  auto prb_ptr = m_field.get_problems();
1171  auto p_miit = prb_ptr->get<Problem_mi_tag>().find(problem_name);
1172  if (p_miit == prb_ptr->get<Problem_mi_tag>().end()) {
1173  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
1174  "problem < %s > is not found (top tip: check spelling)",
1175  problem_name.c_str());
1176  }
1177 
1178  auto nb_rows = p_miit->getNbDofsRow();
1179  auto nb_cols = p_miit->getNbDofsCol();
1180  auto nb_loc_rows = p_miit->getNbLocalDofsRow();
1181  auto nb_loc_cols = p_miit->getNbLocalDofsCol();
1182 
1183  auto row_ptr = p_miit->getNumeredRowDofsPtr();
1184  auto col_ptr = p_miit->getNumeredColDofsPtr();
1185 
1186  BitFieldId fields_ids;
1187  for (auto &c : *col_ptr) {
1188  fields_ids |= c->getId();
1189  }
1190  auto fields = m_field.get_fields();
1191  std::vector<int> fields_bit_numbers;
1192  for (auto &f : *fields) {
1193  if ((fields_ids & f->getId()).any()) {
1194  fields_bit_numbers.push_back(f->getBitNumber());
1195  }
1196  }
1197 
1198  Range adj;
1199  EntityHandle prev_ent = 0;
1200  int prev_dim = -1;
1201 
1202  auto get_adj = [&](auto ent, auto dim) {
1203  if (prev_ent == ent && prev_dim == dim) {
1204  return adj;
1205  } else {
1206  adj.clear();
1207  Range bridge;
1208  CHKERR m_field.get_moab().get_adjacencies(&ent, 1, dim + 1, false,
1209  bridge);
1210  CHKERR m_field.get_moab().get_adjacencies(bridge, dim, false, adj,
1211  moab::Interface::UNION);
1212  prev_ent = ent;
1213  prev_dim = dim;
1214  return adj;
1215  }
1216  };
1217 
1218  int prev_bit_number = -1;
1219  EntityHandle prev_first_ent = 0;
1220  EntityHandle prev_second_ent = 0;
1221  using IT = boost::multi_index::index<NumeredDofEntity_multiIndex,
1222  Unique_mi_tag>::type::iterator;
1223  std::pair<IT, IT> pair_lo_hi;
1224  auto get_col_it = [&](auto bit_number, auto first_ent, auto second_ent) {
1225  if (bit_number == prev_bit_number && first_ent == prev_first_ent &&
1226  second_ent == prev_second_ent) {
1227  return pair_lo_hi;
1228  } else {
1229  auto lo_it = col_ptr->get<Unique_mi_tag>().lower_bound(
1230  DofEntity::getLoFieldEntityUId(bit_number, first_ent));
1231  auto hi_it = col_ptr->get<Unique_mi_tag>().upper_bound(
1232  DofEntity::getHiFieldEntityUId(bit_number, second_ent));
1233  pair_lo_hi = std::make_pair(lo_it, hi_it);
1234  prev_bit_number = bit_number;
1235  prev_first_ent = first_ent;
1236  prev_second_ent = second_ent;
1237  return pair_lo_hi;
1238  }
1239  };
1240 
1241  auto create_ghost_vec = [&]() {
1242  SmartPetscObj<Vec> v;
1243  CHKERR m_field.getInterface<VecManager>()->vecCreateGhost(problem_name, ROW,
1244  v);
1245  CHKERR VecZeroEntries(v);
1246  CHKERR VecGhostUpdateBegin(v, INSERT_VALUES, SCATTER_FORWARD);
1247  CHKERR VecGhostUpdateEnd(v, INSERT_VALUES, SCATTER_FORWARD);
1248  return v;
1249  };
1250 
1251  auto v_o_nnz = create_ghost_vec();
1252  auto v_d_nnz = create_ghost_vec();
1253 
1254 
1255  double *o_nnz_real;
1256  CHKERR VecGetArray(v_o_nnz, &o_nnz_real);
1257  double *d_nnz_real;
1258  CHKERR VecGetArray(v_d_nnz, &d_nnz_real);
1259 
1260  for (auto r = row_ptr->begin(); r != row_ptr->end(); ++r) {
1261 
1262  auto row_loc_idx = (*r)->getPetscLocalDofIdx();
1263  if (row_loc_idx < 0)
1264  continue;
1265 
1266  auto ent = (*r)->getEnt();
1267  auto dim = dimension_from_handle(ent);
1268  auto adj = get_adj(ent, dim);
1269 
1270  for (auto p = adj.pair_begin(); p != adj.pair_end(); ++p) {
1271  auto first_ent = p->first;
1272  auto second_ent = p->second;
1273 
1274  for (auto bit_number : fields_bit_numbers) {
1275 
1276  auto [lo_it, hi_it] = get_col_it(bit_number, first_ent, second_ent);
1277 
1278  for (; lo_it != hi_it; ++lo_it) {
1279  auto col_loc_idx = (*lo_it)->getPetscLocalDofIdx();
1280  if (col_loc_idx < 0)
1281  continue;
1282 
1283  if (
1284 
1285  (*lo_it)->getOwnerProc() == (*r)->getOwnerProc()
1286 
1287  ) {
1288  d_nnz_real[row_loc_idx] += 1;
1289  } else {
1290  o_nnz_real[row_loc_idx] += 1;
1291  }
1292  }
1293  }
1294  }
1295  }
1296 
1297  CHKERR VecRestoreArray(v_o_nnz, &o_nnz_real);
1298  CHKERR VecRestoreArray(v_d_nnz, &d_nnz_real);
1299 
1300  auto update_vec = [&](auto v) {
1302  CHKERR VecGhostUpdateBegin(v, ADD_VALUES, SCATTER_REVERSE);
1303  CHKERR VecGhostUpdateEnd(v, ADD_VALUES, SCATTER_REVERSE);
1304  CHKERR VecAssemblyBegin(v);
1305  CHKERR VecAssemblyEnd(v);
1307  };
1308 
1309  CHKERR update_vec(v_o_nnz);
1310  CHKERR update_vec(v_d_nnz);
1311 
1312  int o_nz = 0;
1313  int d_nz = 0;
1314 
1315  CHKERR VecGetArray(v_d_nnz, &d_nnz_real);
1316  CHKERR VecGetArray(v_o_nnz, &o_nnz_real);
1317 
1318  std::vector<int> d_nnz(nb_loc_rows, 0);
1319  std::vector<int> o_nnz(nb_loc_rows, 0);
1320  for (auto r = row_ptr->begin(); r != row_ptr->end(); ++r) {
1321 
1322  auto row_loc_idx = (*r)->getPetscLocalDofIdx();
1323 
1324  if (
1325 
1326  row_loc_idx >= 0 && row_loc_idx < nb_loc_rows
1327 
1328  ) {
1329  auto row_loc_idx = (*r)->getPetscLocalDofIdx();
1330  d_nz += o_nnz_real[row_loc_idx];
1331  d_nnz[row_loc_idx] = d_nnz_real[row_loc_idx];
1332  o_nz += o_nnz_real[row_loc_idx];
1333  o_nnz[row_loc_idx] = o_nnz_real[row_loc_idx];
1334  }
1335 
1336  }
1337 
1338  CHKERR VecRestoreArray(v_o_nnz, &o_nnz_real);
1339  CHKERR VecRestoreArray(v_d_nnz, &d_nnz_real);
1340 
1341  if (verb >= QUIET) {
1342  MOFEM_LOG("SYNC", Sev::verbose)
1343  << "Hybrid L2 matrix d_nz: " << d_nz << " o_nz: " << o_nz;
1344  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
1345  }
1346 
1347  if (verb >= VERBOSE) {
1348  MOFEM_LOG("SYNC", Sev::noisy) << "Hybrid L2 matrix";
1349  int idx = 0;
1350  for (auto &d : d_nnz) {
1351  MOFEM_LOG("SYNC", Sev::noisy) << idx << ": " << d;
1352  ++idx;
1353  }
1354  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::noisy);
1355  }
1356 
1357  Mat a_raw;
1358  CHKERR MatCreateAIJ(m_field.get_comm(), nb_loc_rows, nb_loc_cols, nb_rows,
1359  nb_cols, 0, &*d_nnz.begin(), 0, &*o_nnz.begin(), &a_raw);
1360  CHKERR MatSetOption(a_raw, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1361  aij_ptr = SmartPetscObj<Mat>(a_raw);
1362 
1363  MOFEM_LOG_CHANNEL("WORLD");
1364  MOFEM_LOG_CHANNEL("SYNC");
1365 
1367 }

◆ createHybridL2MPIAIJ() [3/3]

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

Create a Hybrid MPIAij object.

It assumes that L2 field is on skeleton, and dim+1 is bridge adjacency.

Template Parameters
Tag
Parameters
problem_name
aij_ptr
verb
Returns
MoFEMErrorCode

Definition at line 305 of file MatrixManager.hpp.

307  {
308  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
309  return 0;
310  }

◆ createMPIAdjWithArrays() [1/3]

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

◆ createMPIAdjWithArrays() [2/3]

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

Definition at line 821 of file MatrixManager.cpp.

822  {
823  MoFEM::CoreInterface &m_field = cOre;
824  CreateRowCompressedADJMatrix *core_ptr =
825  static_cast<CreateRowCompressedADJMatrix *>(&cOre);
827  PetscLogEventBegin(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
828 
829  auto problems_ptr = m_field.get_problems();
830  auto &prb = problems_ptr->get<Problem_mi_tag>();
831  auto p_miit = prb.find(name);
832  if (p_miit == prb.end()) {
833  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
834  "problem < %s > is not found (top tip: check spelling)",
835  name.c_str());
836  }
837 
838  std::vector<int> i_vec, j_vec;
839  j_vec.reserve(10000);
840  CHKERR core_ptr->createMatArrays<Idx_mi_tag>(p_miit, MATMPIADJ, i_vec, j_vec,
841  true, verb);
842  int *_i, *_j;
843  CHKERR PetscMalloc(i_vec.size() * sizeof(int), &_i);
844  CHKERR PetscMalloc(j_vec.size() * sizeof(int), &_j);
845  copy(i_vec.begin(), i_vec.end(), _i);
846  copy(j_vec.begin(), j_vec.end(), _j);
847 
848  int nb_col_dofs = p_miit->getNbDofsCol();
849  CHKERR MatCreateMPIAdj(m_field.get_comm(), i_vec.size() - 1, nb_col_dofs, _i,
850  _j, PETSC_NULL, Adj);
851  CHKERR MatSetOption(*Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
852 
853  PetscLogEventEnd(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
855 }

◆ createMPIAdjWithArrays() [3/3]

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

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 166 of file MatrixManager.hpp.

167  {
168  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
169  return MOFEM_NOT_IMPLEMENTED;
170  }

◆ createMPIAIJ() [1/4]

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

◆ createMPIAIJ() [2/4]

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

Definition at line 751 of file MatrixManager.cpp.

752  {
753  MoFEM::CoreInterface &m_field = cOre;
754  CreateRowCompressedADJMatrix *core_ptr =
755  static_cast<CreateRowCompressedADJMatrix *>(&cOre);
757  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJ, 0, 0, 0, 0);
758 
759  auto problems_ptr = m_field.get_problems();
760  auto &prb = problems_ptr->get<Problem_mi_tag>();
761  auto p_miit = prb.find(name);
762  if (p_miit == prb.end()) {
763  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
764  "problem < %s > is not found (top tip: check spelling)",
765  name.c_str());
766  }
767 
768  std::vector<int> i_vec, j_vec;
769  j_vec.reserve(10000);
770  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
771  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
772 
773  int nb_row_dofs = p_miit->getNbDofsRow();
774  int nb_col_dofs = p_miit->getNbDofsCol();
775  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
776  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
777 
778  auto get_layout = [&]() {
779  int start_ranges, end_ranges;
780  PetscLayout layout;
781  CHKERR PetscLayoutCreate(m_field.get_comm(), &layout);
782  CHKERR PetscLayoutSetBlockSize(layout, 1);
783  CHKERR PetscLayoutSetLocalSize(layout, nb_local_dofs_col);
784  CHKERR PetscLayoutSetUp(layout);
785  CHKERR PetscLayoutGetRange(layout, &start_ranges, &end_ranges);
786  CHKERR PetscLayoutDestroy(&layout);
787  return std::make_pair(start_ranges, end_ranges);
788  };
789 
790  auto get_nnz = [&](auto &d_nnz, auto &o_nnz) {
792  auto layout = get_layout();
793  int j = 0;
794  for (int i = 0; i != nb_local_dofs_row; ++i) {
795  for (; j != i_vec[i + 1]; ++j) {
796  if (j_vec[j] < layout.second && j_vec[j] >= layout.first)
797  ++(d_nnz[i]);
798  else
799  ++(o_nnz[i]);
800  }
801  }
803  };
804 
805  std::vector<int> d_nnz(nb_local_dofs_row, 0), o_nnz(nb_local_dofs_row, 0);
806  CHKERR get_nnz(d_nnz, o_nnz);
807 
808  CHKERR MatCreate(m_field.get_comm(), Aij);
809  CHKERR MatSetSizes(*Aij, nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
810  nb_col_dofs);
811  CHKERR MatSetType(*Aij, MATMPIAIJ);
812  CHKERR MatMPIAIJSetPreallocation(*Aij, 0, &*d_nnz.begin(), 0,
813  &*o_nnz.begin());
814 
815  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJ, 0, 0, 0, 0);
817 }

◆ createMPIAIJ() [3/4]

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

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 123 of file MatrixManager.hpp.

124  {
125  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
126  return 0;
127  }

◆ createMPIAIJ() [4/4]

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

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 132 of file MatrixManager.hpp.

134  {
136  Mat aij;
137  CHKERR createMPIAIJ<Tag>(name, &aij, verb);
138  aij_ptr.reset(aij, false);
140  }

◆ createMPIAIJCUSPARSEWithArrays() [1/4]

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

◆ createMPIAIJCUSPARSEWithArrays() [2/4]

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

Definition at line 662 of file MatrixManager.cpp.

663  {
665 
666  MoFEM::CoreInterface &m_field = cOre;
667  CreateRowCompressedADJMatrix *core_ptr =
668  static_cast<CreateRowCompressedADJMatrix *>(&cOre);
669  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays, 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  j_vec.reserve(10000);
682  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
683  p_miit, MATAIJCUSPARSE, i_vec, j_vec, false, verb);
684 
685  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
686  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
687 
688  auto get_layout = [&]() {
689  int start_ranges, end_ranges;
690  PetscLayout layout;
691  CHKERR PetscLayoutCreate(m_field.get_comm(), &layout);
692  CHKERR PetscLayoutSetBlockSize(layout, 1);
693  CHKERR PetscLayoutSetLocalSize(layout, nb_local_dofs_col);
694  CHKERR PetscLayoutSetUp(layout);
695  CHKERR PetscLayoutGetRange(layout, &start_ranges, &end_ranges);
696  CHKERR PetscLayoutDestroy(&layout);
697  return std::make_pair(start_ranges, end_ranges);
698  };
699 
700  auto get_nnz = [&](auto &d_nnz, auto &o_nnz) {
702  auto layout = get_layout();
703  int j = 0;
704  for (int i = 0; i != nb_local_dofs_row; ++i) {
705  for (; j != i_vec[i + 1]; ++j) {
706  if (j_vec[j] < layout.second && j_vec[j] >= layout.first)
707  ++(d_nnz[i]);
708  else
709  ++(o_nnz[i]);
710  }
711  }
713  };
714 
715  std::vector<int> d_nnz(nb_local_dofs_row, 0), o_nnz(nb_local_dofs_row, 0);
716  CHKERR get_nnz(d_nnz, o_nnz);
717 
718 #ifdef PETSC_HAVE_CUDA
719  CHKERR ::MatCreateAIJCUSPARSE(m_field.get_comm(), nb_local_dofs_row,
720  nb_local_dofs_col, nb_row_dofs, nb_col_dofs, 0,
721  &*d_nnz.begin(), 0, &*o_nnz.begin(), Aij);
722 #else
723  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
724  "Error: To use this matrix type compile PETSc with CUDA.");
725 #endif
726 
727  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays, 0, 0, 0, 0);
728 
730 }

◆ createMPIAIJCUSPARSEWithArrays() [3/4]

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

Definition at line 69 of file MatrixManager.hpp.

70  {
71  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
72  return 0;
73  }

◆ createMPIAIJCUSPARSEWithArrays() [4/4]

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

Definition at line 78 of file MatrixManager.hpp.

80  {
82  Mat aij;
83  CHKERR createMPIAIJCUSPARSEWithArrays<Tag>(name, &aij, verb);
84  aij_ptr.reset(aij, false);
86  }

◆ createMPIAIJWithArrays() [1/4]

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

◆ createMPIAIJWithArrays() [2/4]

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

Definition at line 624 of file MatrixManager.cpp.

625  {
627 
628  MoFEM::CoreInterface &m_field = cOre;
629  CreateRowCompressedADJMatrix *core_ptr =
630  static_cast<CreateRowCompressedADJMatrix *>(&cOre);
631  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
632 
633  auto problems_ptr = m_field.get_problems();
634  auto &prb = problems_ptr->get<Problem_mi_tag>();
635  auto p_miit = prb.find(name);
636  if (p_miit == prb.end()) {
637  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
638  "problem < %s > is not found (top tip: check spelling)",
639  name.c_str());
640  }
641 
642  std::vector<int> i_vec, j_vec;
643  j_vec.reserve(10000);
644  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
645  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
646 
647  int nb_row_dofs = p_miit->getNbDofsRow();
648  int nb_col_dofs = p_miit->getNbDofsCol();
649  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
650  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
651 
652  CHKERR ::MatCreateMPIAIJWithArrays(
653  m_field.get_comm(), nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
654  nb_col_dofs, &*i_vec.begin(), &*j_vec.begin(), PETSC_NULL, Aij);
655 
656  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
658 }

◆ createMPIAIJWithArrays() [3/4]

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

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 49 of file MatrixManager.hpp.

50  {
51  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
52  return 0;
53  }

◆ createMPIAIJWithArrays() [4/4]

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

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 58 of file MatrixManager.hpp.

60  {
62  Mat aij;
63  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
64  aij_ptr.reset(aij, false);
66  }

◆ createSeqAIJCUSPARSEWithArrays() [1/4]

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

◆ createSeqAIJCUSPARSEWithArrays() [2/4]

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

Definition at line 734 of file MatrixManager.cpp.

735  {
737 
738 #ifdef PETSC_HAVE_CUDA
739  // CHKERR ::MatCreateSeqAIJCUSPARSE(MPI_Comm comm, PetscInt m, PetscInt n,
740  // PetscInt nz, const PetscInt nnz[], Mat
741  // *A);
742  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
743  "Not implemented type of matrix yet, try MPI version (aijcusparse)");
744 #endif
745 
747 }

◆ createSeqAIJCUSPARSEWithArrays() [3/4]

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

Definition at line 89 of file MatrixManager.hpp.

90  {
91  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
92  return 0;
93  }

◆ createSeqAIJCUSPARSEWithArrays() [4/4]

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

Definition at line 98 of file MatrixManager.hpp.

100  {
102  Mat aij;
103  CHKERR createSeqAIJCUSPARSEWithArrays<Tag>(name, &aij, verb);
104  aij_ptr.reset(aij, false);
106  }

◆ createSeqAIJWithArrays() [1/4]

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

◆ createSeqAIJWithArrays() [2/4]

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

Definition at line 858 of file MatrixManager.cpp.

859  {
860  MoFEM::CoreInterface &m_field = cOre;
861  CreateRowCompressedADJMatrix *core_ptr =
862  static_cast<CreateRowCompressedADJMatrix *>(&cOre);
864  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
865 
866  auto problems_ptr = m_field.get_problems();
867  auto &prb = problems_ptr->get<Problem_mi_tag>();
868  auto p_miit = prb.find(name);
869  if (p_miit == prb.end()) {
870  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
871  "problem < %s > is not found (top tip: check spelling)",
872  name.c_str());
873  }
874 
875  std::vector<int> i_vec, j_vec;
876  j_vec.reserve(10000);
877  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(p_miit, MATAIJ, i_vec,
878  j_vec, false, verb);
879 
880  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
881  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
882 
883  double *_a;
884  CHKERR PetscMalloc(j_vec.size() * sizeof(double), &_a);
885 
886  Mat tmpMat;
887  CHKERR ::MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, nb_local_dofs_row,
888  nb_local_dofs_col, &*i_vec.begin(),
889  &*j_vec.begin(), _a, &tmpMat);
890  CHKERR MatDuplicate(tmpMat, MAT_SHARE_NONZERO_PATTERN, Aij);
891  CHKERR MatDestroy(&tmpMat);
892 
893  CHKERR PetscFree(_a);
894 
895  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
897 }

◆ createSeqAIJWithArrays() [3/4]

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

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 195 of file MatrixManager.hpp.

196  {
197  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
198  return 0;
199  }

◆ createSeqAIJWithArrays() [4/4]

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

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 224 of file MatrixManager.hpp.

226  {
228  Mat aij;
229  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
230  aij_ptr.reset(aij, false);
232  }

◆ query_interface()

MoFEMErrorCode MoFEM::MatrixManager::query_interface ( boost::typeindex::type_index  type_index,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 600 of file MatrixManager.cpp.

601  {
602  *iface = const_cast<MatrixManager *>(this);
603  return 0;
604 }

Member Data Documentation

◆ cOre

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

Definition at line 26 of file MatrixManager.hpp.

◆ MOFEM_EVENT_checkMatrixFillIn

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_checkMatrixFillIn
private

Definition at line 319 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAdjWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
private

Definition at line 315 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJ

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
private

Definition at line 313 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
private

Definition at line 317 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
private

Definition at line 314 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
private

Definition at line 318 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
private

Definition at line 316 of file MatrixManager.hpp.


The documentation for this struct was generated from the following files:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::dimension_from_handle
auto dimension_from_handle(const EntityHandle h)
get entity dimension form handle
Definition: Templates.hpp:1892
MoFEM::CoreInterface::loop_finite_elements
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, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
Definition: MatrixManager.hpp:317
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
MOFEM_LOG_SEVERITY_SYNC
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
MOFEM_LOG_CHANNEL
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
EntityHandle
NOISY
@ NOISY
Definition: definitions.h:224
MatrixManagerFunctionBegin
#define MatrixManagerFunctionBegin
Create adjacent matrices using different indices.
Definition: MatrixManager.cpp:6
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::BitFieldId
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:42
MOFEM_LOG_ATTRIBUTES
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
Definition: LogManager.hpp:296
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
Definition: MatrixManager.hpp:318
MoFEM::CoreInterface::get_finite_elements
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
MoFEM::MatrixManager::MOFEM_EVENT_checkMatrixFillIn
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
Definition: MatrixManager.hpp:319
sdf.r
int r
Definition: sdf.py:8
ROW
@ ROW
Definition: definitions.h:136
NumeredDofEntity_multiIndex
multi_index_container< boost::shared_ptr< NumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, UId, &NumeredDofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Part_mi_tag >, member< NumeredDofEntity, unsigned int, &NumeredDofEntity::pArt > >, ordered_non_unique< tag< Idx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::dofIdx > >, ordered_non_unique< tag< PetscGlobalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscGloablDofIdx > >, ordered_non_unique< tag< PetscLocalIdx_mi_tag >, member< NumeredDofEntity, DofIdx, &NumeredDofEntity::petscLocalDofIdx > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< NumeredDofEntity::interface_type_DofEntity, EntityHandle, &NumeredDofEntity::getEnt > > > > NumeredDofEntity_multiIndex
MultiIndex container keeps NumeredDofEntity.
Definition: DofsMultiIndices.hpp:469
c
const double c
speed of light (cm/ns)
Definition: initial_diffusion.cpp:39
MoFEM::LogManager::BitScope
@ BitScope
Definition: LogManager.hpp:49
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
Definition: MatrixManager.hpp:314
VERBOSE
@ VERBOSE
Definition: definitions.h:222
FEMethod
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
a
constexpr double a
Definition: approx_sphere.cpp:30
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
MoFEM::Types::UId
uint128_t UId
Unique Id.
Definition: Types.hpp:31
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
Definition: MatrixManager.hpp:315
MoFEM::MatrixManager::MatrixManager
MatrixManager(const MoFEM::Core &core)
Definition: MatrixManager.cpp:606
MoFEM::CoreInterface::get_problems
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
MoFEM::CacheTupleSharedPtr
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
Definition: FEMultiIndices.hpp:495
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
PetscLogEvent MOFEM_EVENT_createMPIAIJ
Definition: MatrixManager.hpp:313
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::DofEntity::getLoFieldEntityUId
static UId getLoFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
Definition: DofsMultiIndices.hpp:69
MOFEM_NOT_FOUND
@ MOFEM_NOT_FOUND
Definition: definitions.h:33
MoFEM::CoreInterface::get_fields
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
Range
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::DofEntity::getHiFieldEntityUId
static UId getHiFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
Definition: DofsMultiIndices.hpp:75
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
IT
constexpr IntegrationType IT
Definition: test_broken_space.cpp:24
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
QUIET
@ QUIET
Definition: definitions.h:221
MoFEM::CoreInterface
Interface.
Definition: Interface.hpp:27
MF_EXIST
@ MF_EXIST
Definition: definitions.h:113
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::MatrixManager::cOre
MoFEM::Core & cOre
Definition: MatrixManager.hpp:26
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
VERY_VERBOSE
@ VERY_VERBOSE
Definition: definitions.h:223
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::MatrixManager::checkMatrixFillIn
MoFEMErrorCode checkMatrixFillIn(const std::string problem_name, int row_print, int col_print, Mat A, int verb=QUIET)
check if matrix fill in correspond to finite element indices
Definition: MatrixManager.cpp:899
MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
Definition: MatrixManager.hpp:316