From c14d21292cd7ccbaaee357f39c29e3e8aaf19899 Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Fri, 2 Nov 2012 09:48:52 +0000
Subject: [PATCH] Fixed eps issue in final timestep calculation.
---
AMDiS/src/AdaptInfo.h | 16 ++++-----
AMDiS/src/ElementDofIterator.h | 10 +++---
AMDiS/src/Lagrange.h | 3 +-
AMDiS/src/parallel/PetscSolverFeti.cc | 4 ---
.../src/parallel/PetscSolverFetiOperators.cc | 36 ++-----------------
AMDiS/src/parallel/PetscSolverGlobalMatrix.cc | 13 ++++---
6 files changed, 26 insertions(+), 56 deletions(-)
diff --git a/AMDiS/src/AdaptInfo.h b/AMDiS/src/AdaptInfo.h
index eb095322..8357f513 100644
--- a/AMDiS/src/AdaptInfo.h
+++ b/AMDiS/src/AdaptInfo.h
@@ -511,24 +511,24 @@ namespace AMDiS {
return timestep;
}
- inline void setLastProcessedTimestep(double t){
- lastProcessedTimestep=t;
+ inline void setLastProcessedTimestep(double t)
+ {
+ lastProcessedTimestep = t;
}
- inline double getLastProcessedTimestep(){
+ inline double getLastProcessedTimestep()
+ {
return lastProcessedTimestep;
}
- /** \brief
- * Returns true, if the end time is reached and no more timestep
- * computations must be done.
- */
+ /// Returns true, if the end time is reached and no more timestep
+ /// computations must be done.
inline bool reachedEndTime()
{
if (nTimesteps > 0)
return !(timestepNumber < nTimesteps);
- return !(time < endTime - DBL_TOL);
+ return !(fabs(time - endTime) > DBL_TOL);
}
diff --git a/AMDiS/src/ElementDofIterator.h b/AMDiS/src/ElementDofIterator.h
index 8b02219b..b7579a3d 100644
--- a/AMDiS/src/ElementDofIterator.h
+++ b/AMDiS/src/ElementDofIterator.h
@@ -140,12 +140,10 @@ namespace AMDiS {
/// Dimension dependent geo index of the current position in traverse.
GeoIndex posIndex;
- /** \brief
- * Number of DOFs at the current traverse position. Examples: independent of
- * dimension and degree of basis functions there is only one DOF per vertex.
- * But in 2D and with 3rd degree lagrange basis functions there are two
- * DOFs per edge.
- */
+ /// Number of DOFs at the current traverse position. Examples: independent of
+ /// dimension and degree of basis functions there is only one DOF per vertex.
+ /// But in 2D and with 3rd degree lagrange basis functions there are two
+ /// DOFs per edge.
int nDofs;
/// Displacement of DOF indices. Used if more than one DOF admin is defined
diff --git a/AMDiS/src/Lagrange.h b/AMDiS/src/Lagrange.h
index 310f7ab9..14462f22 100644
--- a/AMDiS/src/Lagrange.h
+++ b/AMDiS/src/Lagrange.h
@@ -154,7 +154,8 @@ namespace AMDiS {
GeoIndex position,
int positionIndex) const;
- /// Calculates the number of DOFs needed for Lagrange of the given dim and degree.
+ /// Calculates the number of DOFs needed for Lagrange of the given dim
+ /// and degree.
static int getNumberOfDofs(int dim, int degree);
private:
diff --git a/AMDiS/src/parallel/PetscSolverFeti.cc b/AMDiS/src/parallel/PetscSolverFeti.cc
index b0863e35..9379f352 100644
--- a/AMDiS/src/parallel/PetscSolverFeti.cc
+++ b/AMDiS/src/parallel/PetscSolverFeti.cc
@@ -1007,12 +1007,9 @@ namespace AMDiS {
Mat qT, jTqT;
MatTranspose(mat_augmented_lagrange, MAT_INITIAL_MATRIX, &qT);
// Mat jT;
- MSG("START COMPUTING MAT TRANS\n");
// MatTranspose(mat_lagrange, MAT_INITIAL_MATRIX, &jT);
- MSG("DONE\n");
MatTransposeMatMult(mat_lagrange, qT, MAT_INITIAL_MATRIX, PETSC_DEFAULT,
&jTqT);
- MSG("DONE WITH THIS!\n");
petsc_helper::blockMatMatSolve(subdomain->getSolver(), jTqT, matTmp);
MatDestroy(&qT);
MatDestroy(&jTqT);
@@ -1041,7 +1038,6 @@ namespace AMDiS {
MatDestroy(&mat10);
MatDestroy(&mat11);
MatDestroy(&matTmp);
- MSG("FINISHED!\n");
} else {
Mat tmp;
diff --git a/AMDiS/src/parallel/PetscSolverFetiOperators.cc b/AMDiS/src/parallel/PetscSolverFetiOperators.cc
index f73c2c06..e1abec27 100644
--- a/AMDiS/src/parallel/PetscSolverFetiOperators.cc
+++ b/AMDiS/src/parallel/PetscSolverFetiOperators.cc
@@ -252,33 +252,12 @@ namespace AMDiS {
double wtime = MPI::Wtime();
- {
- MatNullSpace nullSpace;
- MatGetNullSpace(mat, &nullSpace);
-
- PetscBool hasConst;
- PetscInt nVec;
- const Vec *vecs;
- // MatNullSpaceGetVecs(nullSpace, &hasConst, &nVec, &vecs);
- MatNullSpaceRemove(nullSpace, x, PETSC_NULL);
- }
-
-
Vec x_interface, x_lagrange, y_interface, y_lagrange;
VecNestGetSubVec(x, 0, &x_interface);
VecNestGetSubVec(x, 1, &x_lagrange);
VecNestGetSubVec(y, 0, &y_interface);
VecNestGetSubVec(y, 1, &y_lagrange);
- {
- int n;
- VecGetSize(x_interface, &n);
- double sum;
- VecSum(x_interface, &sum);
- sum = -sum / static_cast<int>(n);
- MSG("xbegin = %e\n", sum);
- }
-
void *ctx;
MatShellGetContext(mat, &ctx);
FetiData* data = static_cast<FetiData*>(ctx);
@@ -337,16 +316,6 @@ namespace AMDiS {
// y_interface = A_{\Gamma B} tmp_vec_b0
MatMult(data->subSolver->getMatCoarseInterior(1), data->tmp_vec_b0, y_interface);
- {
- int n;
- VecGetSize(y_interface, &n);
- double sum;
- VecSum(y_interface, &sum);
- sum = -sum / static_cast<int>(n);
- MSG("yend = %e\n", sum);
- }
-
-
// tmp_vec_primal0 = S_{\Pi\Pi}^{-1} (v_{\Pi} - A_{\Pi B} A_{BB}^{-1} v_{B})
// tmp_vec_interface = A_{\Gamma \Pi} tmp_vec_primal0
MatMult(data->subSolver->getMatCoarse(1, 0), data->tmp_vec_primal0, data->tmp_vec_interface);
@@ -358,6 +327,7 @@ namespace AMDiS {
FetiTimings::fetiSolve += (MPI::Wtime() - wtime);
+
return 0;
}
@@ -497,8 +467,8 @@ namespace AMDiS {
VecNestGetSubVec(yvec, 0, &y_interface);
VecNestGetSubVec(yvec, 1, &y_lagrange);
- VecCopy(x_interface, y_interface);
- // KSPSolve(data->ksp_mass, x_interface, y_interface);
+ //VecCopy(x_interface, y_interface);
+ KSPSolve(data->ksp_mass, x_interface, y_interface);
MatMultTranspose(*(data->mat_lagrange_scaled), x_lagrange, data->tmp_vec_b0);
diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
index 7cb273bf..5d4cf36b 100644
--- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
+++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
@@ -157,7 +157,7 @@ namespace AMDiS {
if (dirichletRows.count(*cursor)) {
if ((!isRowCoarse && !(*interiorMap)[rowComponent].isRankDof(*cursor)) ||
(isRowCoarse && !(*rowCoarseSpace)[rowComponent].isRankDof(*cursor)))
- continue;
+ continue;
}
cols.clear();
@@ -763,11 +763,16 @@ namespace AMDiS {
if (rankOnly && !(*interiorMap)[rowComp].isRankDof(dof))
continue;
+ bool isCoarseDof = isCoarseSpace(rowComp, dof);
+
// Dirichlet rows can be set only be the owner ranks.
- if (dirichletValues.count(dof) && !((*interiorMap)[rowComp].isRankDof(dof)))
- continue;
+ if (dirichletValues.count(dof)) {
+ if ((!isCoarseDof && !((*interiorMap)[rowComp].isRankDof(dof))) ||
+ (isCoarseDof && !((*rowCoarseSpace)[rowComp].isRankDof(dof))))
+ continue;
+ }
- if (isCoarseSpace(rowComp, dof)) {
+ if (isCoarseDof) {
TEST_EXIT_DBG(vecCoarse != PETSC_NULL)("Should not happen!\n");
int index = rowCoarseSpace->getMatIndex(rowComp, dof);
--
GitLab