From 2385dbff25c0770ad612c8ddd57aece5c619e1c8 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Thu, 17 Jan 2013 09:50:03 +0000
Subject: [PATCH] GAMG: Add optional correction interpolate after prolongation
---
src/OpenFOAM/Make/files | 3 +-
.../lduMatrix/solvers/GAMG/GAMGSolver.C | 2 +
.../lduMatrix/solvers/GAMG/GAMGSolver.H | 32 +++---
.../solvers/GAMG/GAMGSolverScalingFactor.C | 88 -----------------
.../lduMatrix/solvers/GAMG/GAMGSolverSolve.C | 99 ++++++++++---------
5 files changed, 76 insertions(+), 148 deletions(-)
delete mode 100644 src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverScalingFactor.C
diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index 03072a6f02c..dc8e946baaa 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -284,7 +284,8 @@ $(lduInterfaceFields)/cyclicLduInterfaceField/cyclicLduInterfaceField.C
GAMG = $(lduMatrix)/solvers/GAMG
$(GAMG)/GAMGSolver.C
$(GAMG)/GAMGSolverAgglomerateMatrix.C
-$(GAMG)/GAMGSolverScalingFactor.C
+$(GAMG)/GAMGSolverInterpolate.C
+$(GAMG)/GAMGSolverScale.C
$(GAMG)/GAMGSolverSolve.C
GAMGInterfaces = $(GAMG)/interfaces
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
index 7ae547c0bb3..a0ed11bb32b 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
@@ -71,6 +71,7 @@ Foam::GAMGSolver::GAMGSolver
postSweepsLevelMultiplier_(1),
maxPostSweeps_(4),
nFinestSweeps_(2),
+ interpolateCorrection_(false),
scaleCorrection_(matrix.symmetric()),
directSolveCoarsest_(false),
agglomeration_(GAMGAgglomeration::New(matrix_, controlDict_)),
@@ -174,6 +175,7 @@ void Foam::GAMGSolver::readControls()
);
controlDict_.readIfPresent("maxPostSweeps", maxPostSweeps_);
controlDict_.readIfPresent("nFinestSweeps", nFinestSweeps_);
+ controlDict_.readIfPresent("interpolateCorrection", interpolateCorrection_);
controlDict_.readIfPresent("scaleCorrection", scaleCorrection_);
controlDict_.readIfPresent("directSolveCoarsest", directSolveCoarsest_);
}
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H
index c2aefbea16a..44e27cc9fc1 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H
@@ -46,6 +46,8 @@ SourceFiles
GAMGSolverCalcAgglomeration.C
GAMGSolverMakeCoarseMatrix.C
GAMGSolverOperations.C
+ GAMGSolverInterpolate.C
+ GAMGSolverScale.C
GAMGSolverSolve.C
\*---------------------------------------------------------------------------*/
@@ -97,6 +99,10 @@ class GAMGSolver
//- Number of smoothing sweeps on finest mesh
label nFinestSweeps_;
+ //- Choose if the corrections should be interpolated after injection.
+ // By default corrections are not interpolated.
+ bool interpolateCorrection_;
+
//- Choose if the corrections should be scaled.
// By default corrections for symmetric matrices are scaled
// but not for asymmetric matrices.
@@ -154,32 +160,34 @@ class GAMGSolver
//- Agglomerate coarse matrix
void agglomerateMatrix(const label fineLevelIndex);
- //- Calculate and return the scaling factor from Acf, coarseSource
+ //- Interpolate the correction after injected prolongation
+ void interpolate
+ (
+ scalarField& psi,
+ scalarField& Apsi,
+ const lduMatrix& m,
+ const FieldField<Field, scalar>& interfaceBouCoeffs,
+ const lduInterfaceFieldPtrsList& interfaces,
+ const scalarField& source,
+ const direction cmpt
+ ) const;
+
+ //- Calculate and apply the scaling factor from Acf, coarseSource
// and coarseField.
// At the same time do a Jacobi iteration on the coarseField using
// the Acf provided after the coarseField values are used for the
// scaling factor.
- scalar scalingFactor
+ void scale
(
scalarField& field,
- const scalarField& source,
- const scalarField& Acf,
- const scalarField& D
- ) const;
-
- //- Calculate Acf and calculate and return the scaling factor.
- scalar scalingFactor
- (
scalarField& Acf,
const lduMatrix& A,
- scalarField& field,
const FieldField<Field, scalar>& interfaceLevelBouCoeffs,
const lduInterfaceFieldPtrsList& interfaceLevel,
const scalarField& source,
const direction cmpt
) const;
-
//- Initialise the data structures for the V-cycle
void initVcycle
(
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverScalingFactor.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverScalingFactor.C
deleted file mode 100644
index 1df912f32bd..00000000000
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverScalingFactor.C
+++ /dev/null
@@ -1,88 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "GAMGSolver.H"
-#include "vector2D.H"
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-Foam::scalar Foam::GAMGSolver::scalingFactor
-(
- scalarField& field,
- const scalarField& source,
- const scalarField& Acf,
- const scalarField& D
-) const
-{
- scalar scalingFactorNum = 0.0;
- scalar scalingFactorDenom = 0.0;
-
- forAll(field, i)
- {
- scalingFactorNum += source[i]*field[i];
- scalingFactorDenom += Acf[i]*field[i];
-
- // While the matrix-multiply done for the scaling it is
- // possible to perform a point-Jacobi smoothing operation cheaply
- field[i] += (source[i] - Acf[i])/D[i];
- }
-
- vector2D scalingVector(scalingFactorNum, scalingFactorDenom);
- reduce(scalingVector, sumOp<vector2D>());
- return scalingVector.x()/stabilise(scalingVector.y(), VSMALL);
-}
-
-
-Foam::scalar Foam::GAMGSolver::scalingFactor
-(
- scalarField& Acf,
- const lduMatrix& A,
- scalarField& field,
- const FieldField<Field, scalar>& interfaceLevelBouCoeffs,
- const lduInterfaceFieldPtrsList& interfaceLevel,
- const scalarField& source,
- const direction cmpt
-) const
-{
- A.Amul
- (
- Acf,
- field,
- interfaceLevelBouCoeffs,
- interfaceLevel,
- cmpt
- );
-
- return scalingFactor
- (
- field,
- source,
- Acf,
- A.diag()
- );
-}
-
-
-// ************************************************************************* //
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
index 25c8e6dbe56..ce2bf39f0cf 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
@@ -28,6 +28,7 @@ License
#include "BICCG.H"
#include "SubField.H"
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::solverPerformance Foam::GAMGSolver::solve
@@ -173,23 +174,16 @@ void Foam::GAMGSolver::Vcycle
// but not on the coarsest level because it evaluates to 1
if (scaleCorrection_ && leveli < coarsestLevel - 1)
{
- scalar sf = scalingFactor
+ scale
(
+ coarseCorrFields[leveli],
const_cast<scalarField&>(ACf.operator const scalarField&()),
matrixLevels_[leveli],
- coarseCorrFields[leveli],
interfaceLevelsBouCoeffs_[leveli],
interfaceLevels_[leveli],
coarseSources[leveli],
cmpt
);
-
- if (debug >= 2)
- {
- Pout<< sf << " ";
- }
-
- coarseCorrFields[leveli] *= sf;
}
// Correct the residual with the new solution
@@ -246,7 +240,7 @@ void Foam::GAMGSolver::Vcycle
coarseCorrFields[leveli].size()
);
- // Only store the preSmoothedCoarseCorrField is pre-smoothing is used
+ // Only store the preSmoothedCoarseCorrField if pre-smoothing is used
if (nPreSweeps_)
{
preSmoothedCoarseCorrField.assign(coarseCorrFields[leveli]);
@@ -259,38 +253,47 @@ void Foam::GAMGSolver::Vcycle
leveli + 1
);
- // Scale coarse-grid correction field
- // but not on the coarsest level because it evaluates to 1
- if (scaleCorrection_ && leveli < coarsestLevel - 1)
+ // Create A.psi for this coarse level as a sub-field of Apsi
+ scalarField::subField ACf
+ (
+ Apsi,
+ coarseCorrFields[leveli].size()
+ );
+
+ scalarField& ACfRef =
+ const_cast<scalarField&>(ACf.operator const scalarField&());
+
+ if (interpolateCorrection_)
{
- // Create A.psi for this coarse level as a sub-field of Apsi
- scalarField::subField ACf
+ interpolate
(
- Apsi,
- coarseCorrFields[leveli].size()
+ coarseCorrFields[leveli],
+ ACfRef,
+ matrixLevels_[leveli],
+ interfaceLevelsBouCoeffs_[leveli],
+ interfaceLevels_[leveli],
+ coarseSources[leveli],
+ cmpt
);
+ }
- scalar sf = scalingFactor
+ // Scale coarse-grid correction field
+ // but not on the coarsest level because it evaluates to 1
+ if (scaleCorrection_ && leveli < coarsestLevel - 1)
+ {
+ scale
(
- const_cast<scalarField&>(ACf.operator const scalarField&()),
- matrixLevels_[leveli],
coarseCorrFields[leveli],
+ ACfRef,
+ matrixLevels_[leveli],
interfaceLevelsBouCoeffs_[leveli],
interfaceLevels_[leveli],
coarseSources[leveli],
cmpt
);
-
-
- if (debug >= 2)
- {
- Pout<< sf << " ";
- }
-
- coarseCorrFields[leveli] *= sf;
}
- // Only add the preSmoothedCoarseCorrField is pre-smoothing is used
+ // Only add the preSmoothedCoarseCorrField if pre-smoothing is used
if (nPreSweeps_)
{
coarseCorrFields[leveli] += preSmoothedCoarseCorrField;
@@ -317,36 +320,38 @@ void Foam::GAMGSolver::Vcycle
0
);
- if (scaleCorrection_)
+ if (interpolateCorrection_)
{
- // Calculate finest level scaling factor
- scalar fsf = scalingFactor
+ interpolate
(
+ finestCorrection,
Apsi,
matrix_,
- finestCorrection,
interfaceBouCoeffs_,
interfaces_,
finestResidual,
cmpt
);
+ }
- if (debug >= 2)
- {
- Pout<< fsf << endl;
- }
-
- forAll(psi, i)
- {
- psi[i] += fsf*finestCorrection[i];
- }
+ if (scaleCorrection_)
+ {
+ // Scale the finest level correction
+ scale
+ (
+ finestCorrection,
+ Apsi,
+ matrix_,
+ interfaceBouCoeffs_,
+ interfaces_,
+ finestResidual,
+ cmpt
+ );
}
- else
+
+ forAll(psi, i)
{
- forAll(psi, i)
- {
- psi[i] += finestCorrection[i];
- }
+ psi[i] += finestCorrection[i];
}
smoothers[0].smooth
--
GitLab