v0.13.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 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 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 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)
 
- 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 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_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.

Definition at line 31 of file MatrixManager.hpp.

Constructor & Destructor Documentation

◆ MatrixManager()

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

Definition at line 568 of file MatrixManager.cpp.

569  : cOre(const_cast<MoFEM::Core &>(core)) {
570  PetscLogEventRegister("MatrixManagerCreateMPIAIJ", 0,
572  PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
574  PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
576  PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
578  PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
580 }
Core (interface) class.
Definition: Core.hpp:92
PetscLogEvent MOFEM_EVENT_createMPIAdjWithArrays
PetscLogEvent MOFEM_EVENT_checkMatrixFillIn
PetscLogEvent MOFEM_EVENT_createMPIAIJ
PetscLogEvent MOFEM_EVENT_createMPIAIJWithArrays
PetscLogEvent MOFEM_EVENT_createSeqAIJWithArrays

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

770  {
771  MoFEM::CoreInterface &m_field = cOre;
773 
774  PetscLogEventBegin(MOFEM_EVENT_checkMatrixFillIn, 0, 0, 0, 0);
775 
776  struct TestMatrixFillIn : public FEMethod {
777  CoreInterface *mFieldPtr;
778 
779  Mat A;
780 
781  int rowPrint, colPrint;
782 
783  TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
784  int col_print)
785  : mFieldPtr(m_field_ptr), A(a), rowPrint(row_print),
786  colPrint(col_print){};
787 
788  MoFEMErrorCode preProcess() {
791  }
792 
793  MoFEMErrorCode operator()() {
795 
796  if (refinedFiniteElementsPtr->find(
797  numeredEntFiniteElementPtr->getEnt()) ==
798  refinedFiniteElementsPtr->end()) {
799  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
800  "data inconsistency");
801  }
802 
803  auto row_dofs = getRowDofsPtr();
804  auto col_dofs = getColDofsPtr();
805 
806  for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
807 
808  if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
809  refinedEntitiesPtr->end()) {
810  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
811  "data inconsistency");
812  }
813  if (!(*cit)->getActive()) {
814  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
815  "data inconsistency");
816  }
817 
819  Composite_Unique_mi_tag>::type::iterator ait;
820  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
821  boost::make_tuple((*cit)->getFieldEntityPtr()->getLocalUniqueId(),
822  numeredEntFiniteElementPtr->getLocalUniqueId()));
823  if (ait == adjacenciesPtr->end()) {
824  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
825  "adjacencies data inconsistency");
826  } else {
827  UId uid = ait->getEntUniqueId();
828  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
829  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
830  "data inconsistency");
831  }
832  if (dofsPtr->find((*cit)->getLocalUniqueId()) == dofsPtr->end()) {
833  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
834  "data inconsistency");
835  }
836  }
837 
838  if ((*cit)->getEntType() != MBVERTEX) {
839 
840  auto range =
841  col_dofs->get<Ent_mi_tag>().equal_range((*cit)->getEnt());
842  int nb_dofs_on_ent = std::distance(range.first, range.second);
843 
844  int max_order = (*cit)->getMaxOrder();
845  if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
846  nb_dofs_on_ent) {
847  MOFEM_LOG("SELF", Sev::warning)
848  << "Warning: Number of Dofs in Col diffrent than number "
849  "of dofs for given entity order "
850  << (*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order)
851  << " " << nb_dofs_on_ent;
852  }
853  }
854  }
855 
856  for (auto rit = row_dofs->begin(); rit != row_dofs->end(); rit++) {
857 
858  if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
859  refinedEntitiesPtr->end()) {
860  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
861  "data inconsistency");
862  }
863  if (!(*rit)->getActive()) {
864  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
865  "data inconsistency");
866  }
867 
869  Composite_Unique_mi_tag>::type::iterator ait;
870  ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
871  boost::make_tuple((*rit)->getFieldEntityPtr()->getLocalUniqueId(),
872  numeredEntFiniteElementPtr->getLocalUniqueId()));
873  if (ait == adjacenciesPtr->end()) {
875  MOFEM_LOG("SELF", Sev::error) << *(*rit);
876  MOFEM_LOG("SELF", Sev::error) << *(*rit);
877  MOFEM_LOG("SELF", Sev::error) << *numeredEntFiniteElementPtr;
878  MOFEM_LOG("SELF", Sev::error) << "dof: " << (*rit)->getBitRefLevel();
879  MOFEM_LOG("SELF", Sev::error)
880  << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel();
881  MOFEM_LOG("SELF", Sev::error)
882  << "problem: " << problemPtr->getBitRefLevel();
883  MOFEM_LOG("SELF", Sev::error)
884  << "problem mask: " << problemPtr->getBitRefLevelMask();
885  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
886  "adjacencies data inconsistency");
887  } else {
888  UId uid = ait->getEntUniqueId();
889  if (entitiesPtr->find(uid) == entitiesPtr->end()) {
890  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
891  "data inconsistency");
892  }
893  if (dofsPtr->find((*rit)->getLocalUniqueId()) == dofsPtr->end()) {
894  SETERRQ(mFieldPtr->get_comm(), MOFEM_DATA_INCONSISTENCY,
895  "data inconsistency");
896  }
897  }
898  int row = (*rit)->getPetscGlobalDofIdx();
899 
900  auto col_dofs = getColDofsPtr();
901  for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
902 
903  int col = (*cit)->getPetscGlobalDofIdx();
904 
905  if (row == rowPrint && col == colPrint) {
906  MOFEM_LOG("SELF", Sev::noisy) << "fe:\n"
907  << *numeredEntFiniteElementPtr;
908  MOFEM_LOG("SELF", Sev::noisy) << "row:\n" << *(*rit);
909  MOFEM_LOG("SELF", Sev::noisy) << "col:\n" << *(*cit);
910  MOFEM_LOG("SELF", Sev::noisy)
911  << "fe:\n"
912  << numeredEntFiniteElementPtr->getBitRefLevel();
913  MOFEM_LOG("SELF", Sev::noisy) << "row:\n"
914  << (*rit)->getBitRefLevel();
915  MOFEM_LOG("SELF", Sev::noisy) << "col:\n"
916  << (*cit)->getBitRefLevel();
917  }
918 
919  CHKERR MatSetValue(A, row, col, 1, INSERT_VALUES);
920  }
921 
922  if ((*rit)->getEntType() != MBVERTEX) {
923 
924  auto range =
925  row_dofs->get<Ent_mi_tag>().equal_range((*rit)->getEnt());
926  int nb_dofs_on_ent = std::distance(range.first, range.second);
927 
928  int max_order = (*rit)->getMaxOrder();
929  if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
930  nb_dofs_on_ent) {
931  MOFEM_LOG("SELF", Sev::warning)
932  << "Warning: Number of Dofs in Row diffrent than number "
933  "of dofs for given entity order "
934  << (*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order)
935  << " " << nb_dofs_on_ent;
936  }
937  }
938  }
939 
941  }
942 
943  MoFEMErrorCode postProcess() {
945 
946  // cerr << mFieldPtr->get_comm_rank() << endl;
947  CHKERR MatAssemblyBegin(A, MAT_FLUSH_ASSEMBLY);
948  CHKERR MatAssemblyEnd(A, MAT_FLUSH_ASSEMBLY);
949 
951  }
952  };
953 
954  // create matrix
955  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
956 
957  if (verb >= VERY_VERBOSE) {
958  MatView(A, PETSC_VIEWER_STDOUT_WORLD);
959  }
960 
961  if (verb >= NOISY) {
962  MatView(A, PETSC_VIEWER_DRAW_WORLD);
963  std::string wait;
964  std::cin >> wait;
965  }
966 
967  TestMatrixFillIn method(&m_field, A, row_print, col_print);
968 
969  // get problem
970  auto problems_ptr = m_field.get_problems();
971  auto &prb_set = problems_ptr->get<Problem_mi_tag>();
972  auto p_miit = prb_set.find(problem_name);
973  if (p_miit == prb_set.end())
974  SETERRQ1(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
975  "problem < %s > not found (top tip: check spelling)",
976  problem_name.c_str());
977  MOFEM_LOG_C("WORLD", Sev::inform, "check problem < %s >",
978  problem_name.c_str());
979 
980  // loop all elements in problem and check if assemble is without error
981  auto fe_ptr = m_field.get_finite_elements();
982  for (auto &fe : *fe_ptr) {
983  MOFEM_LOG_C("WORLD", Sev::verbose, "\tcheck element %s",
984  fe->getName().c_str());
985  CHKERR m_field.loop_finite_elements(problem_name, fe->getName(), method,
986  nullptr, MF_EXIST,
987  CacheTupleSharedPtr(), verb);
988  }
989 
990  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
991  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
992 
993  PetscLogEventEnd(MOFEM_EVENT_checkMatrixFillIn, 0, 0, 0, 0);
994 
996 }
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
#define MatrixManagerFunctionBegin
Create adjacent matrices using different indices.
constexpr double a
@ VERY_VERBOSE
Definition: definitions.h:223
@ NOISY
Definition: definitions.h:224
@ MF_EXIST
Definition: definitions.h:113
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
virtual const Problem_multiIndex * get_problems() const =0
Get the problems object.
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
Definition: LogManager.hpp:299
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:67
uint128_t UId
Unique Id.
Definition: Types.hpp:42
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
double 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 1013 of file MatrixManager.cpp.

