From 218abc718515722773aefd29b88e0328a394cda9 Mon Sep 17 00:00:00 2001
From: Stefano Zampini <stefano.zampini@gmail.com>
Date: Thu, 8 Oct 2020 22:13:58 +0300
Subject: [PATCH] Preliminary support for geometric multigrid from PC_ML
---
src/petsc4Foam/solvers/petscSolver.C | 36 +++++++++++++++++++++++++++-
1 file changed, 35 insertions(+), 1 deletion(-)
diff --git a/src/petsc4Foam/solvers/petscSolver.C b/src/petsc4Foam/solvers/petscSolver.C
index 77fb3d6..900c96f 100644
--- a/src/petsc4Foam/solvers/petscSolver.C
+++ b/src/petsc4Foam/solvers/petscSolver.C
@@ -146,6 +146,7 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
petscLinearSolverContexts::New(fvm);
petscLinearSolverContext& ctx = contexts.getContext(eqName_);
+ const bool firsttimein = !ctx.initialized();
dictionary petscDictCaching = petscDict_.subOrEmptyDict("caching");
ctx.caching.init(petscDictCaching);
@@ -157,7 +158,7 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
List<label>& lowNonZero = ctx.lowNonZero;
label& maxLowNonZeroPerRow = ctx.maxLowNonZeroPerRow;
- if (!ctx.initialized())
+ if (firsttimein)
{
DebugInfo<< "Initializing PETSc Linear Solver " << eqName_ << nl;
@@ -289,6 +290,39 @@ Foam::solverPerformance Foam::petscSolver::scalarSolve
);
}
+ // Geometric multigrid support
+ // (pc_type ml or pc_type gamg pc_gamg_type geo for 2d only as for 3.14)
+ if (firsttimein)
+ {
+ PC pc;
+ void (*f)(void) = NULL;
+
+ KSPGetPC(ksp, &pc);
+ PetscObjectQueryFunction((PetscObject)pc, "PCSetCoordinates_C", &f);
+
+ // Cell-centres are contiguous in memory (array of structures)
+ // Not sure about const-ness here
+ // Could also use the PrecisionAdaptor
+ if (f)
+ {
+ const vectorField& cc = fvm.C().primitiveField();
+ const PetscInt n = cc.size();
+ const PetscInt sdim = vector::nComponents;
+
+ List<PetscReal> ccPoints(cc.size()*vector::nComponents);
+
+ auto iter = ccPoints.data();
+ for (const vector& v : cc)
+ {
+ *(iter++) = v.x();
+ *(iter++) = v.y();
+ *(iter++) = v.z();
+ }
+
+ PCSetCoordinates(pc, sdim, n, ccPoints.data());
+ }
+ }
+
// 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);
--
GitLab