From e0fde3a2dca6c4db337996e2ea8f7d4c470d34e4 Mon Sep 17 00:00:00 2001 From: Thomas Witkowski <thomas.witkowski@gmx.de> Date: Mon, 31 Aug 2009 14:28:23 +0000 Subject: [PATCH] Changes in ProblemVec matrix debugging functions. --- AMDiS/src/ProblemVec.cc | 74 ++++++++++++++++++++++++----------------- 1 file changed, 44 insertions(+), 30 deletions(-) diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc index 49478ab5..7fe5b59a 100644 --- a/AMDiS/src/ProblemVec.cc +++ b/AMDiS/src/ProblemVec.cc @@ -660,8 +660,6 @@ namespace AMDiS { // will be set to false). DOFMatrix *matrix = (*systemMatrix)[i][j]; - std::cout << " i = " << i << " j = " << j << " matrix " << matrix << std::endl; - if (matrix) matrix->calculateNnz(); @@ -765,29 +763,6 @@ namespace AMDiS { if (matrix) nnz += matrix->getBaseMatrix().nnz(); - - std::cout << " i = " << i << " j = " << j << " matrix " << matrix << std::endl; - - if (i == 0 && matrix) { - using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end; - namespace traits= mtl::traits; - typedef DOFMatrix::base_matrix_type Matrix; - - traits::row<Matrix>::type row(matrix->getBaseMatrix()); - traits::col<Matrix>::type col(matrix->getBaseMatrix()); - traits::const_value<Matrix>::type value(matrix->getBaseMatrix()); - - typedef traits::range_generator<major, Matrix>::type cursor_type; - typedef traits::range_generator<nz, cursor_type>::type icursor_type; - - std::cout.precision(10); - for (cursor_type cursor = begin<major>(matrix->getBaseMatrix()), cend = end<major>(matrix->getBaseMatrix()); cursor != cend; ++cursor) - for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) - if (row(*icursor) == 0) - std::cout << "C = " << col(*icursor) << " : " << value(*icursor) << std::endl; - std::cout << "============================" << std::endl; - } - } // And now assemble boundary conditions on the vectors @@ -1464,19 +1439,27 @@ namespace AMDiS { typedef traits::range_generator<major, Matrix>::type cursor_type; typedef traits::range_generator<nz, cursor_type>::type icursor_type; + // Create map with all world coordinates of all DOFs. DofToCoord dofToCoords; createDofToCoordMap(dofToCoords); + // Start writing output file. std::ofstream out(filename.c_str()); + out.precision(15); + + // First, we write the number of DOFs the file contains. out << dofToCoords.size() << std::endl; int dimOfWorld = Global::getGeo(WORLD); + // And now write all the DOF's coords. The first line contains the dof index, all + // the other lines contain one component of the coordinates. for (DofToCoord::iterator it = dofToCoords.begin(); it != dofToCoords.end(); ++it) { out << it->first << std::endl; for (int j = 0; j < dimOfWorld; j++) out << (it->second)[j] << std::endl; } + // Write the matrices. for (int i = 0; i < nComponents; i++) { for (int j = 0; j < nComponents; j++) { DOFMatrix *dofmatrix = (*systemMatrix)[i][j]; @@ -1488,6 +1471,7 @@ namespace AMDiS { int nnz = matrix.nnz(); int testNnz = 0; + // Write to file, how many entries the matrix conatins. out << nnz << std::endl; traits::row<Matrix>::type row(matrix); @@ -1503,6 +1487,7 @@ namespace AMDiS { } } + // Just to check, if nnz() is correct. TEST_EXIT(testNnz == nnz)("NNZ does not fit together!\n"); } } @@ -1519,6 +1504,7 @@ namespace AMDiS { typedef traits::range_generator<major, Matrix>::type cursor_type; typedef traits::range_generator<nz, cursor_type>::type icursor_type; + // Create a map from coords of all DOFs, to the DOF indices in this problem. CoordToDof coordToDof; createCoordToDofMap(coordToDof); @@ -1527,9 +1513,12 @@ namespace AMDiS { DegreeOfFreedom dof; int dimOfWorld = Global::getGeo(WORLD); + // Open file and read the number of DOFs the file contains. std::ifstream in(filename.c_str()); int nReadDofs; in >> nReadDofs; + + // Read all DOF indices and their world coordinates. for (int i = 0; i < nReadDofs; i++) { in >> dof; for (int j = 0; j < dimOfWorld; j++) @@ -1537,6 +1526,7 @@ namespace AMDiS { dofToCoord[dof] = coords; } + // Now we traverse all matrices and check them. for (int i = 0; i < nComponents; i++) { for (int j = 0; j < nComponents; j++) { DOFMatrix *dofmatrix = (*systemMatrix)[i][j]; @@ -1547,6 +1537,8 @@ namespace AMDiS { int readNnz; in >> readNnz; + // Read each entry in file and check it with the corresponding matrix value + // in this problem. for (int k = 0; k < readNnz; k++) { DegreeOfFreedom row, col; double value; @@ -1558,12 +1550,34 @@ namespace AMDiS { WorldVector<double> rowCoords = dofToCoord[row]; WorldVector<double> colCoords = dofToCoord[col]; - if (coordToDof.find(rowCoords) == coordToDof.end()) - std::cout << "CANNOT FIND" << std::endl; + if (coordToDof.find(rowCoords) == coordToDof.end()) { + std::cout << "CANNOT ROW FIND" << std::endl; + rowCoords.print(); + exit(0); + } + + if (coordToDof.find(colCoords) == coordToDof.end()) { + std::cout << "CANNOT COL FIND" << std::endl; + rowCoords.print(); + exit(0); + } + + DegreeOfFreedom rowHere = coordToDof[rowCoords]; + DegreeOfFreedom colHere = coordToDof[colCoords]; + double valueHere = (dofmatrix->getBaseMatrix())[rowHere][colHere]; + double diff = fabs(valueHere - value); + + // And so we check the difference between the value read from file and the + // corresponding value in the matrix. + if (diff > 1e-10) { + std::cout << "WRONG value in matrix[" << i << "][" << j << "]!" << std::endl + << " DOF in other matrix [" << row << "][" << col << "] = " << value << std::endl + << " DOF in this matrix [" << rowHere << "][" << colHere << "] = " << valueHere << std::endl; + + + exit(0); + } - if (coordToDof.find(colCoords) == coordToDof.end()) - std::cout << "CANNOT FIND" << std::endl; - } } } -- GitLab