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
892 {
895
897
898 struct TestMatrixFillIn :
public FEMethod {
899 CoreInterface *mFieldPtr;
900
902
903 int rowPrint, colPrint;
904
905 TestMatrixFillIn(CoreInterface *m_field_ptr, Mat a, int row_print,
906 int col_print)
907 : mFieldPtr(m_field_ptr),
A(
a), rowPrint(row_print),
908 colPrint(col_print){};
909
913 }
914
917
918 if (refinedFiniteElementsPtr->find(
919 numeredEntFiniteElementPtr->getEnt()) ==
920 refinedFiniteElementsPtr->end()) {
922 "data inconsistency");
923 }
924
925 auto row_dofs = getRowDofsPtr();
926 auto col_dofs = getColDofsPtr();
927
928 for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
929
930 if (refinedEntitiesPtr->find((*cit)->getEnt()) ==
931 refinedEntitiesPtr->end()) {
933 "data inconsistency");
934 }
935 if (!(*cit)->getActive()) {
937 "data inconsistency");
938 }
939
940 FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
941 Composite_Unique_mi_tag>::type::iterator ait;
942 ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
943 boost::make_tuple((*cit)->getFieldEntityPtr()->getLocalUniqueId(),
944 numeredEntFiniteElementPtr->getLocalUniqueId()));
945 if (ait == adjacenciesPtr->end()) {
947 "adjacencies data inconsistency");
948 } else {
949 UId uid = ait->getEntUniqueId();
950 if (entitiesPtr->find(uid) == entitiesPtr->end()) {
952 "data inconsistency");
953 }
954 if (dofsPtr->find((*cit)->getLocalUniqueId()) == dofsPtr->end()) {
956 "data inconsistency");
957 }
958 }
959
960 if ((*cit)->getEntType() != MBVERTEX) {
961
962 auto range =
963 col_dofs->get<Ent_mi_tag>().equal_range((*cit)->getEnt());
964 int nb_dofs_on_ent = std::distance(range.first, range.second);
965
966 int max_order = (*cit)->getMaxOrder();
967 if ((*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order) !=
968 nb_dofs_on_ent) {
969
970
971
972
973
975 << "Warning: Number of Dofs in Col different than number "
976 "of dofs for given entity order "
977 << (*cit)->getNbOfCoeffs() * (*cit)->getOrderNbDofs(max_order)
978 << " " << nb_dofs_on_ent;
979 }
980 }
981 }
982
983 for (auto rit = row_dofs->begin(); rit != row_dofs->end(); rit++) {
984
985 if (refinedEntitiesPtr->find((*rit)->getEnt()) ==
986 refinedEntitiesPtr->end()) {
988 "data inconsistency");
989 }
990 if (!(*rit)->getActive()) {
992 "data inconsistency");
993 }
994
995 FieldEntityEntFiniteElementAdjacencyMap_multiIndex::index<
996 Composite_Unique_mi_tag>::type::iterator ait;
997 ait = adjacenciesPtr->get<Composite_Unique_mi_tag>().find(
998 boost::make_tuple((*rit)->getFieldEntityPtr()->getLocalUniqueId(),
999 numeredEntFiniteElementPtr->getLocalUniqueId()));
1000 if (ait == adjacenciesPtr->end()) {
1002 MOFEM_LOG(
"SELF", Sev::error) << *(*rit);
1003 MOFEM_LOG(
"SELF", Sev::error) << *(*rit);
1004 MOFEM_LOG(
"SELF", Sev::error) << *numeredEntFiniteElementPtr;
1005 MOFEM_LOG(
"SELF", Sev::error) <<
"dof: " << (*rit)->getBitRefLevel();
1007 << "fe: " << numeredEntFiniteElementPtr->getBitRefLevel();
1009 << "problem: " << problemPtr->getBitRefLevel();
1011 << "problem mask: " << problemPtr->getBitRefLevelMask();
1013 "adjacencies data inconsistency");
1014 } else {
1015 UId uid = ait->getEntUniqueId();
1016 if (entitiesPtr->find(uid) == entitiesPtr->end()) {
1018 "data inconsistency");
1019 }
1020 if (dofsPtr->find((*rit)->getLocalUniqueId()) == dofsPtr->end()) {
1022 "data inconsistency");
1023 }
1024 }
1025 int row = (*rit)->getPetscGlobalDofIdx();
1026
1027 auto col_dofs = getColDofsPtr();
1028 for (auto cit = col_dofs->begin(); cit != col_dofs->end(); cit++) {
1029
1030 int col = (*cit)->getPetscGlobalDofIdx();
1031
1032 if (row == rowPrint && col == colPrint) {
1033 MOFEM_LOG(
"SELF", Sev::noisy) <<
"fe:\n"
1034 << *numeredEntFiniteElementPtr;
1035 MOFEM_LOG(
"SELF", Sev::noisy) <<
"row:\n" << *(*rit);
1036 MOFEM_LOG(
"SELF", Sev::noisy) <<
"col:\n" << *(*cit);
1038 << "fe:\n"
1039 << numeredEntFiniteElementPtr->getBitRefLevel();
1040 MOFEM_LOG(
"SELF", Sev::noisy) <<
"row:\n"
1041 << (*rit)->getBitRefLevel();
1042 MOFEM_LOG(
"SELF", Sev::noisy) <<
"col:\n"
1043 << (*cit)->getBitRefLevel();
1044 }
1045
1046 CHKERR MatSetValue(
A, row, col, 1, INSERT_VALUES);
1047 }
1048
1049 if ((*rit)->getEntType() != MBVERTEX) {
1050
1051 auto range =
1052 row_dofs->get<Ent_mi_tag>().equal_range((*rit)->getEnt());
1053 int nb_dofs_on_ent = std::distance(range.first, range.second);
1054
1055 int max_order = (*rit)->getMaxOrder();
1056 if ((*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order) !=
1057 nb_dofs_on_ent) {
1058
1059
1060
1061
1062
1064 << "Warning: Number of Dofs in Row different than number "
1065 "of dofs for given entity order "
1066 << (*rit)->getNbOfCoeffs() * (*rit)->getOrderNbDofs(max_order)
1067 << " " << nb_dofs_on_ent;
1068 }
1069 }
1070 }
1071
1073 }
1074
1077
1078
1079 CHKERR MatAssemblyBegin(
A, MAT_FLUSH_ASSEMBLY);
1080 CHKERR MatAssemblyEnd(
A, MAT_FLUSH_ASSEMBLY);
1081
1083 }
1084 };
1085
1086
1087 CHKERR MatSetOption(
A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE);
1088
1090 MatView(
A, PETSC_VIEWER_STDOUT_WORLD);
1091 }
1092
1093 if (verb >=
NOISY) {
1094 MatView(
A, PETSC_VIEWER_DRAW_WORLD);
1095 std::string wait;
1096 std::cin >> wait;
1097 }
1098
1099 TestMatrixFillIn method(&m_field,
A, row_print, col_print);
1100
1101
1103 auto &prb_set = problems_ptr->get<Problem_mi_tag>();
1104 auto p_miit = prb_set.find(problem_name);
1105 if (p_miit == prb_set.end())
1107 "problem < %s > not found (top tip: check spelling)",
1108 problem_name.c_str());
1109 MOFEM_LOG_C(
"WORLD", Sev::inform,
"check problem < %s >",
1110 problem_name.c_str());
1111
1112
1114 for (auto &fe : *fe_ptr) {
1115 MOFEM_LOG_C(
"WORLD", Sev::verbose,
"\tcheck element %s",
1116 fe->getName().c_str());
1120 }
1121
1122 CHKERR MatAssemblyBegin(
A, MAT_FINAL_ASSEMBLY);
1123 CHKERR MatAssemblyEnd(
A, MAT_FINAL_ASSEMBLY);
1124
1126
1128}
#define MOFEM_LOG_C(channel, severity, format,...)
#define MatrixManagerFunctionBegin
Create adjacent matrices using different indices.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
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.
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
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.
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
virtual MPI_Comm & get_comm() const =0
Public Member Functions inherited from MoFEM::UnknownInterface