v0.12.1
Public Member Functions | Public Attributes | Private Attributes | List of all members
MoFEM::MatrixManager Struct Reference

Matrix manager is used to build and partition problems. More...

#include <src/interfaces/MatrixManager.hpp>

Inheritance diagram for MoFEM::MatrixManager:
[legend]
Collaboration diagram for MoFEM::MatrixManager:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 MatrixManager (const MoFEM::Core &core)
 
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 (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

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...
 
- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Matrix manager is used to build and partition problems.

Definition at line 34 of file MatrixManager.hpp.

Constructor & Destructor Documentation

◆ MatrixManager()

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

Definition at line 557 of file MatrixManager.cpp.

558  : cOre(const_cast<MoFEM::Core &>(core)) {
559  PetscLogEventRegister("MatrixManagerCreateMPIAIJ", 0,
561  PetscLogEventRegister("MatrixManagerCreateMPIAIJWithArrays", 0,
563  PetscLogEventRegister("MatrixManagerCreateMPIAdjWithArrays", 0,
565  PetscLogEventRegister("MatrixManagerCreateSeqAIJWithArrays", 0,
567  PetscLogEventRegister("MatrixManagerCheckMPIAIJWithArraysMatrixFillIn", 0,
569 }
Core (interface) class.
Definition: Core.hpp:93
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 757 of file MatrixManager.cpp.

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

1003  {
1004  MoFEM::CoreInterface &m_field = cOre;
1006  // create matrix
1007  SmartPetscObj<Mat> A;
1008  CHKERR createMPIAIJ<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
1009  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1010  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
1012 }
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 259 of file MatrixManager.hpp.

261  {
262  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
263  return 0;
264  }

◆ checkMPIAIJWithArraysMatrixFillIn() [1/3]

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

Definition at line 988 of file MatrixManager.cpp.

990  {
991  MoFEM::CoreInterface &m_field = cOre;
993  // create matrix
994  SmartPetscObj<Mat> A;
995  CHKERR createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(problem_name, A, verb);
996  CHKERR MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
997  CHKERR checkMatrixFillIn(problem_name, row_print, col_print, A, verb);
999 }

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

240  {
241  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
242  return 0;
243  }

◆ createMPIAdjWithArrays() [1/3]

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

Definition at line 678 of file MatrixManager.cpp.

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

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

140  {
141  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
142  return MOFEM_NOT_IMPLEMENTED;
143  }
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:125

◆ createMPIAIJ() [1/4]

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

Definition at line 608 of file MatrixManager.cpp.

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

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

97  {
98  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
99  return 0;
100  }

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

107  {
109  Mat aij;
110  CHKERR createMPIAIJ<Tag>(name, &aij, verb);
111  aij_ptr.reset(aij, false);
113  }

◆ createMPIAIJWithArrays() [1/4]

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

Definition at line 572 of file MatrixManager.cpp.

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

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

63  {
64  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
65  return 0;
66  }

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

73  {
75  Mat aij;
76  CHKERR createMPIAIJWithArrays<Tag>(name, &aij, verb);
77  aij_ptr.reset(aij, false);
79  }

◆ createSeqAIJWithArrays() [1/4]

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

Definition at line 716 of file MatrixManager.cpp.

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

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

169  {
170  static_assert(!std::is_same<Tag, Tag>::value, "not implemented");
171  return 0;
172  }

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

199  {
201  Mat aij;
202  CHKERR createSeqAIJWithArrays<Tag>(name, &aij, verb);
203  aij_ptr.reset(aij, false);
205  }

◆ query_interface()

MoFEMErrorCode MoFEM::MatrixManager::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 545 of file MatrixManager.cpp.

546  {
548  *iface = NULL;
549  if (uuid == IDD_MOFEMMatrixManager) {
550  *iface = const_cast<MatrixManager *>(this);
552  }
553  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
555 }
static const MOFEMuuid IDD_MOFEMMatrixManager
MatrixManager(const MoFEM::Core &core)

Member Data Documentation

◆ cOre

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

Definition at line 39 of file MatrixManager.hpp.

◆ MOFEM_EVENT_checkMatrixFillIn

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_checkMatrixFillIn
private

Definition at line 271 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAdjWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAdjWithArrays
private

Definition at line 269 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJ

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJ
private

Definition at line 267 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createMPIAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createMPIAIJWithArrays
private

Definition at line 268 of file MatrixManager.hpp.

◆ MOFEM_EVENT_createSeqAIJWithArrays

PetscLogEvent MoFEM::MatrixManager::MOFEM_EVENT_createSeqAIJWithArrays
private

Definition at line 270 of file MatrixManager.hpp.


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