1014  {
1015  MoFEM::CoreInterface &m_field = cOre;
1017  // create matrix
1018  SmartPetscObj<Mat> A;
1019  CHKERR createMPIAIJ<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1020  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1021  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1023 }
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 
)

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

258  {
259  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
260  return 0;
261  }

◆ checkMPIAIJWithArraysMatrixFillIn() [1/3]

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

Definition at line 999 of file MatrixManager.cpp.

1001  {
1002  MoFEM::CoreInterface &m_field = cOre;
1004  // create matrix
1005  SmartPetscObj<Mat> A;
1006  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1007  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1008  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1010 }

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

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

237  {
238  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
239  return 0;
240  }

◆ createMPIAdjWithArrays() [1/3]

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

Definition at line 689 of file MatrixManager.cpp.

691  {
692  MoFEM::CoreInterface &m_field = cOre;
693  CreateRowComressedADJMatrix *core_ptr =
694  static_cast<CreateRowComressedADJMatrix *>(&cOre);
696  PetscLogEventBegin(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
697 
698  auto problems_ptr = m_field.get_problems();
699  auto &prb = problems_ptr->get<Problem_mi_tag>();
700  auto p_miit = prb.find(name);
701  if (p_miit == prb.end()) {
702  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
703  "problem < %s > is not found (top tip: check spelling)",
704  name.c_str());
705  }
706 
707  std::vector<int> i_vec, j_vec;
708  j_vec.reserve(10000);
709  CHKERR core_ptr->createMatArrays<Idx_mi_tag>(p_miit, MATMPIADJ, i_vec, j_vec,
710  true, verb);
711  int *_i, *_j;
712  CHKERR PetscMalloc(i_vec.size() * sizeof(int), &_i);
713  CHKERR PetscMalloc(j_vec.size() * sizeof(int), &_j);
714  copy(i_vec.begin(), i_vec.end(), _i);
715  copy(j_vec.begin(), j_vec.end(), _j);
716 
717  int nb_col_dofs = p_miit->getNbDofsCol();
718  CHKERR MatCreateMPIAdj(m_field.get_comm(), i_vec.size() - 1, nb_col_dofs, _i,
719  _j, PETSC_NULL, Adj);
720  CHKERR MatSetOption(*Adj, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
721 
722  PetscLogEventEnd(MOFEM_EVENT_createMPIAdjWithArrays, 0, 0, 0, 0);
724 }
@ MOFEM_NOT_FOUND
Definition: definitions.h:46

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

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

137  {
138  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
139  return MOFEM_NOT_IMPLEMENTED;
140  }
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45

◆ createMPIAIJ() [1/4]

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

Definition at line 619 of file MatrixManager.cpp.

621  {
622  MoFEM::CoreInterface &m_field = cOre;
623  CreateRowComressedADJMatrix *core_ptr =
624  static_cast<CreateRowComressedADJMatrix *>(&cOre);
626  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJ, 0, 0, 0, 0);
627 
628  auto problems_ptr = m_field.get_problems();
629  auto &prb = problems_ptr->get<Problem_mi_tag>();
630  auto p_miit = prb.find(name);
631  if (p_miit == prb.end()) {
632  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
633  "problem < %s > is not found (top tip: check spelling)",
634  name.c_str());
635  }
636 
637  std::vector<int> i_vec, j_vec;
638  j_vec.reserve(10000);
639  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
640  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
641 
642  int nb_row_dofs = p_miit->getNbDofsRow();
643  int nb_col_dofs = p_miit->getNbDofsCol();
644  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
645  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
646 
647  auto get_layout = [&]() {
648  int start_ranges, end_ranges;
649  PetscLayout layout;
650  CHKERR PetscLayoutCreate(m_field.get_comm(), &layout);
651  CHKERR PetscLayoutSetBlockSize(layout, 1);
652  CHKERR PetscLayoutSetLocalSize(layout, nb_local_dofs_col);
653  CHKERR PetscLayoutSetUp(layout);
654  CHKERR PetscLayoutGetRange(layout, &start_ranges, &end_ranges);
655  CHKERR PetscLayoutDestroy(&layout);
656  return std::make_pair(start_ranges, end_ranges);
657  };
658 
659  auto get_nnz = [&](auto &d_nnz, auto &o_nnz) {
661  auto layout = get_layout();
662  int j = 0;
663  for (int i = 0; i != nb_local_dofs_row; ++i) {
664  for (; j != i_vec[i + 1]; ++j) {
665  if (j_vec[j] < layout.second && j_vec[j] >= layout.first)
666  ++(d_nnz[i]);
667  else
668  ++(o_nnz[i]);
669  }
670  }
672  };
673 
674  std::vector<int> d_nnz(nb_local_dofs_row, 0), o_nnz(nb_local_dofs_row, 0);
675  CHKERR get_nnz(d_nnz, o_nnz);
676 
677  CHKERR MatCreate(m_field.get_comm(), Aij);
678  CHKERR MatSetSizes(*Aij, nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
679  nb_col_dofs);
680  CHKERR MatSetType(*Aij, MATMPIAIJ);
681  CHKERR MatMPIAIJSetPreallocation(*Aij, 0, &*d_nnz.begin(), 0,
682  &*o_nnz.begin());
683 
684  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJ, 0, 0, 0, 0);
686 }
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 
)

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

