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

598  : cOre(const_cast<MoFEM::Core &>(core)) {
599  PetscLogEventRegister("MatMngCrtMPIAIJ", 0, &MOFEM_EVENT_createMPIAIJ);
600  PetscLogEventRegister("MatMngCrtMPIAIJWthArr", 0,
602  PetscLogEventRegister("MatMngCrtMPIAdjWithArr", 0,
604  PetscLogEventRegister("MatMngCrtMPIAIJCUSPARSEWthArr", 0,
606  PetscLogEventRegister("MatMngCrtSeqAIJCUSPARSEWthArrs", 0,
608  PetscLogEventRegister("MatMngCrtSeqAIJWthArrs", 0,
610  PetscLogEventRegister("MatMngCrtCheckMPIAIJWthArrFillIn", 0,
612 }

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

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

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

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

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

1132  {
1134  // create matrix
1135  SmartPetscObj<Mat> A;
1136  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1137  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1138  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1140 }

◆ 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  }

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

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

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

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

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

654  {
656 
657  MoFEM::CoreInterface &m_field = cOre;
658  CreateRowComressedADJMatrix *core_ptr =
659  static_cast<CreateRowComressedADJMatrix *>(&cOre);
660  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays, 0, 0, 0, 0);
661 
662  auto problems_ptr = m_field.get_problems();
663  auto &prb = problems_ptr->get<Problem_mi_tag>();
664  auto p_miit = prb.find(name);
665  if (p_miit == prb.end()) {
666  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
667  "problem < %s > is not found (top tip: check spelling)",
668  name.c_str());
669  }
670 
671  std::vector<int> i_vec, j_vec;
672  j_vec.reserve(10000);
673  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
674  p_miit, MATAIJCUSPARSE, i_vec, j_vec, false, verb);
675 
676  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
677  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
678 
679  auto get_layout = [&]() {
680  int start_ranges, end_ranges;
681  PetscLayout layout;
682  CHKERR PetscLayoutCreate(m_field.get_comm(), &layout);
683  CHKERR PetscLayoutSetBlockSize(layout, 1);
684  CHKERR PetscLayoutSetLocalSize(layout, nb_local_dofs_col);
685  CHKERR PetscLayoutSetUp(layout);
686  CHKERR PetscLayoutGetRange(layout, &start_ranges, &end_ranges);
687  CHKERR PetscLayoutDestroy(&layout);
688  return std::make_pair(start_ranges, end_ranges);
689  };
690 
691  auto get_nnz = [&](auto &d_nnz, auto &o_nnz) {
693  auto layout = get_layout();
694  int j = 0;
695  for (int i = 0; i != nb_local_dofs_row; ++i) {
696  for (; j != i_vec[i + 1]; ++j) {
697  if (j_vec[j] < layout.second && j_vec[j] >= layout.first)
698  ++(d_nnz[i]);
699  else
700  ++(o_nnz[i]);
701  }
702  }
704  };
705 
706  std::vector<int> d_nnz(nb_local_dofs_row, 0), o_nnz(nb_local_dofs_row, 0);
707  CHKERR get_nnz(d_nnz, o_nnz);
708 
709 #ifdef PETSC_HAVE_CUDA
710  CHKERR ::MatCreateAIJCUSPARSE(m_field.get_comm(), nb_local_dofs_row,
711  nb_local_dofs_col, nb_row_dofs, nb_col_dofs, 0,
712  &*d_nnz.begin(), 0, &*o_nnz.begin(), Aij);
713 #else
714  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
715  "Error: To use this matrix type compile PETSc with CUDA.");
716 #endif
717 
718  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays, 0, 0, 0, 0);
719 
721 }

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

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

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

726  {
728 
729 #ifdef PETSC_HAVE_CUDA
730  // CHKERR ::MatCreateSeqAIJCUSPARSE(MPI_Comm comm, PetscInt m, PetscInt n,
731  // PetscInt nz, const PetscInt nnz[], Mat *A);
732  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
733  "Not implemented type of matrix yet, try MPI version (aijcusparse)");
734 #endif
735 
737 }

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

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

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

592  {
593  *iface = const_cast<MatrixManager *>(this);
594  return 0;
595 }

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

◆ MOFEM_EVENT_createMPIAdjWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
private

Definition at line 296 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJ

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
private

Definition at line 294 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays
private

Definition at line 298 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
private

Definition at line 295 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
private

Definition at line 299 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
private

Definition at line 297 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:447
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:298
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
NOISY
@ NOISY
Definition: definitions.h:211
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_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:299
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:300
MoFEM::LogManager::BitScope
@ BitScope
Definition: LogManager.hpp:49
MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
Definition: MatrixManager.hpp:295
FEMethod
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
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:296
MoFEM::MatrixManager::MatrixManager
MatrixManager(const MoFEM::Core &core)
Definition: MatrixManager.cpp:597
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:294
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MOFEM_NOT_FOUND
@ MOFEM_NOT_FOUND
Definition: definitions.h:33
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
MoFEM::CoreInterface
Interface.
Definition: Interface.hpp:27
MF_EXIST
@ MF_EXIST
Definition: definitions.h:100
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
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:210
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
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:889
MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
Definition: MatrixManager.hpp:297