diff --git a/src/petsc4Foam/solvers/petscSolver.C b/src/petsc4Foam/solvers/petscSolver.C
index 56e9e5b7e620c27ffea3c3f52dc8f94f51fb5d07..4d9d81164c8d903d5c62f1aac885f284877c9f77 100644
--- a/src/petsc4Foam/solvers/petscSolver.C
+++ b/src/petsc4Foam/solvers/petscSolver.C
@@ -175,17 +175,38 @@ Foam::solverPerformance Foam::petscSolver::solve
updateKsp(ksp, Amat, !ctx.caching.needsPrecondUpdate());
}
- // Solver name from petsc
- KSPType ksptype;
- KSPGetType(ksp, &ksptype);
- const word solverName(ksptype);
-
// Use built-in residual norm computation
const bool usePetscResidualNorm
(
petscDict_.getOrDefault("use_petsc_residual_norm", false)
);
+ // Monitor built-in residual
+ const bool monitorFoamResidualNorm
+ (
+ petscDict_.getOrDefault("monitor_foam_residual_norm", false)
+ );
+
+
+ // Disable KSP default computation of residual norm if we are
+ // monitoring convergence a-la OpenFOAM
+ if (!usePetscResidualNorm)
+ {
+ KSPSetNormType(ksp, KSP_NORM_NONE);
+ }
+ else
+ {
+ KSPSetNormType(ksp, KSP_NORM_DEFAULT);
+ }
+
+ // ksp set options from db (may change norm type here if needed)
+ KSPSetFromOptions(ksp);
+
+ // Solver name from petsc (may have been changed from option database)
+ KSPType ksptype;
+ KSPGetType(ksp, &ksptype);
+ const word solverName(ksptype);
+
// Setup class containing solver performance data
solverPerformance solverPerf
(
@@ -200,59 +221,61 @@ Foam::solverPerformance Foam::petscSolver::solve
PetscWrappedVector petsc_psi(psi, Amat);
PetscWrappedVector petsc_source(source, Amat);
- // Initial residual (scaled)
+ // Initialize data to compute L1 residual or to monitor it while solving
+ if (!usePetscResidualNorm || monitorFoamResidualNorm)
+ {
+ ctx.createAuxVecs(petsc_psi);
+ ctx.initArowsSumVec(Amat);
+ ctx.computeNormFactor(Amat,petsc_psi,petsc_source);
+ }
+
+ // Convergence testing
if (usePetscResidualNorm)
{
+ // Add monitor to record initial residual
solverPerf.initialResidual() = 0;
KSPMonitorSet
(
ksp,
PetscUtils::foamKSPMonitorRecordInit,
&ctx,
- nullptr
+ NULL
);
}
else
{
- ctx.createAuxVecs(petsc_psi);
- ctx.initArowsSumVec(Amat);
- ctx.computeNormFactor(Amat,petsc_psi,petsc_source);
-
KSPSetConvergenceTest
(
ksp,
PetscUtils::foamKSPConverge,
&ctx,
- nullptr
- );
-
- // Monitor built-in residual
- const bool monitorFoamResidualNorm
- (
- petscDict_.getOrDefault("monitor_foam_residual_norm", false)
+ NULL
);
+ }
- if (monitorFoamResidualNorm)
- {
- KSPMonitorSet
- (
- ksp,
- PetscUtils::foamKSPMonitorFoam,
- &ctx,
- nullptr
- );
- }
+ // Monitoring (print to stdout) residual reduction in OpenFOAM norm
+ if (monitorFoamResidualNorm)
+ {
+ KSPMonitorSet
+ (
+ ksp,
+ PetscUtils::foamKSPMonitorFoam,
+ &ctx,
+ NULL
+ );
}
+ // Setup KSP (this is not explicitly needed, but we call it to separate
+ // PCSetUp from KSPSolve timings when requesting -log_view from PETSc)
+ KSPSetUp(ksp);
+
// Solve A x = b
KSPSolve(ksp, petsc_source, petsc_psi);
- // Set nIterations and compute final residual
- {
- PetscInt nIters;
- KSPGetIterationNumber(ksp, &nIters);
- solverPerf.nIterations() = nIters;
- }
+ // Set nIterations and final residual
+ PetscInt nIters;
+ KSPGetIterationNumber(ksp, &nIters);
+ solverPerf.nIterations() = nIters;
if (usePetscResidualNorm)
{
@@ -329,11 +352,6 @@ void Foam::petscSolver::buildKsp
// Use solution from the previous timestep as initial guess
KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
-
- // ksp set options from db
- KSPSetFromOptions(ksp);
-
- KSPSetUp(ksp);
}
diff --git a/src/petsc4Foam/utils/petscControls.C b/src/petsc4Foam/utils/petscControls.C
index 963fcbac2ebe483121ea893a4ab17ff82c8b1c87..ad5834d6b33940187bfc9b684ece6bb5ccdaed80 100644
--- a/src/petsc4Foam/utils/petscControls.C
+++ b/src/petsc4Foam/utils/petscControls.C
@@ -67,8 +67,10 @@ void Foam::petscControls::start(const fileName& optionsFile)
{
err = PetscInitialize
(
- nullptr, nullptr,
- optionsFile.c_str(), PETSC_NULL
+ NULL,
+ NULL,
+ optionsFile.c_str(),
+ NULL
);
}
else
diff --git a/src/petsc4Foam/utils/petscUtils.C b/src/petsc4Foam/utils/petscUtils.C
index 4dcd0af45057084a50870f199bfe7f059dac7ce0..8e27a4474020877fd746a80ae727ac894bd2da36 100644
--- a/src/petsc4Foam/utils/petscUtils.C
+++ b/src/petsc4Foam/utils/petscUtils.C
@@ -208,11 +208,32 @@ PetscErrorCode Foam::PetscUtils::foamKSPMonitorFoam
ierr = PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);CHKERRQ(ierr);
ierr = PetscViewerASCIIAddTab(viewer, tablevel);CHKERRQ(ierr);
ierr = KSPGetOptionsPrefix(ksp, &prefix);CHKERRQ(ierr);
- if (it == 0) {
- ierr = PetscViewerASCIIPrintf(viewer, " Residual norms for %s solve.\n", prefix);CHKERRQ(ierr);
+ if (it == 0)
+ {
+ ierr = PetscViewerASCIIPrintf(viewer, " Residual norms for %s solve.\n", prefix);CHKERRQ(ierr);
}
ierr = KSPGetNormType(ksp, &ntype);CHKERRQ(ierr);
- ierr = PetscViewerASCIIPrintf(viewer, "%3D KSP Residual foamnorm %14.12e (PETSc %s norm %14.12e)\n", it, (double)fnorm, KSPNormTypes[ntype], (double)rnorm);CHKERRQ(ierr);
+ if (ntype != KSP_NORM_NONE) // Print both norms, KSP built-in and OpenFOAM built-in
+ {
+ char normtype[256];
+ PetscErrorCode (*converge)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*);
+
+ ierr = PetscStrncpy(normtype,KSPNormTypes[ntype],sizeof(normtype));CHKERRQ(ierr);
+ ierr = PetscStrtolower(normtype);CHKERRQ(ierr);
+ ierr = KSPGetConvergenceTest(ksp, &converge, NULL, NULL);CHKERRQ(ierr);
+ if (converge == Foam::PetscUtils::foamKSPConverge) // we are using foam convergence testing, list foam norm first
+ {
+ ierr = PetscViewerASCIIPrintf(viewer, "%3D KSP Residual foam norm %14.12e (PETSc %s norm %14.12e)\n", it, (double)fnorm, normtype, (double)rnorm);CHKERRQ(ierr);
+ }
+ else // we are using KSP default convergence testing, list KSP norm first
+ {
+ ierr = PetscViewerASCIIPrintf(viewer, "%3D KSP Residual %s norm %14.12e (foam norm %14.12e)\n", it, normtype, (double)rnorm, (double)fnorm);CHKERRQ(ierr);
+ }
+ }
+ else // KSP has no norm, list foam norm only
+ {
+ ierr = PetscViewerASCIIPrintf(viewer, "%3D KSP Residual foam norm %14.12e\n", it, (double)fnorm);CHKERRQ(ierr);
+ }
ierr = PetscViewerASCIISubtractTab(viewer, tablevel);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
@@ -315,7 +336,7 @@ void Foam::PetscUtils::setFlag
Info<< key << ' ' << val << nl;
}
- PetscOptionsSetValue(nullptr, key.c_str(), val.c_str());
+ PetscOptionsSetValue(NULL, key.c_str(), val.c_str());
}
@@ -336,7 +357,7 @@ void Foam::PetscUtils::setFlags
Info<< key << ' ' << val << nl;
}
- PetscOptionsSetValue(nullptr, key.c_str(), val.c_str());
+ PetscOptionsSetValue(NULL, key.c_str(), val.c_str());
}
}
diff --git a/src/petsc4Foam/utils/petscWrappedVector.H b/src/petsc4Foam/utils/petscWrappedVector.H
index 3a9a4a6d1d51dd1ecf676e6e82914825a82f5b0d..97f12124f0b66a009690bf6faea07be868aff021 100644
--- a/src/petsc4Foam/utils/petscWrappedVector.H
+++ b/src/petsc4Foam/utils/petscWrappedVector.H
@@ -74,7 +74,7 @@ class PetscWrappedVector
Mat mat
)
{
- MatCreateVecs(mat, &vec_, PETSC_NULL);
+ MatCreateVecs(mat, &vec_, NULL);
VecPlaceArray
(
vec_,
diff --git a/tutorials/basic/laplacianFoam/flange/system/fvSolution-petsc b/tutorials/basic/laplacianFoam/flange/system/fvSolution-petsc
index aeb983b521503f5735b7382af598a752c190a3f1..dc373cb2fee90ff8b491e755f3dba77d54c2388c 100644
--- a/tutorials/basic/laplacianFoam/flange/system/fvSolution-petsc
+++ b/tutorials/basic/laplacianFoam/flange/system/fvSolution-petsc
@@ -30,10 +30,16 @@ solvers
pc_type bjacobi;
sub_pc_type icc;
ksp_cg_single_reduction true;
+ // uncomment the line below to obtain a textbook CG with use_petsc_residual_norm == true
+ // ksp_norm_type natural;
mat_type mpiaij;
}
+ // use petsc default convergence testing instead of OpenFOAM's L1 based criterion
+ // default is false
use_petsc_residual_norm false;
- monitor_foam_residual_norm true;
+ // monitor (print to stdout) residual reduction in OpenFOAM norm
+ // default is false
+ monitor_foam_residual_norm false;
}
tolerance 1e-06;