94  {
95  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
96  return 0;
97  }

◆ createMPIAIJ() [4/4]

template<class Tag >
MoFEMErrorCode MoFEM::MatrixManager::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.

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij
verb
Returns
MoFEMErrorCode

Definition at line 102 of file MatrixManager.hpp.

104  {
106  Mat aij;
107  CHKERR createMPIAIJ<Tag>(name, &aij, verb);
108  aij_ptr.reset(aij, false);
110  }

◆ createMPIAIJWithArrays() [1/4]

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

Definition at line 583 of file MatrixManager.cpp.

584  {
585  MoFEM::CoreInterface &m_field = cOre;
586  CreateRowComressedADJMatrix *core_ptr =
587  static_cast<CreateRowComressedADJMatrix *>(&cOre);
589  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
590 
591  auto problems_ptr = m_field.get_problems();
592  auto &prb = problems_ptr->get<Problem_mi_tag>();
593  auto p_miit = prb.find(name);
594  if (p_miit == prb.end()) {
595  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
596  "problem < %s > is not found (top tip: check spelling)",
597  name.c_str());
598  }
599 
600  std::vector<int> i_vec, j_vec;
601  j_vec.reserve(10000);
602  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(
603  p_miit, MATMPIAIJ, i_vec, j_vec, false, verb);
604 
605  int nb_row_dofs = p_miit->getNbDofsRow();
606  int nb_col_dofs = p_miit->getNbDofsCol();
607  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
608  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
609 
610  CHKERR ::MatCreateMPIAIJWithArrays(
611  m_field.get_comm(), nb_local_dofs_row, nb_local_dofs_col, nb_row_dofs,
612  nb_col_dofs, &*i_vec.begin(), &*j_vec.begin(), PETSC_NULL, Aij);
613 
614  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
616 }

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

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

