diff --git a/MathLib/LinAlg/PETSc/PETScLinearSolver.cpp b/MathLib/LinAlg/PETSc/PETScLinearSolver.cpp
index e280ef38e9f305b2723184e6733a277fbee630df..6687a67bbe344c59e29383452a23e6c7fda39495 100644
--- a/MathLib/LinAlg/PETSc/PETScLinearSolver.cpp
+++ b/MathLib/LinAlg/PETSc/PETScLinearSolver.cpp
@@ -25,15 +25,15 @@ PETScLinearSolver::PETScLinearSolver(PETScMatrix &A, const std::string &prefix)
KSPGetPC(_solver, &_pc);
//
- if ( !prefix.empty() )
- {
- KSPSetOptionsPrefix(_solver, prefix.c_str());
+ if ( !prefix.empty() )
+ {
+ KSPSetOptionsPrefix(_solver, prefix.c_str());
}
-
+
KSPSetFromOptions(_solver); // set running time option
}
-void PETScLinearSolver::solve(const PETScVector &b, PETScVector &x)
+bool PETScLinearSolver::solve(const PETScVector &b, PETScVector &x)
{
// define TEST_MEM_PETSC
#ifdef TEST_MEM_PETSC
@@ -48,16 +48,8 @@ void PETScLinearSolver::solve(const PETScVector &b, PETScVector &x)
KSPConvergedReason reason;
KSPGetConvergedReason(_solver, &reason);
- if(reason == KSP_DIVERGED_INDEFINITE_PC)
- {
- PetscPrintf(PETSC_COMM_WORLD,"\nDivergence because of indefinite preconditioner;\n");
- PetscPrintf(PETSC_COMM_WORLD,"Run the executable again but with -pc_factor_shift_positive_definite option.\n");
- }
- else if(reason < 0)
- {
- PetscPrintf(PETSC_COMM_WORLD,"\nOther kind of divergence: this should not happen.\n");
- }
- else
+ bool converged = true;
+ if(reason > 0)
{
const char *slv_type;
const char *prc_type;
@@ -71,14 +63,47 @@ void PETScLinearSolver::solve(const PETScVector &b, PETScVector &x)
int its;
KSPGetIterationNumber(_solver, &its);
PetscPrintf(PETSC_COMM_WORLD,"\nConvergence in %d iterations.\n", its);
- PetscPrintf(PETSC_COMM_WORLD,"\n================================================");
+ PetscPrintf(PETSC_COMM_WORLD,"\n================================================\n");
+ }
+ else if(reason == KSP_DIVERGED_ITS)
+ {
+ PetscPrintf(PETSC_COMM_WORLD, "\nWaning: maximum iteration reached.\n");
+ converged = false;
+ }
+ else
+ {
+ if(reason == KSP_DIVERGED_INDEFINITE_PC)
+ {
+ PetscPrintf(PETSC_COMM_WORLD, "\nDivergence because of indefinite preconditioner,");
+ PetscPrintf(PETSC_COMM_WORLD, "\nTry to run again with -pc_factor_shift_positive_definite option.\n");
+ }
+ else if(reason == KSP_DIVERGED_BREAKDOWN_BICG)
+ {
+ PetscPrintf(PETSC_COMM_WORLD, "\nKSPBICG method was detected so the method could not continue to enlarge the Krylov space.");
+ PetscPrintf(PETSC_COMM_WORLD, "\nTry to run again with another solver.\n");
+ }
+
+ else if(reason == KSP_DIVERGED_NONSYMMETRIC)
+ {
+ PetscPrintf(PETSC_COMM_WORLD, "\nMatrx or preconditioner is unsymmetric but KSP requires symmetric.\n");
+ }
+ else
+ {
+ PetscPrintf(PETSC_COMM_WORLD, "\nOther kind divergence, use command option -ksp_monitor or -log_summary to check the details.\n");
+ }
+
+ PetscPrintf(PETSC_COMM_WORLD, "\nLinear solver (PETSc KSP) failed, quit now.\n");
+ KSPDestroy(&_solver);
+ PetscFinalize();
+ exit(EXIT_FAILURE);
}
- PetscPrintf(PETSC_COMM_WORLD,"\n");
#ifdef TEST_MEM_PETSC
PetscMemoryGetCurrentUsage(&mem2);
PetscPrintf(PETSC_COMM_WORLD, "###Memory usage by solver. Before :%f After:%f Increase:%d\n", mem1, mem2, (int)(mem2 - mem1));
#endif
+
+ return converged;
}
} //end of namespace
diff --git a/MathLib/LinAlg/PETSc/PETScLinearSolver.h b/MathLib/LinAlg/PETSc/PETScLinearSolver.h
index 430374c5bc364034f18c8733b4a29e5f2670a5a6..4e8aab6ced1a0b23f30b5e0658969203a784059d 100644
--- a/MathLib/LinAlg/PETSc/PETScLinearSolver.h
+++ b/MathLib/LinAlg/PETSc/PETScLinearSolver.h
@@ -54,16 +54,30 @@ class PETScLinearSolver
Solve a system of equations.
\param b The right hand of the equations.
\param x The solutions to be solved.
+ \return true: converged, false: diverged due to exeeding
+ the maximum iterations.
*/
- void solve(const PETScVector &b, PETScVector &x);
+ bool solve(const PETScVector &b, PETScVector &x);
+
+ /*!
+ \brief Get number of iterations.
+ If function solve(...) returns false, it returns the maximum
+ iterations.
+ */
+ PetscInt getIterations() const
+ {
+ PetscInt its = 0;
+ KSPGetIterationNumber(_solver, &its);
+ return its;
+ }
private:
/// Matrix, kept as a member only for solving successive linear equation
- /// that preconditioner matrix may vary.
- PETScMatrix &_A;
-
+ /// that preconditioner matrix may vary.
+ PETScMatrix &_A;
+
KSP _solver; ///< Slover type.
- PC _pc; ///< Preconditioner type.
+ PC _pc; ///< Preconditioner type.
};
} // end namespace
diff --git a/MathLib/LinAlg/PETSc/PETScTools.cpp b/MathLib/LinAlg/PETSc/PETScTools.cpp
index 3d9e0e46213826642dc81e707a2cd3d4e16fe53e..aea99c360604df6bd2d443733072e37f69b99010 100644
--- a/MathLib/LinAlg/PETSc/PETScTools.cpp
+++ b/MathLib/LinAlg/PETSc/PETScTools.cpp
@@ -28,10 +28,10 @@ void applyKnownSolution(PETScMatrix &A, PETScVector &b, PETScVector &x,
if(vec_knownX_id.size() > 0)
{
- x.set(vec_knownX_id, vec_knownX_x);
- b.set(vec_knownX_id, vec_knownX_x);
+ x.set(vec_knownX_id, vec_knownX_x);
+ b.set(vec_knownX_id, vec_knownX_x);
}
-
+
x.finalizeAssembly();
b.finalizeAssembly();
}
diff --git a/MathLib/LinAlg/PETSc/PETScTools.h b/MathLib/LinAlg/PETSc/PETScTools.h
index e417adbf1f30dfff23a1d4700bf66fb3f2cc68e0..aa9011719c2e01918368a4bfee54a7feb526ce3d 100644
--- a/MathLib/LinAlg/PETSc/PETScTools.h
+++ b/MathLib/LinAlg/PETSc/PETScTools.h
@@ -26,7 +26,7 @@ namespace MathLib
{
/*!
\brief apply known solutions to a system of linear equations
-
+
\param A Coefficient matrix
\param b RHS vector
\param vec_knownX_id a vector of known solution entry IDs
diff --git a/Tests/MathLib/TestLinearSolver.cpp b/Tests/MathLib/TestLinearSolver.cpp
index 0b5ff0d799de5223f57f317fdaf341b9d6949eb4..2881cc6ed3a79ed4d63b74e58c97543a406a1ede 100644
--- a/Tests/MathLib/TestLinearSolver.cpp
+++ b/Tests/MathLib/TestLinearSolver.cpp
@@ -155,7 +155,7 @@ void checkLinearSolverInterface(T_MATRIX &A, T_VECTOR &b, const std::string &pre
std::vector<double> local_vec(2);
local_vec[0] = mrank+1;
local_vec[1] = 2. * (mrank+1);
- x.set(col_pos, local_vec);
+ x.set(row_pos, local_vec);
double x0[6];
double x1[6];
@@ -164,7 +164,7 @@ void checkLinearSolverInterface(T_MATRIX &A, T_VECTOR &b, const std::string &pre
A.multiply(x, b);
// apply BC
- std::vector<int> bc_id; /// Type must be int to match Petsc_Int
+ std::vector<int> bc_id; // Type must be int to match Petsc_Int
std::vector<double> bc_value;
if(mrank == 1)
@@ -182,9 +182,11 @@ void checkLinearSolverInterface(T_MATRIX &A, T_VECTOR &b, const std::string &pre
// solve
T_LINEAR_SOVLER ls(A, prefix_name);
ls.solve(b, x);
+
+ EXPECT_GT(ls.getIterations(), 0u);
+
x.getGlobalVector(x1);
-
- ASSERT_ARRAY_NEAR(x0, x1, 6, 1e-5);
+ ASSERT_ARRAY_NEAR(x0, x1, 6, 1e-5);
}
#endif