diff --git a/src/petsc4Foam/solvers/petscSolver.C b/src/petsc4Foam/solvers/petscSolver.C
index 6718ff6deb026573c3e276bd26809c2563e347a9..690d45a4405f5f421744af83182795a23e2784bc 100644
--- a/src/petsc4Foam/solvers/petscSolver.C
+++ b/src/petsc4Foam/solvers/petscSolver.C
@@ -157,14 +157,11 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
List<label>& lowNonZero = ctx.lowNonZero;
label& maxLowNonZeroPerRow = ctx.maxLowNonZeroPerRow;
- bool needsCacheUpdate = true;
-
if (!ctx.initialized())
{
DebugInfo<< "Initializing PETSc Linear Solver " << eqName_ << nl;
ctx.initialized() = true;
- needsCacheUpdate = false;
PetscLogStageRegister
(
@@ -172,20 +169,30 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
&ctx.matstage
);
PetscLogStageRegister
+ (
+ ("foam_" + eqName_ + "_pc").c_str(),
+ &ctx.pcstage
+ );
+ PetscLogStageRegister
(
("foam_" + eqName_ + "_ksp").c_str(),
&ctx.kspstage
);
- computeMatAllocation(Amat, lowNonZero, maxLowNonZeroPerRow);
- buildMat(Amat, lowNonZero, maxLowNonZeroPerRow, ctx.matstage);
+ PetscLogStagePush(ctx.matstage);
+ buildMat(Amat, lowNonZero, maxLowNonZeroPerRow);
+ PetscLogStagePop();
buildKsp(Amat, ksp);
}
- if (ctx.caching.needsMatrixUpdate() && needsCacheUpdate)
+ bool matup = ctx.caching.needsMatrixUpdate();
+ bool pcup = ctx.caching.needsPrecondUpdate();
+ if (matup)
{
- buildMat(Amat, lowNonZero, maxLowNonZeroPerRow, ctx.matstage);
- updateKsp(ksp, Amat, !ctx.caching.needsPrecondUpdate());
+ PetscLogStagePush(ctx.matstage);
+ updateMat(Amat, lowNonZero, maxLowNonZeroPerRow);
+ PetscLogStagePop();
}
+ updateKsp(ksp, Amat, !pcup);
// Use built-in residual norm computation
const bool usePetscResidualNorm
@@ -273,24 +280,24 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
// Monitoring (print to stdout) residual reduction in OpenFOAM norm
if (monitorFoamResidualNorm)
{
- KSPMonitorSet
- (
- ksp,
- PetscUtils::foamKSPMonitorFoam,
- &ctx,
- NULL
- );
+ KSPMonitorSet
+ (
+ ksp,
+ PetscUtils::foamKSPMonitorFoam,
+ &ctx,
+ NULL
+ );
}
- PetscLogStagePush(ctx.kspstage);
-
// Setup KSP (this is not explicitly needed, but we call it to separate
// PCSetUp from KSPSolve timings when requesting -log_view from PETSc)
+ PetscLogStagePush(ctx.pcstage);
KSPSetUp(ksp);
+ PetscLogStagePop();
// Solve A x = b
+ PetscLogStagePush(ctx.kspstage);
KSPSolve(ksp, petsc_source, petsc_psi);
-
PetscLogStagePop();
// Set nIterations and final residual
@@ -336,15 +343,15 @@ void Foam::petscSolver::updateKsp
const bool reuse
) const
{
- DebugInfo<< "Cache-Hit: updating PETSc Amat " << eqName_ << nl;
-
if (reuse)
{
+ DebugInfo<< "Cache-Hit: reuse preconditioner " << eqName_ << nl;
+
KSPSetReusePreconditioner(ksp, PETSC_TRUE);
}
else
{
- DebugInfo<< "Cache-Hit: updating PETSc Pmat " << eqName_ << nl;
+ DebugInfo<< "Cache-Hit: rebuild preconditioner " << eqName_ << nl;
KSPSetReusePreconditioner(ksp, PETSC_FALSE);
}
@@ -393,7 +400,7 @@ void Foam::petscSolver::buildKsp
}
-void Foam::petscSolver::computeMatAllocation
+void Foam::petscSolver::buildMat
(
Mat& Amat,
List<label>& lowNonZero,
@@ -505,7 +512,7 @@ void Foam::petscSolver::computeMatAllocation
nonZeroNeigh_.data()
);
- // The general AIJ preallocator do not handle SELL
+ // The general AIJ preallocator does not handle SELL
MatSeqSELLSetPreallocation
(
Amat,
@@ -537,16 +544,13 @@ void Foam::petscSolver::computeMatAllocation
}
-void Foam::petscSolver::buildMat
+void Foam::petscSolver::updateMat
(
Mat& Amat,
List<label>& lowNonZero,
- label& maxLowNonZeroPerRow,
- PetscLogStage stage
+ label& maxLowNonZeroPerRow
) const
{
- PetscLogStagePush(stage);
-
const lduAddressing& lduAddr = matrix_.mesh().lduAddr();
const lduInterfacePtrsList interfaces(matrix_.mesh().interfaces());
@@ -727,7 +731,6 @@ void Foam::petscSolver::buildMat
MatAssemblyBegin(Amat, MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(Amat, MAT_FINAL_ASSEMBLY);
- PetscLogStagePop();
}
diff --git a/src/petsc4Foam/solvers/petscSolver.H b/src/petsc4Foam/solvers/petscSolver.H
index d0287a55c4e463661fa72e2aa1c099238bc088fa..983ba27b452a152d87315e0701b743c8fb44c1ac 100644
--- a/src/petsc4Foam/solvers/petscSolver.H
+++ b/src/petsc4Foam/solvers/petscSolver.H
@@ -93,20 +93,19 @@ class petscSolver
//- Compute allocation (proc and off-proc)
//- for the petsc matrix Amat_
- void computeMatAllocation
+ void buildMat
(
Mat& Amat,
List<label>& lowNonZero,
label& maxLowNonZeroPerRow
) const;
- //- Build Amat_ matrix by inserting new values
- void buildMat
+ //- Update Amat_ matrix by inserting new values
+ void updateMat
(
Mat& Amat,
List<label>& lowNonZero,
- label& maxLowNonZeroPerRow,
- PetscLogStage stage
+ label& maxLowNonZeroPerRow
) const;
//- Build the ksp_ krylov solver context
diff --git a/src/petsc4Foam/utils/petscCacheManager.H b/src/petsc4Foam/utils/petscCacheManager.H
index 701562822db3eea81bdad4048be4bf10fc30cfa4..f6d5cdd05305c023d3c4225b09328d35a76debb7 100644
--- a/src/petsc4Foam/utils/petscCacheManager.H
+++ b/src/petsc4Foam/utils/petscCacheManager.H
@@ -162,8 +162,9 @@ class petscCacheManager
{
// Private Data
- //- The current (relative) iteration
- label iter;
+ //- The current (relative) iterations for matrix and preconditioner
+ label miter;
+ label piter;
PetscUtils::Caching matrixCaching;
PetscUtils::Caching precondCaching;
@@ -176,7 +177,8 @@ public:
//- Default construct
petscCacheManager()
:
- iter{0}
+ miter{0},
+ piter{0}
{}
@@ -195,12 +197,12 @@ public:
bool needsMatrixUpdate()
{
- return needsUpdate(matrixCaching, iter);
+ return needsUpdate(matrixCaching, miter);
}
bool needsPrecondUpdate()
{
- return needsUpdate(precondCaching, iter);
+ return needsUpdate(precondCaching, piter);
}
void eventBegin()
@@ -219,11 +221,11 @@ public:
// Elapsed timing interval (s)
precondCaching.timeIter = precondCaching.timer_.elapsed();
- if (iter == 0)
+ if (piter == 0)
{
precondCaching.time0 = precondCaching.timeIter;
}
- else if (iter == 1)
+ else if (piter == 1)
{
precondCaching.time1 = precondCaching.timeIter;
}
@@ -235,8 +237,16 @@ public:
|| precondCaching.updateType_ == PetscUtils::Caching::Adaptive
)
{
- ++iter;
+ ++piter;
}
+ if
+ (
+ matrixCaching.updateType_ == PetscUtils::Caching::Periodic
+ )
+ {
+ ++miter;
+ }
+
}
@@ -276,14 +286,14 @@ private:
case PetscUtils::Caching::Adaptive:
{
- need = false;
if (iter > 3) // we need at least three times
{
+ need = false;
const double ratio0 =
- precondCaching.time0 / precondCaching.timeIter;
+ caching.time0 / caching.timeIter;
const double ratio1 =
- precondCaching.time1 / precondCaching.timeIter;
+ caching.time1 / caching.timeIter;
const int nsteps =
min(1e5, ratio0 * (1. / mag(1. - ratio1 + 1e-6)));
diff --git a/src/petsc4Foam/utils/petscLinearSolverContext.H b/src/petsc4Foam/utils/petscLinearSolverContext.H
index 7eb25b7e040e244509f701d7f3c36883241ca5fa..21b4873ebe85d82277923da8861daf13e761c35d 100644
--- a/src/petsc4Foam/utils/petscLinearSolverContext.H
+++ b/src/petsc4Foam/utils/petscLinearSolverContext.H
@@ -70,6 +70,7 @@ public:
KSP ksp;
PetscLogStage matstage;
+ PetscLogStage pcstage;
PetscLogStage kspstage;
solverPerformance performance;
@@ -92,23 +93,21 @@ public:
:
init_(false),
init_aux_vectors_(false)
- {}
+ {
+ Amat = NULL;
+ ksp = NULL;
+ ArowsSum = NULL;
+ res_l1_w = NULL;
+ }
//- Destructor
virtual ~petscLinearSolverContext()
{
- if (init_)
- {
- MatDestroy(&Amat);
- KSPDestroy(&ksp);
- }
-
- if (init_aux_vectors_)
- {
- VecDestroy(&ArowsSum);
- VecDestroyVecs(2,&res_l1_w);
- }
+ MatDestroy(&Amat);
+ KSPDestroy(&ksp);
+ VecDestroy(&ArowsSum);
+ if (res_l1_w) VecDestroyVecs(2,&res_l1_w);
}