60  {
61  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
62  return 0;
63  }

◆ createMPIAIJWithArrays() [4/4]

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

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

See PETSc for details

Template Parameters
Tag
Parameters
name
Aij
verb
Returns
MoFEMErrorCode

Definition at line 68 of file MatrixManager.hpp.

70  {
72  Mat aij;
73  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
74  aij_ptr.reset(aij, false);
76  }

◆ createSeqAIJWithArrays() [1/4]

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

Definition at line 727 of file MatrixManager.cpp.

728  {
729  MoFEM::CoreInterface &m_field = cOre;
730  CreateRowComressedADJMatrix *core_ptr =
731  static_cast<CreateRowComressedADJMatrix *>(&cOre);
733  PetscLogEventBegin(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
734 
735  auto problems_ptr = m_field.get_problems();
736  auto &prb = problems_ptr->get<Problem_mi_tag>();
737  auto p_miit = prb.find(name);
738  if (p_miit == prb.end()) {
739  SETERRQ1(m_field.get_comm(), MOFEM_NOT_FOUND,
740  "problem < %s > is not found (top tip: check spelling)",
741  name.c_str());
742  }
743 
744  std::vector<int> i_vec, j_vec;
745  j_vec.reserve(10000);
746  CHKERR core_ptr->createMatArrays<PetscGlobalIdx_mi_tag>(p_miit, MATAIJ, i_vec,
747  j_vec, false, verb);
748 
749  int nb_local_dofs_row = p_miit->getNbLocalDofsRow();
750  int nb_local_dofs_col = p_miit->getNbLocalDofsCol();
751 
752  double *_a;
753  CHKERR PetscMalloc(j_vec.size() * sizeof(double), &_a);
754 
755  Mat tmpMat;
756  CHKERR ::MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, nb_local_dofs_row,
757  nb_local_dofs_col, &*i_vec.begin(),
758  &*j_vec.begin(), _a, &tmpMat);
759  CHKERR MatDuplicate(tmpMat, MAT_SHARE_NONZERO_PATTERN, Aij);
760  CHKERR MatDestroy(&tmpMat);
761 
762  CHKERR PetscFree(_a);
763 
764  PetscLogEventEnd(MOFEM_EVENT_createMPIAIJWithArrays, 0, 0, 0, 0);
766 }

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

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

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

◆ createSeqAIJWithArrays() [4/4]

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

Create sequencial matrix.

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

See PETSc for details

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

Definition at line 194 of file MatrixManager.hpp.

196  {
198  Mat aij;
199  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
200  aij_ptr.reset(aij, false);
202  }

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 562 of file MatrixManager.cpp.

563  {
564  *iface = const_cast<MatrixManager *>(this);
565  return 0;
566 }
MatrixManager(const MoFEM::Core &core)

Member Data Documentation

◆ cOre

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

Definition at line 36 of file MatrixManager.hpp.

◆ MOFEM_EVENT_checkMatrixFillIn

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_checkMatrixFillIn
private

Definition at line 268 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAdjWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
private

Definition at line 266 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJ

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
private

Definition at line 264 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
private

Definition at line 265 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
private

Definition at line 267 of file MatrixManager.hpp.


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