v0.13.2
Loading...
Searching...
No Matches
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
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

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

593 : cOre(const_cast<MoFEM::Core &>(core)) {
594 PetscLogEventRegister("MatrixManagerCreateMPIAIJ", 0,
596 PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
598 PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
600 PetscLogEventRegister("MatrixManagerCreateMPIAIJCUSPARSEWithArrays", 0,
602 PetscLogEventRegister("MatrixManagerCreateSeqAIJCUSPARSEWithArrays", 0,
604 PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
606 PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
608}
Core (interface) class.
Definition: Core.hpp:82
PetscLogEvent MOFEM_EVENT_createSeqAIJCUSPARSEWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
PetscLogEvent MOFEM_EVENT_createMPIAIJ
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays
PetscLogEvent MOFEM_EVENT_createMPIAIJCUSPARSEWithArrays

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

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

◆ checkMPIAIJMatrixFillIn() [1/3]

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

Definition at line 1142 of file MatrixManager.cpp.

1143 {
1144 MoFEM::CoreInterface &m_field = cOre;
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}
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

◆ checkMPIAIJMatrixFillIn() [2/3]

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

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

Definition at line 1128 of file MatrixManager.cpp.

1130 {
1131 MoFEM::CoreInterface &m_field = cOre;
1133 // create matrix
1134 SmartPetscObj<Mat> A;
1135 CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1136 CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1137 CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1139}

◆ checkMPIAIJWithArraysMatrixFillIn() [2/3]

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

◆ checkMPIAIJWithArraysMatrixFillIn() [3/3]

template<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 
)

Definition at line 808 of file MatrixManager.cpp.

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

◆ createMPIAdjWithArrays() [2/3]

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

◆ 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");
170 }
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32

◆ createMPIAIJ() [1/4]

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

Definition at line 738 of file MatrixManager.cpp.

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

◆ createMPIAIJ() [2/4]

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

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

Definition at line 648 of file MatrixManager.cpp.

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

◆ createMPIAIJCUSPARSEWithArrays() [2/4]

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

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

Definition at line 611 of file MatrixManager.cpp.

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

◆ createMPIAIJWithArrays() [2/4]

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

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

Definition at line 722 of file MatrixManager.cpp.

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

◆ createSeqAIJCUSPARSEWithArrays() [2/4]

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

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

Definition at line 846 of file MatrixManager.cpp.

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

◆ createSeqAIJWithArrays() [2/4]

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

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

587 {
588 *iface = const_cast<MatrixManager *>(this);
589 return 0;
590}
MatrixManager(const MoFEM::Core &core)

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: