diff --git a/AMDiS/bin/Makefile.am b/AMDiS/bin/Makefile.am
index bae4bb4c9d954dd943215db50b51248a1aa6b8ea..943384e385e1effa353f678b23b536a8fa951391 100644
--- a/AMDiS/bin/Makefile.am
+++ b/AMDiS/bin/Makefile.am
@@ -82,6 +82,7 @@ $(SOURCE_DIR)/ProblemIterationInterface.h \
$(SOURCE_DIR)/StandardProblemIteration.h $(SOURCE_DIR)/StandardProblemIteration.cc \
$(SOURCE_DIR)/ProblemScal.h $(SOURCE_DIR)/ProblemScal.cc \
$(SOURCE_DIR)/ProblemVec.h $(SOURCE_DIR)/ProblemVec.cc \
+$(SOURCE_DIR)/ProblemVecDbg.h $(SOURCE_DIR)/ProblemVecDbg.cc \
$(SOURCE_DIR)/DualTraverse.h $(SOURCE_DIR)/DualTraverse.cc \
$(SOURCE_DIR)/ElementPartition_ED.h $(SOURCE_DIR)/SurfacePartition_ED.h \
$(SOURCE_DIR)/ElementData.h $(SOURCE_DIR)/ElementData.cc \
diff --git a/AMDiS/bin/Makefile.in b/AMDiS/bin/Makefile.in
index 06d1aa21ef859dc758ab347b158a971f0a2a4335..ccfca10ca623980bdbc15b01afe65a09830cea9d 100644
--- a/AMDiS/bin/Makefile.in
+++ b/AMDiS/bin/Makefile.in
@@ -109,6 +109,7 @@ am__libamdis_la_SOURCES_DIST = $(PARALLEL_DIR)/ParallelDomainBase.h \
$(SOURCE_DIR)/StandardProblemIteration.cc \
$(SOURCE_DIR)/ProblemScal.h $(SOURCE_DIR)/ProblemScal.cc \
$(SOURCE_DIR)/ProblemVec.h $(SOURCE_DIR)/ProblemVec.cc \
+ $(SOURCE_DIR)/ProblemVecDbg.h $(SOURCE_DIR)/ProblemVecDbg.cc \
$(SOURCE_DIR)/DualTraverse.h $(SOURCE_DIR)/DualTraverse.cc \
$(SOURCE_DIR)/ElementPartition_ED.h \
$(SOURCE_DIR)/SurfacePartition_ED.h \
@@ -248,7 +249,8 @@ am_libamdis_la_OBJECTS = $(am__objects_2) libamdis_la-DOFIndexed.lo \
libamdis_la-AdaptBase.lo \
libamdis_la-StandardProblemIteration.lo \
libamdis_la-ProblemScal.lo libamdis_la-ProblemVec.lo \
- libamdis_la-DualTraverse.lo libamdis_la-ElementData.lo \
+ libamdis_la-ProblemVecDbg.lo libamdis_la-DualTraverse.lo \
+ libamdis_la-ElementData.lo \
libamdis_la-ComponentTraverseInfo.lo libamdis_la-CreatorMap.lo \
libamdis_la-ProblemInterpolScal.lo \
libamdis_la-ProblemInterpolVec.lo libamdis_la-MacroReader.lo \
@@ -493,6 +495,7 @@ $(SOURCE_DIR)/ProblemIterationInterface.h \
$(SOURCE_DIR)/StandardProblemIteration.h $(SOURCE_DIR)/StandardProblemIteration.cc \
$(SOURCE_DIR)/ProblemScal.h $(SOURCE_DIR)/ProblemScal.cc \
$(SOURCE_DIR)/ProblemVec.h $(SOURCE_DIR)/ProblemVec.cc \
+$(SOURCE_DIR)/ProblemVecDbg.h $(SOURCE_DIR)/ProblemVecDbg.cc \
$(SOURCE_DIR)/DualTraverse.h $(SOURCE_DIR)/DualTraverse.cc \
$(SOURCE_DIR)/ElementPartition_ED.h $(SOURCE_DIR)/SurfacePartition_ED.h \
$(SOURCE_DIR)/ElementData.h $(SOURCE_DIR)/ElementData.cc \
@@ -777,6 +780,7 @@ distclean-compile:
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-ProblemNonLin.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-ProblemScal.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-ProblemVec.Plo@am__quote@
+@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-ProblemVecDbg.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-Projection.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-QPsiPhi.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libamdis_la-Quadrature.Plo@am__quote@
@@ -969,6 +973,13 @@ libamdis_la-ProblemVec.lo: $(SOURCE_DIR)/ProblemVec.cc
@AMDEP_TRUE@@am__fastdepCXX_FALSE@ DEPDIR=$(DEPDIR) $(CXXDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCXX_FALSE@ $(LIBTOOL) --tag=CXX --mode=compile $(CXX) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(libamdis_la_CXXFLAGS) $(CXXFLAGS) -c -o libamdis_la-ProblemVec.lo `test -f '$(SOURCE_DIR)/ProblemVec.cc' || echo '$(srcdir)/'`$(SOURCE_DIR)/ProblemVec.cc
+libamdis_la-ProblemVecDbg.lo: $(SOURCE_DIR)/ProblemVecDbg.cc
+@am__fastdepCXX_TRUE@ if $(LIBTOOL) --tag=CXX --mode=compile $(CXX) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(libamdis_la_CXXFLAGS) $(CXXFLAGS) -MT libamdis_la-ProblemVecDbg.lo -MD -MP -MF "$(DEPDIR)/libamdis_la-ProblemVecDbg.Tpo" -c -o libamdis_la-ProblemVecDbg.lo `test -f '$(SOURCE_DIR)/ProblemVecDbg.cc' || echo '$(srcdir)/'`$(SOURCE_DIR)/ProblemVecDbg.cc; \
+@am__fastdepCXX_TRUE@ then mv -f "$(DEPDIR)/libamdis_la-ProblemVecDbg.Tpo" "$(DEPDIR)/libamdis_la-ProblemVecDbg.Plo"; else rm -f "$(DEPDIR)/libamdis_la-ProblemVecDbg.Tpo"; exit 1; fi
+@AMDEP_TRUE@@am__fastdepCXX_FALSE@ source='$(SOURCE_DIR)/ProblemVecDbg.cc' object='libamdis_la-ProblemVecDbg.lo' libtool=yes @AMDEPBACKSLASH@
+@AMDEP_TRUE@@am__fastdepCXX_FALSE@ DEPDIR=$(DEPDIR) $(CXXDEPMODE) $(depcomp) @AMDEPBACKSLASH@
+@am__fastdepCXX_FALSE@ $(LIBTOOL) --tag=CXX --mode=compile $(CXX) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(libamdis_la_CXXFLAGS) $(CXXFLAGS) -c -o libamdis_la-ProblemVecDbg.lo `test -f '$(SOURCE_DIR)/ProblemVecDbg.cc' || echo '$(srcdir)/'`$(SOURCE_DIR)/ProblemVecDbg.cc
+
libamdis_la-DualTraverse.lo: $(SOURCE_DIR)/DualTraverse.cc
@am__fastdepCXX_TRUE@ if $(LIBTOOL) --tag=CXX --mode=compile $(CXX) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(libamdis_la_CXXFLAGS) $(CXXFLAGS) -MT libamdis_la-DualTraverse.lo -MD -MP -MF "$(DEPDIR)/libamdis_la-DualTraverse.Tpo" -c -o libamdis_la-DualTraverse.lo `test -f '$(SOURCE_DIR)/DualTraverse.cc' || echo '$(srcdir)/'`$(SOURCE_DIR)/DualTraverse.cc; \
@am__fastdepCXX_TRUE@ then mv -f "$(DEPDIR)/libamdis_la-DualTraverse.Tpo" "$(DEPDIR)/libamdis_la-DualTraverse.Plo"; else rm -f "$(DEPDIR)/libamdis_la-DualTraverse.Tpo"; exit 1; fi
diff --git a/AMDiS/libtool b/AMDiS/libtool
index e5172c05a262fdd3d3a02b791a0fdd52443eb55b..5984c8db56b2733b2919e078dda5b71f4a5a1124 100755
--- a/AMDiS/libtool
+++ b/AMDiS/libtool
@@ -44,7 +44,7 @@ available_tags=" CXX F77"
# ### BEGIN LIBTOOL CONFIG
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p2d125:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
@@ -6760,7 +6760,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac`
# End:
# ### BEGIN LIBTOOL TAG CONFIG: CXX
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p2d125:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
@@ -7065,7 +7065,7 @@ include_expsyms=""
# ### BEGIN LIBTOOL TAG CONFIG: F77
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p2d125:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
diff --git a/AMDiS/src/ParallelDomainBase.cc b/AMDiS/src/ParallelDomainBase.cc
index bb0f6930643375d654f434bd4023cc2010a18790..b7964d6edc3d391fa8344ed211a5d403e6584cf7 100644
--- a/AMDiS/src/ParallelDomainBase.cc
+++ b/AMDiS/src/ParallelDomainBase.cc
@@ -16,7 +16,6 @@
#include "ProblemStatBase.h"
#include "StandardProblemIteration.h"
#include "ElementFileWriter.h"
-#include "Serializer.h"
#include "petscksp.h"
@@ -35,13 +34,13 @@ namespace AMDiS {
return (*dof1 < *dof2);
}
- ParallelDomainBase::ParallelDomainBase(std::string name,
- ProblemIterationInterface *iIF,
+ ParallelDomainBase::ParallelDomainBase(ProblemIterationInterface *iIF,
ProblemTimeInterface *tIF,
FiniteElemSpace *fe,
RefinementManager *refineManager)
: iterationIF(iIF),
timeIF(tIF),
+ name(iIF->getName()),
feSpace(fe),
mesh(fe->getMesh()),
refinementManager(refineManager),
@@ -1605,7 +1604,7 @@ namespace AMDiS {
SerUtil::serialize(out, mapLocalToDofIndex);
SerUtil::serialize(out, isRankDof);
- // SerUtil::serialize(out, vertexDof);
+ serialize(out, vertexDof);
SerUtil::serialize(out, &rstart);
SerUtil::serialize(out, &nRankRows);
@@ -1649,7 +1648,7 @@ namespace AMDiS {
SerUtil::deserialize(in, mapLocalToDofIndex);
SerUtil::deserialize(in, isRankDof);
- // SerUtil::deserialize(in, vertexDof);
+ deserialize(in, vertexDof, dofMap);
SerUtil::deserialize(in, &rstart);
SerUtil::deserialize(in, &nRankRows);
diff --git a/AMDiS/src/ParallelDomainBase.h b/AMDiS/src/ParallelDomainBase.h
index 60bfd19e1b4182ad1c5fcb5f7bc60a822b875fcc..88ea2565cd61971d373cceaffccfb258804de483 100644
--- a/AMDiS/src/ParallelDomainBase.h
+++ b/AMDiS/src/ParallelDomainBase.h
@@ -32,6 +32,7 @@
#include "FiniteElemSpace.h"
#include "AdaptInfo.h"
#include "InteriorBoundary.h"
+#include "Serializer.h"
#include "AMDiS_fwd.h"
#include "petsc.h"
@@ -81,8 +82,7 @@ namespace AMDiS {
typedef std::map<const DegreeOfFreedom*, DegreeOfFreedom> DofIndexMap;
public:
- ParallelDomainBase(std::string name,
- ProblemIterationInterface *iterationIF,
+ ParallelDomainBase(ProblemIterationInterface *iterationIF,
ProblemTimeInterface *timeIF,
FiniteElemSpace *feSpace,
RefinementManager *refineManager);
@@ -305,7 +305,38 @@ namespace AMDiS {
/// Reads a \ref RankToDofContainer from an input stream.
void deserialize(std::istream &in, RankToDofContainer &data,
std::map<int, const DegreeOfFreedom*> &dofMap);
-
+
+ /// Writes a mapping from dof pointers to some values to an output stream.
+ template<typename T>
+ void serialize(std::ostream &out, std::map<const DegreeOfFreedom*, T> &data)
+ {
+ int mapSize = data.size();
+ SerUtil::serialize(out, &mapSize);
+ for (typename std::map<const DegreeOfFreedom*, T>::iterator it = data.begin();
+ it != data.end(); ++it) {
+ int v1 = (*(it->first));
+ T v2 = it->second;
+ SerUtil::serialize(out, &v1);
+ SerUtil::serialize(out, &v2);
+ }
+ }
+
+ /// Reads a mapping from dof pointer to some values from an input stream.
+ template<typename T>
+ void deserialize(std::istream &in, std::map<const DegreeOfFreedom*, T> &data,
+ std::map<int, const DegreeOfFreedom*> &dofMap)
+ {
+ int mapSize = 0;
+ SerUtil::deserialize(in, &mapSize);
+ for (int i = 0; i < mapSize; i++) {
+ int v1 = 0;
+ T v2;
+ SerUtil::deserialize(in, &v1);
+ SerUtil::deserialize(in, &v2);
+ data[dofMap[v1]] = v2;
+ }
+ }
+
inline void orderDOFs(const DegreeOfFreedom* dof1,
const DegreeOfFreedom* dof2,
const DegreeOfFreedom* dof3,
diff --git a/AMDiS/src/ParallelDomainScal.cc b/AMDiS/src/ParallelDomainScal.cc
index a9f877556e8bde63f29e0d16fc57e8c6fa2a4f4b..4335aef95976fa835bfafa72c28e19412f5dc139 100644
--- a/AMDiS/src/ParallelDomainScal.cc
+++ b/AMDiS/src/ParallelDomainScal.cc
@@ -5,11 +5,9 @@
namespace AMDiS {
- ParallelDomainScal::ParallelDomainScal(std::string name,
- ProblemScal *problem,
+ ParallelDomainScal::ParallelDomainScal(ProblemScal *problem,
ProblemInstatScal *problemInstat)
- : ParallelDomainBase(name,
- problem,
+ : ParallelDomainBase(problem,
problemInstat,
problem->getFESpace(),
problem->getRefinementManager()),
diff --git a/AMDiS/src/ParallelDomainScal.h b/AMDiS/src/ParallelDomainScal.h
index 226e32a8f67b63c26f30207397ec3802c0916a40..99f571a3155c0f8e295aedf92af565abf123052f 100644
--- a/AMDiS/src/ParallelDomainScal.h
+++ b/AMDiS/src/ParallelDomainScal.h
@@ -29,8 +29,7 @@ namespace AMDiS {
class ParallelDomainScal : public ParallelDomainBase
{
public:
- ParallelDomainScal(std::string name,
- ProblemScal *problem,
+ ParallelDomainScal(ProblemScal *problem,
ProblemInstatScal *problemInstat);
void initParallelization(AdaptInfo *adaptInfo);
diff --git a/AMDiS/src/ParallelDomainVec.cc b/AMDiS/src/ParallelDomainVec.cc
index cd5d47ca7bb490a9089d87c8274781d643d1dd0f..c3ff8b5cdd2a16a507925de7119a64b1c4d7c2e3 100644
--- a/AMDiS/src/ParallelDomainVec.cc
+++ b/AMDiS/src/ParallelDomainVec.cc
@@ -5,11 +5,9 @@
namespace AMDiS {
- ParallelDomainVec::ParallelDomainVec(std::string name,
- ProblemVec *problem,
+ ParallelDomainVec::ParallelDomainVec(ProblemVec *problem,
ProblemInstatVec *problemInstat)
- : ParallelDomainBase(name,
- problem,
+ : ParallelDomainBase(problem,
problemInstat,
problem->getFESpace(0),
problem->getRefinementManager(0)),
diff --git a/AMDiS/src/ParallelDomainVec.h b/AMDiS/src/ParallelDomainVec.h
index 0d81d547bffc7d9c64c101945839f93cc2025964..bea1cf9ab12a737f564c3fce97130feee15c805c 100644
--- a/AMDiS/src/ParallelDomainVec.h
+++ b/AMDiS/src/ParallelDomainVec.h
@@ -30,8 +30,7 @@ namespace AMDiS {
class ParallelDomainVec : public ParallelDomainBase
{
public:
- ParallelDomainVec(std::string name,
- ProblemVec *problem,
+ ParallelDomainVec(ProblemVec *problem,
ProblemInstatVec *problemInstat);
void initParallelization(AdaptInfo *adaptInfo);
diff --git a/AMDiS/src/ProblemScal.cc b/AMDiS/src/ProblemScal.cc
index c8cf73031794db705a5a07bb3976433ff37128ba..da30f3541a2474d25b8f4742c94591565fe2c000 100644
--- a/AMDiS/src/ProblemScal.cc
+++ b/AMDiS/src/ProblemScal.cc
@@ -485,10 +485,14 @@ namespace AMDiS {
{
fileWriters.push_back(new FileWriter(name + "->output", mesh, solution));
int writeSerialization = 0;
- GET_PARAMETER(0, name + "->output->write serialization", "%d",
- &writeSerialization);
+ GET_PARAMETER(0, name + "->output->write serialization", "%d", &writeSerialization);
+
+ // Serialization is not allowed to be done by the problem, if its part of a parallel
+ // problem definition. Than, the parallel problem object must be serialized.
+#ifndef HAVE_PARALLEL_DOMAIN_AMDIS
if (writeSerialization)
fileWriters.push_back(new Serializer<ProblemScal>(this));
+#endif
}
void ProblemScal::estimate(AdaptInfo *adaptInfo)
diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc
index c3474cfbfe463f74c35b3d5c1c36f52a2a301433..3cee67a6f99abb65cb293519c9cfb2f15429eb47 100644
--- a/AMDiS/src/ProblemVec.cc
+++ b/AMDiS/src/ProblemVec.cc
@@ -470,8 +470,13 @@ namespace AMDiS {
}
GET_PARAMETER(0, name + "->output->write serialization", "%d", &writeSerialization);
+
+ // Serialization is not allowed to be done by the problem, if its part of a parallel
+ // problem definition. Than, the parallel problem object must be serialized.
+#ifndef HAVE_PARALLEL_DOMAIN_AMDIS
if (writeSerialization)
fileWriters.push_back(new Serializer<ProblemVec>(this));
+#endif
}
void ProblemVec::doOtherStuff()
@@ -1385,404 +1390,6 @@ namespace AMDiS {
delete sol;
}
}
-
-
- void ProblemVec::createDofToCoordMap(DofToCoord &dofMap)
- {
- const BasisFunction *basisFcts = componentSpaces[0]->getBasisFcts();
- WorldVector<double> coords;
- ElementDofIterator elDofIter(componentSpaces[0], true);
-
- TraverseStack stack;
- ElInfo *elInfo = stack.traverseFirst(componentMeshes[0], -1,
- Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
- while (elInfo) {
- elDofIter.reset(elInfo->getElement());
- int i = 0;
- do {
- DimVec<double> *baryCoord = basisFcts->getCoords(i);
- elInfo->coordToWorld(*baryCoord, coords);
- i++;
- dofMap[elDofIter.getDof()] = coords;
- } while (elDofIter.next());
-
- elInfo = stack.traverseNext(elInfo);
- }
- }
-
-
- void ProblemVec::createCoordToDofMap(CoordToDof &dofMap)
- {
- const BasisFunction *basisFcts = componentSpaces[0]->getBasisFcts();
- WorldVector<double> coords;
- ElementDofIterator elDofIter(componentSpaces[0], true);
-
- TraverseStack stack;
- ElInfo *elInfo = stack.traverseFirst(componentMeshes[0], -1,
- Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
- while (elInfo) {
- elDofIter.reset(elInfo->getElement());
- int i = 0;
- do {
- DimVec<double> *baryCoord = basisFcts->getCoords(i);
- elInfo->coordToWorld(*baryCoord, coords);
- i++;
- dofMap[coords] = elDofIter.getDof();
- } while (elDofIter.next());
-
- elInfo = stack.traverseNext(elInfo);
- }
- }
-
-
- void ProblemVec::writeDbgMatrix(std::string filename)
- {
- FUNCNAME("ProblemVec::writeDbMatrix()");
-
- using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end;
- namespace traits = mtl::traits;
- typedef DOFMatrix::base_matrix_type Matrix;
- 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];
-
- if (!dofmatrix)
- continue;
-
- Matrix& matrix = dofmatrix->getBaseMatrix();
- 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);
- traits::col<Matrix>::type col(matrix);
- traits::const_value<Matrix>::type value(matrix);
-
- for (cursor_type cursor = begin<major>(matrix), cend = end<major>(matrix); cursor != cend; ++cursor) {
- for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) {
- out << row(*icursor) << " "
- << col(*icursor) << " "
- << value(*icursor) << std::endl;
- testNnz++;
- }
- }
-
- // Just to check, if nnz() is correct.
- TEST_EXIT(testNnz == nnz)("NNZ does not fit together!\n");
- }
- }
-
-
- for (int i = 0; i < nComponents; i++) {
- out << rhs->getDOFVector(i)->getUsedSize() << std::endl;
- DOFIterator<double> it(rhs->getDOFVector(i), USED_DOFS);
- int c = 0;
- for (it.reset(); !it.end(); ++it) {
- out << c << " " << *it << std::endl;
- c++;
- }
-
- TEST_EXIT(c == rhs->getDOFVector(i)->getUsedSize())
- ("RHS NNZ does not fit together!\n");
- }
-
-
- for (int i = 0; i < nComponents; i++) {
- out << solution->getDOFVector(i)->getUsedSize() << std::endl;
- DOFIterator<double> it(solution->getDOFVector(i), USED_DOFS);
- int c = 0;
- for (it.reset(); !it.end(); ++it) {
- out << c << " " << *it << std::endl;
- c++;
- }
-
- TEST_EXIT(c == rhs->getDOFVector(i)->getUsedSize())
- ("RHS NNZ does not fit together!\n");
- }
-
- out.close();
- }
-
-
- void ProblemVec::readAndCompareDbgMatrix(std::vector<std::string> filenames)
- {
- using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end;
- namespace traits = mtl::traits;
- typedef DOFMatrix::base_matrix_type Matrix;
- 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);
-
- int dimOfWorld = Global::getGeo(WORLD);
- std::vector<std::vector<std::map<std::pair<int, int>, double> > > nnzValues(nComponents);
- std::vector<std::map<int, double> > rhsValues(nComponents);
- std::vector<std::map<int, double> > solValues(nComponents);
- for (int i = 0; i < nComponents; i++)
- nnzValues[i].resize(nComponents);
-
- // Stores to each dof index of this problem a map from rank indices (of each rank
- // that also has this dof) to the corresponding local dof index.
- std::map<DegreeOfFreedom, std::vector<DegreeOfFreedom> > dofMapHereToFiles;
-
- for (std::vector<std::string>::iterator fileIt = filenames.begin();
- fileIt != filenames.end(); ++fileIt) {
-
- // Open file and read the number of DOFs the file contains.
- std::ifstream in(fileIt->c_str());
- int nReadDofs;
- in >> nReadDofs;
-
- // Is used to map the dof indices in the files to the global coordinates.
- DofToCoord dofToCoord;
-
- // Read all DOF indices and their world coordinates.
- for (int i = 0; i < nReadDofs; i++) {
- DegreeOfFreedom dof;
- WorldVector<double> coords;
-
- in >> dof;
- for (int j = 0; j < dimOfWorld; j++)
- in >> coords[j];
- dofToCoord[dof] = coords;
- }
-
- std::map<DegreeOfFreedom, DegreeOfFreedom> dofMapFileToHere;
-
- for (DofToCoord::iterator it = dofToCoord.begin(); it != dofToCoord.end(); ++it) {
- DegreeOfFreedom dofIndexInFile = it->first;
- WorldVector<double> &dofCoords = it->second;
-
- if (coordToDof.find(dofCoords) == coordToDof.end()) {
- std::cout << "Cannot find dof index for coords: " << std::endl;
- dofCoords.print();
- exit(0);
- }
-
- DegreeOfFreedom dofIndexHere = coordToDof[dofCoords];
- dofMapHereToFiles[dofIndexHere].push_back(dofIndexInFile);
- dofMapFileToHere[dofIndexInFile] = dofIndexHere;
- }
-
- // 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];
-
- if (!dofmatrix)
- continue;
-
- 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;
-
- in >> row;
- in >> col;
- in >> value;
-
- if (dofMapFileToHere.count(row) == 0) {
- std::cout << "Cannot find row index for: " << row << std::endl;
- exit(0);
- }
-
- if (dofMapFileToHere.count(col) == 0) {
- std::cout << "Cannot find col index for: " << col << std::endl;
- exit(0);
- }
-
- // Get dof indices for row and col of this problem matrix.
- DegreeOfFreedom rowHere = dofMapFileToHere[row];
- DegreeOfFreedom colHere = dofMapFileToHere[col];
- std::pair<int, int> rowcol = make_pair(rowHere, colHere);
-
- if (nnzValues[i][j].count(rowcol) == 0)
- nnzValues[i][j][rowcol] = value;
- else
- nnzValues[i][j][rowcol] += value;
- }
- }
- }
-
- for (int i = 0; i < nComponents; i++) {
- int readNnz;
- in >> readNnz;
- for (int k = 0; k < readNnz; k++) {
- int row;
- double value;
-
- in >> row;
- in >> value;
-
- WorldVector<double> rowCoords = dofToCoord[row];
-
- if (coordToDof.find(rowCoords) == coordToDof.end()) {
- std::cout << "Cannot find row index for coords: " << std::endl;
- rowCoords.print();
- exit(0);
- }
-
- DegreeOfFreedom rowHere = coordToDof[rowCoords];
-
- if (rhsValues[i].count(rowHere) == 0)
- rhsValues[i][rowHere] = value;
- else
- rhsValues[i][rowHere] += value;
- }
- }
-
- for (int i = 0; i < nComponents; i++) {
- int readNnz;
- in >> readNnz;
- for (int k = 0; k < readNnz; k++) {
- int row;
- double value;
-
- in >> row;
- in >> value;
-
- WorldVector<double> rowCoords = dofToCoord[row];
-
- if (coordToDof.find(rowCoords) == coordToDof.end()) {
- std::cout << "Cannot find row index for coords: " << std::endl;
- rowCoords.print();
- exit(0);
- }
-
- DegreeOfFreedom rowHere = coordToDof[rowCoords];
-
- if (solValues[i].count(rowHere) == 0) {
- solValues[i][rowHere] = value;
- } else {
- double diff = fabs(solValues[i][rowHere] - value);
- if (diff > 1e-8) {
- std::cout << "DIFFERENT values in solution vector!" << std::endl;
- exit(0);
- }
- }
- }
- }
-
- in.close();
-
- std::cout << "File read!" << std::endl;
- }
-
- std::cout << "Now checking ..." << std::endl;
-
- std::cout.precision(15);
-
- for (int i = 0; i < nComponents; i++) {
- for (int j = 0; j < nComponents; j++) {
- DOFMatrix *dofmatrix = (*systemMatrix)[i][j];
-
- if (!dofmatrix)
- continue;
-
- for (std::map<std::pair<int, int>, double>::iterator nnzIt =
- nnzValues[i][j].begin(); nnzIt != nnzValues[i][j].end(); ++nnzIt) {
-
- int row = nnzIt->first.first;
- int col = nnzIt->first.second;
- double value = nnzIt->second;
-
- double valueHere = (dofmatrix->getBaseMatrix())[row][col];
- 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-8) {
- std::cout << "Wrong value in component matrix " << i
- << "/" << j << ": " << std::endl
- << " DOF in this matrix[" << row << "][" << col << "] = "
- << valueHere << std::endl
- << " DOF in other matrix[" << dofMapHereToFiles[row][0] << "]["
- << dofMapHereToFiles[col][0] << "] = "
- << value << std::endl;
-
- exit(0);
- }
- }
- }
- }
-
- for (int i = 0; i < nComponents; i++) {
- for (std::map<int, double>::iterator rhsIt = rhsValues[i].begin();
- rhsIt != rhsValues[i].end(); ++rhsIt) {
-
- int row = rhsIt->first;
- double value = rhsIt->second;
-
- double valueHere = (*(rhs->getDOFVector(i)))[row];
- double diff = fabs(valueHere - value);
-
- if (diff > 1e-8) {
- std::cout << "WRONG value in rhs[" << i << "]!" << std::endl
- << " DOF in other rhs[" << row << "] = " << value << std::endl
- << " DOF in this rhs[" << row << "] = " << valueHere << std::endl;
-
- exit(0);
- }
- }
- }
-
- for (int i = 0; i < nComponents; i++) {
- for (std::map<int, double>::iterator solIt = solValues[i].begin();
- solIt != solValues[i].end(); ++solIt) {
-
- int row = solIt->first;
- double value = solIt->second;
-
- double valueHere = (*(solution->getDOFVector(i)))[row];
- double diff = fabs(valueHere - value);
-
- if (diff > 1e-8) {
- std::cout << "WRONG value in sol[" << i << "]!" << std::endl
- << " DOF in other sol[" << row << "] = " << value << std::endl
- << " DOF in this sol[" << row << "] = " << valueHere << std::endl;
-
- exit(0);
- }
- }
- }
-
- std::cout << "FINISHED COMPARING!" << std::endl;
-
- exit(0);
- }
-
}
diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h
index 3445d0b3f2da8d58e84cf94c026942398f4133ca..451d8b111dd7e2b6dc6402dfa47bf70b7af470b2 100644
--- a/AMDiS/src/ProblemVec.h
+++ b/AMDiS/src/ProblemVec.h
@@ -40,7 +40,7 @@ namespace AMDiS {
class ProblemVec : public ProblemStatBase,
public StandardProblemIteration
{
- private:
+ protected:
// Defines a mapping type from dof indices to world coordinates.
typedef std::map<DegreeOfFreedom, WorldVector<double> > DofToCoord;
@@ -473,21 +473,6 @@ namespace AMDiS {
return fileWriters;
}
- /** \brief
- * This function writes the system matrix and the system vector to a file
- * for debugging. The entries of both, the matrix and the vector, are not
- * indicated by the dof indices, but by the world coordinates of the dofs.
- * This makes it possible to compare dof matrices, where the dof indices
- * are different, e.g., when using domain decomposition parallelization.
- */
- void writeDbgMatrix(std::string filename);
-
- /** \brief
- * Reads a file, that was created using the function \ref writeDbgMatrix, and
- * compares the date with the system matrix of this problem.
- */
- void readAndCompareDbgMatrix(std::vector<std::string> filenames);
-
protected:
/** \brief
* If the exact solution is known, the problem can compute the exact
@@ -495,20 +480,6 @@ namespace AMDiS {
*/
void computeError(AdaptInfo *adaptInfo);
- /** \brief
- * Create from the current problem a map from dof indices to world coordinates.
- * This function is used by the debugging function \ref writeDbgMatrix and \ref
- * readAndCompareDbgMatrix.
- */
- void createDofToCoordMap(DofToCoord &dofMap);
-
- /** \brief
- * Create from the current problem a map from world coordinates of all dofs to
- * the corresponding dof indices. This function is used by the debugging function
- * \ref writeDbgMatrix and \ref readAndCompareDbgMatrix.
- */
- void createCoordToDofMap(CoordToDof &dofMap);
-
protected:
/// Name of this problem.
diff --git a/AMDiS/src/ProblemVecDbg.cc b/AMDiS/src/ProblemVecDbg.cc
new file mode 100644
index 0000000000000000000000000000000000000000..4702df5c9986da7127f141c8e9c4df213659b875
--- /dev/null
+++ b/AMDiS/src/ProblemVecDbg.cc
@@ -0,0 +1,405 @@
+#include <map>
+#include <utility>
+#include "ProblemVecDbg.h"
+#include "SystemVector.h"
+#include "ElementDofIterator.h"
+
+namespace AMDiS {
+
+ void ProblemVecDbg::writeDbgMatrix(std::string filename)
+ {
+ FUNCNAME("ProblemVec::writeDbMatrix()");
+
+ using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end;
+ namespace traits = mtl::traits;
+ typedef DOFMatrix::base_matrix_type Matrix;
+ 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];
+
+ if (!dofmatrix)
+ continue;
+
+ Matrix& matrix = dofmatrix->getBaseMatrix();
+ 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);
+ traits::col<Matrix>::type col(matrix);
+ traits::const_value<Matrix>::type value(matrix);
+
+ for (cursor_type cursor = begin<major>(matrix), cend = end<major>(matrix); cursor != cend; ++cursor) {
+ for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) {
+ out << row(*icursor) << " "
+ << col(*icursor) << " "
+ << value(*icursor) << std::endl;
+ testNnz++;
+ }
+ }
+
+ // Just to check, if nnz() is correct.
+ TEST_EXIT(testNnz == nnz)("NNZ does not fit together!\n");
+ }
+ }
+
+
+ for (int i = 0; i < nComponents; i++) {
+ out << rhs->getDOFVector(i)->getUsedSize() << std::endl;
+ DOFIterator<double> it(rhs->getDOFVector(i), USED_DOFS);
+ int c = 0;
+ for (it.reset(); !it.end(); ++it) {
+ out << c << " " << *it << std::endl;
+ c++;
+ }
+
+ TEST_EXIT(c == rhs->getDOFVector(i)->getUsedSize())
+ ("RHS NNZ does not fit together!\n");
+ }
+
+
+ for (int i = 0; i < nComponents; i++) {
+ out << solution->getDOFVector(i)->getUsedSize() << std::endl;
+ DOFIterator<double> it(solution->getDOFVector(i), USED_DOFS);
+ int c = 0;
+ for (it.reset(); !it.end(); ++it) {
+ out << c << " " << *it << std::endl;
+ c++;
+ }
+
+ TEST_EXIT(c == rhs->getDOFVector(i)->getUsedSize())
+ ("RHS NNZ does not fit together!\n");
+ }
+
+ out.close();
+ }
+
+
+ void ProblemVecDbg::readAndCompareDbgMatrix(std::vector<std::string> filenames)
+ {
+ using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end;
+ namespace traits = mtl::traits;
+ typedef DOFMatrix::base_matrix_type Matrix;
+ 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);
+
+ int dimOfWorld = Global::getGeo(WORLD);
+ std::vector<std::vector<std::map<std::pair<int, int>, double> > > nnzValues(nComponents);
+ std::vector<std::map<int, double> > rhsValues(nComponents);
+ std::vector<std::map<int, double> > solValues(nComponents);
+ for (int i = 0; i < nComponents; i++)
+ nnzValues[i].resize(nComponents);
+
+ // Stores to each dof index of this problem a map from rank indices (of each rank
+ // that also has this dof) to the corresponding local dof index.
+ std::map<DegreeOfFreedom, std::vector<DegreeOfFreedom> > dofMapHereToFiles;
+
+ for (std::vector<std::string>::iterator fileIt = filenames.begin();
+ fileIt != filenames.end(); ++fileIt) {
+
+ // Open file and read the number of DOFs the file contains.
+ std::ifstream in(fileIt->c_str());
+ int nReadDofs;
+ in >> nReadDofs;
+
+ // Is used to map the dof indices in the files to the global coordinates.
+ DofToCoord dofToCoord;
+
+ // Read all DOF indices and their world coordinates.
+ for (int i = 0; i < nReadDofs; i++) {
+ DegreeOfFreedom dof;
+ WorldVector<double> coords;
+
+ in >> dof;
+ for (int j = 0; j < dimOfWorld; j++)
+ in >> coords[j];
+ dofToCoord[dof] = coords;
+ }
+
+ std::map<DegreeOfFreedom, DegreeOfFreedom> dofMapFileToHere;
+
+ for (DofToCoord::iterator it = dofToCoord.begin(); it != dofToCoord.end(); ++it) {
+ DegreeOfFreedom dofIndexInFile = it->first;
+ WorldVector<double> &dofCoords = it->second;
+
+ if (coordToDof.find(dofCoords) == coordToDof.end()) {
+ std::cout << "Cannot find dof index for coords: " << std::endl;
+ dofCoords.print();
+ exit(0);
+ }
+
+ DegreeOfFreedom dofIndexHere = coordToDof[dofCoords];
+ dofMapHereToFiles[dofIndexHere].push_back(dofIndexInFile);
+ dofMapFileToHere[dofIndexInFile] = dofIndexHere;
+ }
+
+ // 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];
+
+ if (!dofmatrix)
+ continue;
+
+ 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;
+
+ in >> row;
+ in >> col;
+ in >> value;
+
+ if (dofMapFileToHere.count(row) == 0) {
+ std::cout << "Cannot find row index for: " << row << std::endl;
+ exit(0);
+ }
+
+ if (dofMapFileToHere.count(col) == 0) {
+ std::cout << "Cannot find col index for: " << col << std::endl;
+ exit(0);
+ }
+
+ // Get dof indices for row and col of this problem matrix.
+ DegreeOfFreedom rowHere = dofMapFileToHere[row];
+ DegreeOfFreedom colHere = dofMapFileToHere[col];
+ std::pair<int, int> rowcol = std::make_pair(rowHere, colHere);
+
+ if (nnzValues[i][j].count(rowcol) == 0)
+ nnzValues[i][j][rowcol] = value;
+ else
+ nnzValues[i][j][rowcol] += value;
+ }
+ }
+ }
+
+ for (int i = 0; i < nComponents; i++) {
+ int readNnz;
+ in >> readNnz;
+ for (int k = 0; k < readNnz; k++) {
+ int row;
+ double value;
+
+ in >> row;
+ in >> value;
+
+ WorldVector<double> rowCoords = dofToCoord[row];
+
+ if (coordToDof.find(rowCoords) == coordToDof.end()) {
+ std::cout << "Cannot find row index for coords: " << std::endl;
+ rowCoords.print();
+ exit(0);
+ }
+
+ DegreeOfFreedom rowHere = coordToDof[rowCoords];
+
+ if (rhsValues[i].count(rowHere) == 0)
+ rhsValues[i][rowHere] = value;
+ else
+ rhsValues[i][rowHere] += value;
+ }
+ }
+
+ for (int i = 0; i < nComponents; i++) {
+ int readNnz;
+ in >> readNnz;
+ for (int k = 0; k < readNnz; k++) {
+ int row;
+ double value;
+
+ in >> row;
+ in >> value;
+
+ WorldVector<double> rowCoords = dofToCoord[row];
+
+ if (coordToDof.find(rowCoords) == coordToDof.end()) {
+ std::cout << "Cannot find row index for coords: " << std::endl;
+ rowCoords.print();
+ exit(0);
+ }
+
+ DegreeOfFreedom rowHere = coordToDof[rowCoords];
+
+ if (solValues[i].count(rowHere) == 0) {
+ solValues[i][rowHere] = value;
+ } else {
+ double diff = fabs(solValues[i][rowHere] - value);
+ if (diff > 1e-8) {
+ std::cout << "DIFFERENT values in solution vector!" << std::endl;
+ exit(0);
+ }
+ }
+ }
+ }
+
+ in.close();
+
+ std::cout << "File read!" << std::endl;
+ }
+
+ std::cout << "Now checking ..." << std::endl;
+
+ std::cout.precision(15);
+
+ for (int i = 0; i < nComponents; i++) {
+ for (int j = 0; j < nComponents; j++) {
+ DOFMatrix *dofmatrix = (*systemMatrix)[i][j];
+
+ if (!dofmatrix)
+ continue;
+
+ for (std::map<std::pair<int, int>, double>::iterator nnzIt =
+ nnzValues[i][j].begin(); nnzIt != nnzValues[i][j].end(); ++nnzIt) {
+
+ int row = nnzIt->first.first;
+ int col = nnzIt->first.second;
+ double value = nnzIt->second;
+
+ double valueHere = (dofmatrix->getBaseMatrix())[row][col];
+ 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-8) {
+ std::cout << "Wrong value in component matrix " << i
+ << "/" << j << ": " << std::endl
+ << " DOF in this matrix[" << row << "][" << col << "] = "
+ << valueHere << std::endl
+ << " DOF in other matrix[" << dofMapHereToFiles[row][0] << "]["
+ << dofMapHereToFiles[col][0] << "] = "
+ << value << std::endl;
+
+ exit(0);
+ }
+ }
+ }
+ }
+
+ for (int i = 0; i < nComponents; i++) {
+ for (std::map<int, double>::iterator rhsIt = rhsValues[i].begin();
+ rhsIt != rhsValues[i].end(); ++rhsIt) {
+
+ int row = rhsIt->first;
+ double value = rhsIt->second;
+
+ double valueHere = (*(rhs->getDOFVector(i)))[row];
+ double diff = fabs(valueHere - value);
+
+ if (diff > 1e-8) {
+ std::cout << "WRONG value in rhs[" << i << "]!" << std::endl
+ << " DOF in other rhs[" << row << "] = " << value << std::endl
+ << " DOF in this rhs[" << row << "] = " << valueHere << std::endl;
+
+ exit(0);
+ }
+ }
+ }
+
+ for (int i = 0; i < nComponents; i++) {
+ for (std::map<int, double>::iterator solIt = solValues[i].begin();
+ solIt != solValues[i].end(); ++solIt) {
+
+ int row = solIt->first;
+ double value = solIt->second;
+
+ double valueHere = (*(solution->getDOFVector(i)))[row];
+ double diff = fabs(valueHere - value);
+
+ if (diff > 1e-8) {
+ std::cout << "WRONG value in sol[" << i << "]!" << std::endl
+ << " DOF in other sol[" << row << "] = " << value << std::endl
+ << " DOF in this sol[" << row << "] = " << valueHere << std::endl;
+
+ exit(0);
+ }
+ }
+ }
+
+ std::cout << "FINISHED COMPARING!" << std::endl;
+
+ exit(0);
+ }
+
+
+ void ProblemVecDbg::createDofToCoordMap(DofToCoord &dofMap)
+ {
+ const BasisFunction *basisFcts = componentSpaces[0]->getBasisFcts();
+ WorldVector<double> coords;
+ ElementDofIterator elDofIter(componentSpaces[0], true);
+
+ TraverseStack stack;
+ ElInfo *elInfo = stack.traverseFirst(componentMeshes[0], -1,
+ Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
+ while (elInfo) {
+ elDofIter.reset(elInfo->getElement());
+ int i = 0;
+ do {
+ DimVec<double> *baryCoord = basisFcts->getCoords(i);
+ elInfo->coordToWorld(*baryCoord, coords);
+ i++;
+ dofMap[elDofIter.getDof()] = coords;
+ } while (elDofIter.next());
+
+ elInfo = stack.traverseNext(elInfo);
+ }
+ }
+
+
+ void ProblemVecDbg::createCoordToDofMap(CoordToDof &dofMap)
+ {
+ const BasisFunction *basisFcts = componentSpaces[0]->getBasisFcts();
+ WorldVector<double> coords;
+ ElementDofIterator elDofIter(componentSpaces[0], true);
+
+ TraverseStack stack;
+ ElInfo *elInfo = stack.traverseFirst(componentMeshes[0], -1,
+ Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
+ while (elInfo) {
+ elDofIter.reset(elInfo->getElement());
+ int i = 0;
+ do {
+ DimVec<double> *baryCoord = basisFcts->getCoords(i);
+ elInfo->coordToWorld(*baryCoord, coords);
+ i++;
+ dofMap[coords] = elDofIter.getDof();
+ } while (elDofIter.next());
+
+ elInfo = stack.traverseNext(elInfo);
+ }
+ }
+
+}
diff --git a/AMDiS/src/ProblemVecDbg.h b/AMDiS/src/ProblemVecDbg.h
new file mode 100644
index 0000000000000000000000000000000000000000..665ca261779507c3c2012244d0363fb1c9daf326
--- /dev/null
+++ b/AMDiS/src/ProblemVecDbg.h
@@ -0,0 +1,70 @@
+// ============================================================================
+// == ==
+// == AMDiS - Adaptive multidimensional simulations ==
+// == ==
+// ============================================================================
+// == ==
+// == TU Dresden ==
+// == ==
+// == Institut für Wissenschaftliches Rechnen ==
+// == Zellescher Weg 12-14 ==
+// == 01069 Dresden ==
+// == germany ==
+// == ==
+// ============================================================================
+// == ==
+// == https://gforge.zih.tu-dresden.de/projects/amdis/ ==
+// == ==
+// ============================================================================
+
+/** \file ProblemVecDbg.h */
+
+#ifndef AMDIS_PROBLEMVECDBG_H
+#define AMDIS_PROBLEMVECDBG_H
+
+#include "ProblemVec.h"
+
+namespace AMDiS {
+
+ class ProblemVecDbg : public ProblemVec
+ {
+ public:
+ ProblemVecDbg(std::string nameStr,
+ ProblemIterationInterface *problemIteration = NULL)
+ : ProblemVec(nameStr, problemIteration)
+ {}
+
+ /** \brief
+ * This function writes the system matrix and the system vector to a file
+ * for debugging. The entries of both, the matrix and the vector, are not
+ * indicated by the dof indices, but by the world coordinates of the dofs.
+ * This makes it possible to compare dof matrices, where the dof indices
+ * are different, e.g., when using domain decomposition parallelization.
+ */
+ void writeDbgMatrix(std::string filename);
+
+ /** \brief
+ * Reads a file, that was created using the function \ref writeDbgMatrix, and
+ * compares the date with the system matrix of this problem.
+ */
+ void readAndCompareDbgMatrix(std::vector<std::string> filenames);
+
+ protected:
+ /** \brief
+ * Create from the current problem a map from dof indices to world coordinates.
+ * This function is used by the debugging function \ref writeDbgMatrix and \ref
+ * readAndCompareDbgMatrix.
+ */
+ void createDofToCoordMap(DofToCoord &dofMap);
+
+ /** \brief
+ * Create from the current problem a map from world coordinates of all dofs to
+ * the corresponding dof indices. This function is used by the debugging function
+ * \ref writeDbgMatrix and \ref readAndCompareDbgMatrix.
+ */
+ void createCoordToDofMap(CoordToDof &dofMap);
+ };
+
+}
+
+#endif