From e10fa30c3c7bbc3e333da9277b47120f4efa546e Mon Sep 17 00:00:00 2001 From: Thomas Witkowski <thomas.witkowski@gmx.de> Date: Wed, 13 Apr 2011 07:55:32 +0000 Subject: [PATCH] And some more bug fixes for parallel code. --- AMDiS/libtool | 6 +- AMDiS/src/MeshStructure.cc | 30 +++++ AMDiS/src/MeshStructure.h | 2 + AMDiS/src/RefinementManager3d.cc | 8 +- AMDiS/src/io/ElementFileWriter.cc | 172 ++++++++++++++----------- AMDiS/src/io/VtkWriter.cc | 14 +- AMDiS/src/parallel/MeshDistributor.cc | 9 +- AMDiS/src/parallel/MeshManipulation.cc | 157 +++++++++++++--------- 8 files changed, 238 insertions(+), 160 deletions(-) diff --git a/AMDiS/libtool b/AMDiS/libtool index 20ca6d75..a502627c 100755 --- a/AMDiS/libtool +++ b/AMDiS/libtool @@ -44,7 +44,7 @@ available_tags=" CXX F77" # ### BEGIN LIBTOOL CONFIG -# Libtool was configured on host deimos102: +# Libtool was configured on host p1q024: # 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 deimos102: +# Libtool was configured on host p1q024: # 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 deimos102: +# Libtool was configured on host p1q024: # Shell to use when invoking shell scripts. SHELL="/bin/sh" diff --git a/AMDiS/src/MeshStructure.cc b/AMDiS/src/MeshStructure.cc index dc54ab8a..0d72c7b4 100644 --- a/AMDiS/src/MeshStructure.cc +++ b/AMDiS/src/MeshStructure.cc @@ -190,6 +190,36 @@ namespace AMDiS { } + int MeshStructure::lookAhead(unsigned int n) + { + FUNCNAME("MeshStructure::lookAhead()"); + + int returnValue = 0; + + int tmp_pos = pos; + int tmp_currentElement = currentElement; + int tmp_currentIndex = currentIndex; + uint64_t tmp_currentCode = currentCode; + + for (unsigned int i = 0; i < n; i++) { + if (nextElement() == false) { + returnValue = -1; + break; + } + } + + if (returnValue != -1) + returnValue = static_cast<int>(!isLeafElement()); + + pos = tmp_pos; + currentElement = tmp_currentElement; + currentIndex = tmp_currentIndex; + currentCode = tmp_currentCode; + + return returnValue; + } + + bool MeshStructure::skipBranch(MeshStructure *insert) { FUNCNAME("MeshStructure::skipBranch()"); diff --git a/AMDiS/src/MeshStructure.h b/AMDiS/src/MeshStructure.h index f062fcee..97ff767b 100644 --- a/AMDiS/src/MeshStructure.h +++ b/AMDiS/src/MeshStructure.h @@ -91,6 +91,8 @@ namespace AMDiS { bool nextElement(MeshStructure *insert = NULL); + int lookAhead(unsigned int n = 1); + inline bool isLeafElement() { return (currentCode & 1) == 0; diff --git a/AMDiS/src/RefinementManager3d.cc b/AMDiS/src/RefinementManager3d.cc index 80927697..7a34ab91 100644 --- a/AMDiS/src/RefinementManager3d.cc +++ b/AMDiS/src/RefinementManager3d.cc @@ -675,6 +675,8 @@ namespace AMDiS { bool foundEdge = false; while (elInfo2) { + MSG("TRY TO FIND ON EL %d\n", elInfo2->getElement()->getIndex()); + for (int i = 0; i < 6; i++) { DofEdge edge2 = elInfo2->getElement()->getEdge(i); if (edge2.first == *(edge[0]) && @@ -688,13 +690,14 @@ namespace AMDiS { // must be refined at least once to get a refinement edge. if (i == 0) { + MSG("AND FOUND!\n"); + // Edge is refinement edge, so add it to refine list. refineList.setElType(n_neigh, elInfo2->getType()); refineList.setElement(n_neigh, elInfo2->getElement()); n_neigh++; - foundEdge = true; TraverseStack *tmpStack = stack; stack = &stack2; @@ -704,11 +707,12 @@ namespace AMDiS { } stack = tmpStack; + foundEdge = true; break; } else { // Edge i not refinement edge, so refine the element further. - Element *el2 = elInfo->getElement(); + Element *el2 = elInfo2->getElement(); el2->setMark(std::max(el2->getMark(), 1)); TraverseStack *tmpStack = stack; diff --git a/AMDiS/src/io/ElementFileWriter.cc b/AMDiS/src/io/ElementFileWriter.cc index d1262447..05bf0789 100644 --- a/AMDiS/src/io/ElementFileWriter.cc +++ b/AMDiS/src/io/ElementFileWriter.cc @@ -10,6 +10,11 @@ // See also license.opensource.txt in the distribution. +#include <boost/iostreams/filtering_stream.hpp> +#include <boost/iostreams/device/file_descriptor.hpp> +#include <boost/filesystem/operations.hpp> +#include <boost/filesystem/convenience.hpp> + #include "ElementFileWriter.h" #include "BasisFunction.h" #include "Parameters.h" @@ -18,9 +23,11 @@ namespace AMDiS { - ElementFileWriter::ElementFileWriter(std::string name_, + using namespace std; + + ElementFileWriter::ElementFileWriter(string name_, Mesh *mesh_, - std::map<int, double> &mapvec) + map<int, double> &mapvec) : name(name_), amdisMeshDatExt(".elem.mesh"), vtkExt(".vtu"), @@ -59,7 +66,7 @@ namespace AMDiS { if (timestepNumber != 0 && !force) return; - std::string fn = filename; + string fn = filename; if (appendIndex) { TEST_EXIT(indexLength <= 99)("index lenght > 99\n"); @@ -92,9 +99,9 @@ namespace AMDiS { } - void ElementFileWriter::writeFile(std::map<int, double> &vec, + void ElementFileWriter::writeFile(map<int, double> &vec, Mesh *mesh, - std::string filename, + string filename, int level) { ElementFileWriter efw("", mesh, vec); @@ -102,27 +109,34 @@ namespace AMDiS { } - void ElementFileWriter::writeMeshDatValues(std::string filename, double time) + void ElementFileWriter::writeMeshDatValues(string filename, double time) { FUNCNAME("ElementFileWriter::writeMeshDatValues()"); - std::ofstream fout(filename.c_str()); - TEST_EXIT(fout)("Could not open file %s !\n", filename.c_str()); + boost::iostreams::filtering_ostream file; + { + //boost::iostreams seems not to truncate the file + ofstream swapfile(filename.c_str(), ios::out | ios::trunc); + TEST_EXIT(swapfile.is_open()) + ("Cannot open file %s for writing!\n", name.c_str()); + swapfile.close(); + } + file.push(boost::iostreams::file_descriptor_sink(filename, ios::trunc)); int dim = mesh->getDim(); double val; // === Write header. === - fout << "mesh name: " << mesh->getName().c_str() << "\n\n"; - fout << "time: " << time << "\n\n"; - fout << "DIM: " << dim << "\n"; - fout << "DIM_OF_WORLD: " << Global::getGeo(WORLD) << "\n\n"; - fout << "number of vertices: " << (dim+1)*mesh->getNumberOfLeaves() << "\n"; - fout << "number of elements: " << mesh->getNumberOfLeaves() << "\n\n"; + file << "mesh name: " << mesh->getName().c_str() << "\n\n"; + file << "time: " << time << "\n\n"; + file << "DIM: " << dim << "\n"; + file << "DIM_OF_WORLD: " << Global::getGeo(WORLD) << "\n\n"; + file << "number of vertices: " << (dim+1)*mesh->getNumberOfLeaves() << "\n"; + file << "number of elements: " << mesh->getNumberOfLeaves() << "\n\n"; // === Write vertex coordinates (every vertex for every element). === - fout << "vertex coordinates:\n"; + file << "vertex coordinates:\n"; TraverseStack stack; ElInfo *elInfo = stack.traverseFirst(mesh, @@ -135,9 +149,9 @@ namespace AMDiS { // Write coordinates of all element vertices. for (int i = 0; i <= dim; ++i) { for (int j = 0; j < dim; ++j) { - fout << elInfo->getCoord(i)[j] << " "; + file << elInfo->getCoord(i)[j] << " "; } - fout << "\n"; + file << "\n"; } elInfo = stack.traverseNext(elInfo); @@ -147,33 +161,33 @@ namespace AMDiS { // === Write elements. === int numLeaves = mesh->getNumberOfLeaves(); int vertCntr = 0; - fout << "\n"; - fout << "element vertices:\n"; + file << "\n"; + file << "element vertices:\n"; for (int i = 0; i < numLeaves; ++i) { for (int j = 0; j <= dim; ++j) { - fout << vertCntr << " "; + file << vertCntr << " "; ++vertCntr; } - fout << "\n"; + file << "\n"; } // === Write values. === // Write values header. - fout << "\n"; - fout << "number of values: 1\n\n"; - fout << "value description: " << name.c_str() << "\n"; - fout << "number of interpolation points: 0" << "\n"; - fout << "type: scalar" << "\n"; - fout << "interpolation type: lagrange" << "\n"; - fout << "interpolation degree: 1" << "\n"; - fout << "end of description: " << name.c_str() << "\n\n"; + file << "\n"; + file << "number of values: 1\n\n"; + file << "value description: " << name.c_str() << "\n"; + file << "number of interpolation points: 0" << "\n"; + file << "type: scalar" << "\n"; + file << "interpolation type: lagrange" << "\n"; + file << "interpolation degree: 1" << "\n"; + file << "end of description: " << name.c_str() << "\n\n"; // Write values. - fout << "vertex values: " << name.c_str() << "\n"; + file << "vertex values: " << name.c_str() << "\n"; - fout.setf(std::ios::scientific,std::ios::floatfield); + file.setf(ios::scientific,ios::floatfield); elInfo = stack.traverseFirst(mesh, -1, @@ -185,7 +199,7 @@ namespace AMDiS { // Write value for each vertex of each element. for (int i = 0; i <= dim; ++i) { - fout << val << "\n"; + file << val << "\n"; } elInfo = stack.traverseNext(elInfo); @@ -193,21 +207,25 @@ namespace AMDiS { // Write values trailor. - fout << "\n"; - fout << "interpolation values: " << name.c_str() << "\n\n\n"; - fout << "element interpolation points: " << name.c_str() << "\n"; - - - fout.close(); + file << "\n"; + file << "interpolation values: " << name.c_str() << "\n\n\n"; + file << "element interpolation points: " << name.c_str() << "\n"; } - void ElementFileWriter::writeVtkValues(std::string filename, int level) + void ElementFileWriter::writeVtkValues(string filename, int level) { FUNCNAME("ElementFileWriter::writeVtkValues()"); - std::ofstream fout(filename.c_str()); - TEST_EXIT(fout)("Could not open file %s !\n", filename.c_str()); + boost::iostreams::filtering_ostream file; + { + //boost::iostreams seems not to truncate the file + ofstream swapfile(filename.c_str(), ios::out | ios::trunc); + TEST_EXIT(swapfile.is_open()) + ("Cannot open file %s for writing!\n", name.c_str()); + swapfile.close(); + } + file.push(boost::iostreams::file_descriptor_sink(filename, ios::trunc)); int dim = mesh->getDim(); int dimOfWorld = Global::getGeo(WORLD); @@ -226,13 +244,13 @@ namespace AMDiS { } } - fout << "<?xml version=\"1.0\"?>" << std::endl; - fout << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl; - fout << " <UnstructuredGrid>" << std::endl; - fout << " <Piece NumberOfPoints=\"" << (dim + 1) * nElements - << "\" NumberOfCells=\"" << nElements << "\">" << std::endl; - fout << " <Points>" << std::endl; - fout << " <DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\"ascii\">" << std::endl; + file << "<?xml version=\"1.0\"?>\n"; + file << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n"; + file << " <UnstructuredGrid>\n"; + file << " <Piece NumberOfPoints=\"" << (dim + 1) * nElements + << "\" NumberOfCells=\"" << nElements << "\">\n"; + file << " <Points>\n"; + file << " <DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\"ascii\">\n"; // === Write vertex coordinates (every vertex for every element). === @@ -244,61 +262,61 @@ namespace AMDiS { // Write coordinates of all element vertices. for (int i = 0; i <= dim; i++) { for (int j = 0; j < dimOfWorld; j++) - fout << elInfo->getCoord(i)[j] << " "; + file << elInfo->getCoord(i)[j] << " "; for (int j = dimOfWorld; j < 3; j++) - fout << "0.0"; + file << "0.0"; - fout << "\n"; + file << "\n"; } elInfo = stack.traverseNext(elInfo); } - fout << " </DataArray>" << std::endl; - fout << " </Points>" << std::endl; - fout << " <Cells>" << std::endl; + file << " </DataArray>\n"; + file << " </Points>\n"; + file << " <Cells>\n"; - fout << " <DataArray type=\"Int32\" Name=\"offsets\">" << std::endl; + file << " <DataArray type=\"Int32\" Name=\"offsets\">\n"; for (int i = 0; i < nElements; i++) - fout << " " << (i + 1) * vertices << std::endl; - fout << " </DataArray>" << std::endl; + file << " " << (i + 1) * vertices << "\n"; + file << " </DataArray>\n"; - fout << " <DataArray type=\"UInt8\" Name=\"types\">" << std::endl; + file << " <DataArray type=\"UInt8\" Name=\"types\">\n"; for (int i = 0; i < nElements; i++) { switch (vertices) { case 2: - fout << " 3" << std::endl; + file << " 3\n"; break; case 3: - fout << " 5" << std::endl; + file << " 5\n"; break; case 4: - fout << " 10" << std::endl; + file << " 10\n"; break; default: break; } } - fout << " </DataArray>" << std::endl; + file << " </DataArray>\n"; - fout << " <DataArray type=\"Int32\" Name=\"connectivity\">" << std::endl; + file << " <DataArray type=\"Int32\" Name=\"connectivity\">\n"; int vertCntr = 0; for (int i = 0; i < nElements; ++i) { for (int j = 0; j <= dim; ++j) { - fout << vertCntr << " "; + file << vertCntr << " "; ++vertCntr; } - fout << std::endl; + file << "\n"; } - fout << " </DataArray>" << std::endl; + file << " </DataArray>\n"; - fout << " </Cells>" << std::endl; - fout << " <PointData>" << std::endl; - fout << " <DataArray type=\"Float32\" Name=\"value\" format=\"ascii\">" << std::endl; + file << " </Cells>\n"; + file << " <PointData>\n"; + file << " <DataArray type=\"Float32\" Name=\"value\" format=\"ascii\">\n"; - fout.setf(std::ios::scientific,std::ios::floatfield); + file.setf(ios::scientific,ios::floatfield); elInfo = stack.traverseFirst(mesh, level, traverseFlag | Mesh::FILL_COORDS); int vc = 0; @@ -308,20 +326,18 @@ namespace AMDiS { // Write value for each vertex of each element. for (int i = 0; i <= dim; i++) - fout << (fabs(val) < 1e-40 ? 0.0 : val) << "\n"; + file << (fabs(val) < 1e-40 ? 0.0 : val) << "\n"; vc++; elInfo = stack.traverseNext(elInfo); } - fout << " </DataArray>" << std::endl; + file << " </DataArray>\n"; - fout << " </PointData>" << std::endl; - fout << " </Piece>" << std::endl; - fout << " </UnstructuredGrid>" << std::endl; - fout << "</VTKFile>" << std::endl; - - fout.close(); + file << " </PointData>\n"; + file << " </Piece>\n"; + file << " </UnstructuredGrid>\n"; + file << "</VTKFile>\n"; } } diff --git a/AMDiS/src/io/VtkWriter.cc b/AMDiS/src/io/VtkWriter.cc index 41db97c3..b1e9008c 100644 --- a/AMDiS/src/io/VtkWriter.cc +++ b/AMDiS/src/io/VtkWriter.cc @@ -53,19 +53,9 @@ namespace AMDiS { std::ofstream swapfile(name.c_str(), std::ios::out | std::ios::trunc); swapfile.close(); } - file.push(boost::iostreams::file_descriptor_sink(name, std::ios::trunc)); + file.push(boost::iostreams::file_descriptor_sink(name, std::ios::trunc)); writeFileToStream(file); - -#if 0 - - std::ofstream file; - file.open(name.c_str()); - TEST_EXIT(file.is_open())("Cannot open file %s for writing\n", name.c_str()); - writeFileToStream(file); - file.close(); - -#endif - + return 0; } diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc index 5a38b94d..bb99f1cb 100644 --- a/AMDiS/src/parallel/MeshDistributor.cc +++ b/AMDiS/src/parallel/MeshDistributor.cc @@ -818,6 +818,10 @@ namespace AMDiS { nMeshChangesAfterLastRepartitioning >= repartitionIthChange) { repartitionMesh(); nMeshChangesAfterLastRepartitioning = 0; + } else { + MSG_DBG("Repartitioning not tried because tryRepartitioning = %d reparttitioningAllowed = %d nMeshChange =%d repartitionIthChange = %d\n", + tryRepartition, repartitioningAllowed, + nMeshChangesAfterLastRepartitioning, repartitionIthChange); } // === Print imbalance factor. === @@ -881,6 +885,10 @@ namespace AMDiS { elCode.init(boundIt->rankObj); if (elCode.getCode() != recvCodes[i].getCode()) { +// MSG("MACRO EL %d %d %d DOES NOT FIT WITH MACRO EL %d %d %d ON RANK %d\n", +// boundIt->rankObj.elIndex, boundIt->rankObj.subObj, boundIt->rankObj.ithObj, +// boundIt->neighObj.elIndex, boundIt->neighObj.subObj, boundIt->neighObj.ithObj, +// it->first); TEST_EXIT_DBG(refineManager)("Refinement manager is not set correctly!\n"); MeshManipulation mm(feSpace); @@ -1051,7 +1059,6 @@ namespace AMDiS { elInfo = stack.traverseNext(elInfo); } - // === Run mesh partitioner to calculate a new mesh partitioning. === partitioner->setLocalGlobalDofMap(&mapLocalGlobalDofs); diff --git a/AMDiS/src/parallel/MeshManipulation.cc b/AMDiS/src/parallel/MeshManipulation.cc index b0439acf..b982e741 100644 --- a/AMDiS/src/parallel/MeshManipulation.cc +++ b/AMDiS/src/parallel/MeshManipulation.cc @@ -374,82 +374,111 @@ namespace AMDiS { bool meshChanged = false; Element *el = elInfo->getElement(); - + int s0 = el->getSubObjOfChild(0, subObj, ithObj, elInfo->getType()); + int s1 = el->getSubObjOfChild(1, subObj, ithObj, elInfo->getType()); // === If element is leaf (and the code is not), refine the element. === + bool elementRefined = true; + if (el->isLeaf()) { TEST_EXIT_DBG(elInfo->getLevel() < 255)("This should not happen!\n"); - el->setMark(1); - refineManager->setMesh(el->getMesh()); - refineManager->setStack(&stack); - refineManager->refineFunction(elInfo); - meshChanged = true; - } - - - // === Continue fitting the mesh structure code to the children of the === - // === current element. === + // In some situations refinement of an element can be omitted, altough + // it is defined in the mesh structure code. One case is when the + // following holds: + // - the code is used to adapt a face of a tetrahedron + // - only one of the children of the current element contains this face + // - the face to be adapted of the current element is either the + // local face 0 or 1 + // - the next structure code is 0, i.e., we would be finished after + // the refinement of this element + // In this scenario, the element must be refined due to the structure + // code, but the refinement does not introduce new DOFs on the face, + // that should be adapted. Thus, we can ommit the refinement. + + if (subObj == FACE) { + if (s0 != -1 && s1 == -1 || s0 == -1 && s1 != -1) { + if (ithObj <= 1 && code.lookAhead() == 0) { + elementRefined = false; + code.nextElement(); + stack.traverseNext(elInfo); + } + } + } - int s0 = el->getSubObjOfChild(0, subObj, ithObj, elInfo->getType()); - int s1 = el->getSubObjOfChild(1, subObj, ithObj, elInfo->getType()); - Element *child0 = el->getFirstChild(); - Element *child1 = el->getSecondChild(); - if (reverseMode) { - swap(s0, s1); - swap(child0, child1); + if (elementRefined) { + MeshStructure t = code; + el->setMark(1); + refineManager->setMesh(el->getMesh()); + refineManager->setStack(&stack); + refineManager->refineFunction(elInfo); + meshChanged = true; + } } - - // === Traverse left child. === - - if (s0 != -1) { - // The edge/face is contained in the left child, therefore fit this - // child to the mesh structure code. - - stack.traverseNext(elInfo); - code.nextElement(); - meshChanged |= fitElementToMeshCode(code, stack, subObj, s0, reverseMode); - elInfo = stack.getElInfo(); - } else { - // The edge/face is not contained in the left child. Hence we need - // to traverse through all subelements of the left child until we get - // the second child of the current element. - do { - elInfo = stack.traverseNext(elInfo); - } while (elInfo && elInfo->getElement() != child1); - - TEST_EXIT_DBG(elInfo != NULL)("This should not happen!\n"); - } - - TEST_EXIT_DBG(elInfo->getElement() == child1) - ("Got wrong child with idx = %d! Searched for child idx = %d\n", - elInfo->getElement()->getIndex(), child1->getIndex()); - - - // === Traverse right child. === - - if (s1 != -1) { - // The edge/face is contained in the right child, therefore fit this - // child to the mesh structure code. - - code.nextElement(); - meshChanged |= fitElementToMeshCode(code, stack, subObj, s1, reverseMode); - } else { - // The edge/face is not contained in the right child. Hence we need - // to traverse through all subelements of the right child until we have - // finished traversing the current element with all its subelements. - - int level = elInfo->getLevel(); + if (elementRefined) { + // === Continue fitting the mesh structure code to the children of the === + // === current element. === + + Element *child0 = el->getFirstChild(); + Element *child1 = el->getSecondChild(); + if (reverseMode) { + swap(s0, s1); + swap(child0, child1); + } - do { - elInfo = stack.traverseNext(elInfo); - } while (elInfo && elInfo->getLevel() > level); + + // === Traverse left child. === + + if (s0 != -1) { + // The edge/face is contained in the left child, therefore fit this + // child to the mesh structure code. + + stack.traverseNext(elInfo); + code.nextElement(); + meshChanged |= fitElementToMeshCode(code, stack, subObj, s0, reverseMode); + elInfo = stack.getElInfo(); + } else { + // The edge/face is not contained in the left child. Hence we need + // to traverse through all subelements of the left child until we get + // the second child of the current element. + + do { + elInfo = stack.traverseNext(elInfo); + } while (elInfo && elInfo->getElement() != child1); + + TEST_EXIT_DBG(elInfo != NULL)("This should not happen!\n"); + } + + TEST_EXIT_DBG(elInfo->getElement() == child1) + ("Got wrong child with idx = %d! Searched for child idx = %d\n", + elInfo->getElement()->getIndex(), child1->getIndex()); + + + // === Traverse right child. === + + if (s1 != -1) { + // The edge/face is contained in the right child, therefore fit this + // child to the mesh structure code. + + code.nextElement(); + meshChanged |= fitElementToMeshCode(code, stack, subObj, s1, reverseMode); + } else { + // The edge/face is not contained in the right child. Hence we need + // to traverse through all subelements of the right child until we have + // finished traversing the current element with all its subelements. + + int level = elInfo->getLevel(); + + do { + elInfo = stack.traverseNext(elInfo); + } while (elInfo && elInfo->getLevel() > level); + } + } - - + return meshChanged; } -- GitLab