From f94e0e7ce93fc585176a5909e712b817cbc6027f Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Fri, 2 Nov 2012 21:16:52 +0000
Subject: [PATCH] Changed coarse space creation for FETI-DP method.
---
AMDiS/src/AMDiS.cc | 6 +-
AMDiS/src/parallel/PetscSolverFeti.cc | 136 +++++++++++++--------
AMDiS/src/parallel/PetscSolverFeti.h | 4 +-
AMDiS/src/parallel/PetscSolverFetiDebug.cc | 35 +++---
4 files changed, 109 insertions(+), 72 deletions(-)
diff --git a/AMDiS/src/AMDiS.cc b/AMDiS/src/AMDiS.cc
index 876a63b2..1237ccc1 100644
--- a/AMDiS/src/AMDiS.cc
+++ b/AMDiS/src/AMDiS.cc
@@ -22,7 +22,7 @@ namespace AMDiS {
using namespace std;
#if defined(HAVE_PARALLEL_DOMAIN_AMDIS) && defined(HAVE_PARALLEL_MTL4)
- mtl::par::environment* mtl_environment=NULL;
+ mtl::par::environment* mtl_environment = NULL;
#endif
void init(int argc, char **argv, std::string initFileName)
@@ -37,8 +37,6 @@ namespace AMDiS {
mpi::startRand();
#endif
- Parameters::get("parallel->log main rank", Msg::outputMainRank);
-
#ifdef HAVE_ZOLTAN
float zoltanVersion = 0.0;
Zoltan_Initialize(argc, argv, &zoltanVersion);
@@ -98,6 +96,8 @@ namespace AMDiS {
Parameters::init(initFileName);
}
+ Parameters::get("parallel->log main rank", Msg::outputMainRank);
+
// reset parameters from command line
bool ignoreCommandline = false;
Parameters::get("ignore commandline options", ignoreCommandline);
diff --git a/AMDiS/src/parallel/PetscSolverFeti.cc b/AMDiS/src/parallel/PetscSolverFeti.cc
index 7190ff4d..cb2b5f75 100644
--- a/AMDiS/src/parallel/PetscSolverFeti.cc
+++ b/AMDiS/src/parallel/PetscSolverFeti.cc
@@ -734,51 +734,78 @@ namespace AMDiS {
}
- std::set<BoundaryObject> PetscSolverFeti::getAugmentedLagrange()
+ vector<vector<BoundaryObject> > PetscSolverFeti::getCoarseEdges()
{
FUNCNAME("PetscSolverFeti::getAugmentedLagrange()");
- int nEdges = 0;
- int nFaces = 0;
-
InteriorBoundary &intBound = meshDistributor->getIntBoundary();
std::set<BoundaryObject> allEdges;
- for (InteriorBoundary::iterator it(intBound.getOwn()); !it.end(); ++it) {
- if ((it->rankObj.subObj == FACE ||
- (it->rankObj.subObj == EDGE &&
- testWirebasketEdge(it->rankObj, feSpaces[0]))) &&
- allEdges.count(it->rankObj) == 0) {
- bool dirichletOnlyEdge = true;
-
- DofContainer edgeDofs;
- it->rankObj.el->getAllDofs(feSpaces[0], it->rankObj, edgeDofs);
-
- for (DofContainer::iterator dit = edgeDofs.begin();
- dit != edgeDofs.end(); ++dit) {
- if (dirichletRows[0].count(**dit) == 0) {
- dirichletOnlyEdge = false;
- break;
- }
- }
+ for (InteriorBoundary::iterator it(intBound.getOwn()); !it.end(); ++it)
+ if (it->rankObj.subObj == EDGE &&
+ testWirebasketEdge(it->rankObj, feSpaces[0]) &&
+ allEdges.count(it->rankObj) == 0)
+ allEdges.insert(it->rankObj);
+
+ int nEmptyEdges = 0;
+ vector<vector<BoundaryObject> > data;
+ for (std::set<BoundaryObject>::iterator it = allEdges.begin();
+ it != allEdges.end(); ++it) {
+ DofContainer edgeDofs;
+ it->el->getAllDofs(feSpaces[0], *it, edgeDofs);
+ if (edgeDofs.size() == 0) {
+ nEmptyEdges++;
+ } else {
+ vector<BoundaryObject> oneBoundary;
+ oneBoundary.push_back(*it);
+ data.push_back(oneBoundary);
+ }
+ }
- if (!dirichletOnlyEdge) {
- if (it->rankObj.subObj == EDGE)
- nEdges++;
-
- if (it->rankObj.subObj == FACE)
- nFaces++;
+ int nEdges = allEdges.size();
+ mpi::globalAdd(nEdges);
+ mpi::globalAdd(nEmptyEdges);
- allEdges.insert(it->rankObj);
- }
+ MSG("Coarse space edges: %d (empty: %d)\n", nEdges, nEmptyEdges);
+
+ return data;
+ }
+
+
+ vector<vector<BoundaryObject> > PetscSolverFeti::getCoarseFaces()
+ {
+ FUNCNAME("PetscSolverFeti::getAugmentedLagrange()");
+
+ InteriorBoundary &intBound = meshDistributor->getIntBoundary();
+ map<int, std::set<BoundaryObject> > allFaces;
+ for (InteriorBoundary::iterator it(intBound.getOwn()); !it.end(); ++it)
+ if (it->rankObj.subObj == FACE)
+ allFaces[it.getRank()].insert(it->rankObj);
+
+ int nEmptyFaces = 0;
+ vector<vector<BoundaryObject> > data;
+ for (map<int, std::set<BoundaryObject> >::iterator it = allFaces.begin();
+ it != allFaces.end(); ++it) {
+ vector<BoundaryObject> oneBoundary;
+ for (std::set<BoundaryObject>::iterator bIt = it->second.begin();
+ bIt != it->second.end(); ++bIt) {
+ DofContainer faceDofs;
+ bIt->el->getAllDofs(feSpaces[0], *bIt, faceDofs);
+ if (faceDofs.size())
+ oneBoundary.push_back(*bIt);
}
+ if (oneBoundary.size())
+ data.push_back(oneBoundary);
+ else
+ nEmptyFaces++;
}
- mpi::globalAdd(nEdges);
+ int nFaces = allFaces.size();
mpi::globalAdd(nFaces);
+ mpi::globalAdd(nEmptyFaces);
- MSG("NEDGES = %d NFACES = %d\n", nEdges, nFaces);
-
- return allEdges;
+ MSG("Coarse space faces: %d (empty: %d)\n", nFaces, nEmptyFaces);
+
+ return data;
}
@@ -791,8 +818,11 @@ namespace AMDiS {
double wtime = MPI::Wtime();
- std::set<BoundaryObject> allEdges = getAugmentedLagrange();
- nRankEdges = allEdges.size();
+ vector<vector<BoundaryObject> > allEdges = getCoarseEdges();
+ vector<vector<BoundaryObject> > allFaces = getCoarseFaces();
+ allEdges.insert(allEdges.end(), allFaces.begin(), allFaces.end());
+
+ nRankEdges = allEdges.size();
int rStartEdges = 0;
mpi::getDofNumbering(mpiCommGlobal, nRankEdges, rStartEdges, nOverallEdges);
@@ -810,27 +840,29 @@ namespace AMDiS {
MatSetOption(mat_augmented_lagrange, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE);
int rowCounter = rStartEdges;
- for (std::set<BoundaryObject>::iterator edgeIt = allEdges.begin();
- edgeIt != allEdges.end(); ++edgeIt) {
+ for (vector<vector<BoundaryObject> >::iterator it = allEdges.begin();
+ it != allEdges.end(); ++it) {
for (int component = 0; component < componentSpaces.size(); component++) {
- DofContainer edgeDofs;
- edgeIt->el->getAllDofs(componentSpaces[component], *edgeIt, edgeDofs);
-
- for (DofContainer::iterator it = edgeDofs.begin();
- it != edgeDofs.end(); ++it) {
- TEST_EXIT_DBG(isPrimal(component, **it) == false)
- ("Should not be primal!\n");
+ for (vector<BoundaryObject>::iterator edgeIt = it->begin();
+ edgeIt != it->end(); ++edgeIt) {
- if (dirichletRows[component].count(**it))
- continue;
+ DofContainer edgeDofs;
+ edgeIt->el->getAllDofs(componentSpaces[component], *edgeIt, edgeDofs);
+
+ TEST_EXIT(edgeDofs.size())("Should not happen!\n");
- int col = lagrangeMap.getMatIndex(component, **it);
- double value = 1.0;
- MatSetValue(mat_augmented_lagrange, rowCounter, col, value, INSERT_VALUES);
+ for (DofContainer::iterator it = edgeDofs.begin();
+ it != edgeDofs.end(); ++it) {
+ TEST_EXIT(isPrimal(component, **it) == false)
+ ("Should not be primal!\n");
+
+ int col = lagrangeMap.getMatIndex(component, **it);
+ double value = 1.0;
+ MatSetValue(mat_augmented_lagrange, rowCounter, col, value, INSERT_VALUES);
+ }
}
-
- rowCounter++;
- }
+ rowCounter++;
+ }
}
MatAssemblyBegin(mat_augmented_lagrange, MAT_FINAL_ASSEMBLY);
diff --git a/AMDiS/src/parallel/PetscSolverFeti.h b/AMDiS/src/parallel/PetscSolverFeti.h
index 9e308d8d..58aab40e 100644
--- a/AMDiS/src/parallel/PetscSolverFeti.h
+++ b/AMDiS/src/parallel/PetscSolverFeti.h
@@ -131,7 +131,9 @@ namespace AMDiS {
/// to \ref mat_lagrange.
void createMatLagrange();
- std::set<BoundaryObject> getAugmentedLagrange();
+ vector<vector<BoundaryObject> > getCoarseEdges();
+
+ vector<vector<BoundaryObject> > getCoarseFaces();
void createMatAugmentedLagrange();
diff --git a/AMDiS/src/parallel/PetscSolverFetiDebug.cc b/AMDiS/src/parallel/PetscSolverFetiDebug.cc
index dc7ec495..04510815 100644
--- a/AMDiS/src/parallel/PetscSolverFetiDebug.cc
+++ b/AMDiS/src/parallel/PetscSolverFetiDebug.cc
@@ -641,23 +641,26 @@ namespace AMDiS {
StdMpi<vector<double> > stdMpi(MPI::COMM_WORLD);
vector<double> sendData;
- std::set<BoundaryObject> coarseSpace = feti.getAugmentedLagrange();
- for (std::set<BoundaryObject>::iterator it = coarseSpace.begin();
+ vector<vector<BoundaryObject> > coarseSpace = feti.getCoarseFaces();
+ for (vector<vector<BoundaryObject> >::iterator it = coarseSpace.begin();
it != coarseSpace.end(); ++it) {
- if (it->subObj == FACE) {
- DofFace face = it->el->getFace(it->ithObj);
- WorldVector<double> c0 = coords[face.get<0>()];
- WorldVector<double> c1 = coords[face.get<1>()];
- WorldVector<double> c2 = coords[face.get<2>()];
- sendData.push_back(c0[0]);
- sendData.push_back(c0[1]);
- sendData.push_back(c0[2]);
- sendData.push_back(c1[0]);
- sendData.push_back(c1[1]);
- sendData.push_back(c1[2]);
- sendData.push_back(c2[0]);
- sendData.push_back(c2[1]);
- sendData.push_back(c2[2]);
+ for (vector<BoundaryObject>::iterator edgeIt = it->begin();
+ edgeIt != it->end(); ++edgeIt) {
+ if (edgeIt->subObj == FACE) {
+ DofFace face = edgeIt->el->getFace(edgeIt->ithObj);
+ WorldVector<double> c0 = coords[face.get<0>()];
+ WorldVector<double> c1 = coords[face.get<1>()];
+ WorldVector<double> c2 = coords[face.get<2>()];
+ sendData.push_back(c0[0]);
+ sendData.push_back(c0[1]);
+ sendData.push_back(c0[2]);
+ sendData.push_back(c1[0]);
+ sendData.push_back(c1[1]);
+ sendData.push_back(c1[2]);
+ sendData.push_back(c2[0]);
+ sendData.push_back(c2[1]);
+ sendData.push_back(c2[2]);
+ }
}
}
--
GitLab