View Issue Details

IDProjectCategoryView StatusLast Update
0000270OpenFOAMBugpublic2013-12-31 17:18
Reporterjherb Assigned Touser4 
PrioritynormalSeverityfeatureReproducibilityalways
Status resolvedResolutionfixed 
PlatformLinuxOSRHESOS Version5.4
Summary0000270: Make normalization factor of residuals accessible to solver applications
DescriptionThe residuals used for the convergence checks are normalized by a factor defined in
https://github.com/OpenCFD/OpenFOAM-2.0.x/blob/master/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C#L179

The (not normalized) residuals can be accessed by
https://github.com/OpenCFD/OpenFOAM-2.0.x/blob/master/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C#L219

If the residuals should be saved to a file by the solver application, there is no way to normalize them, because the normalization factor cannot be accessed in the solver applications.

The attached solver makes the normFactor accessible to solver applications. Unfortunately I did not find a way to do this without solving the matrix, i.e. a call to solve().

Perhaps someone can suggest a different way.
TagsNo tags attached.

Activities

jherb

2011-08-05 14:36

reporter  

residual-patch-110805-hej.patch (13,100 bytes)   
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
index 881f88f..6c2db2c 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
@@ -96,6 +96,8 @@ public:
         label  noIterations_;
         bool   converged_;
         bool   singular_;
+        //- Normalization factor for residuals
+        scalar normFactor_;
 
 
     public:
@@ -109,7 +111,8 @@ public:
                 finalResidual_(0),
                 noIterations_(0),
                 converged_(false),
-                singular_(false)
+                singular_(false),
+                normFactor_(Foam::lduMatrix::small_)
             {}
 
             //- Construct from components
@@ -121,7 +124,8 @@ public:
                 const scalar fRes = 0,
                 const label  nIter = 0,
                 const bool   converged = false,
-                const bool   singular = false
+                const bool   singular = false,
+                const scalar normFactor = Foam::lduMatrix::small_
             )
             :
                 solverName_(solverName),
@@ -130,7 +134,8 @@ public:
                 finalResidual_(fRes),
                 noIterations_(nIter),
                 converged_(converged),
-                singular_(singular)
+                singular_(singular),
+                normFactor_(normFactor)
             {}
 
             //- Construct from Istream
@@ -210,6 +215,12 @@ public:
                 return singular_;
             }
 
+            //- Return final residual
+            scalar& normFactor()
+            {
+                return normFactor_;
+            }
+
             //- Convergence test
             bool checkConvergence
             (
@@ -413,9 +424,9 @@ public:
                 const scalarField& Apsi,
                 scalarField& tmpField
             ) const;
-    };
-
 
+    };
+    
     //- Abstract base-class for lduMatrix smoothers
     class smoother
     {
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
index 3e347e0..d93440e 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
@@ -62,6 +62,7 @@ Foam::lduMatrix::solverPerformance Foam::GAMGSolver::solve
     // Calculate normalised residual for convergence test
     solverPerf.initialResidual() = gSumMag(finestResidual)/normFactor;
     solverPerf.finalResidual() = solverPerf.initialResidual();
+    solverPerf.normFactor() = normFactor;
 
 
     // Check convergence, solve if not converged
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
index c433a8d..821c455 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
@@ -116,6 +116,7 @@ Foam::lduMatrix::solverPerformance Foam::PBiCG::solve
     // --- Calculate normalised residual norm
     solverPerf.initialResidual() = gSumMag(rA)/normFactor;
     solverPerf.finalResidual() = solverPerf.initialResidual();
+    solverPerf.normFactor() = normFactor;
 
     // --- Check convergence, solve if not converged
     if (!solverPerf.checkConvergence(tolerance_, relTol_))
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
index 965c072..6b57aea 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
@@ -107,6 +107,7 @@ Foam::lduMatrix::solverPerformance Foam::PCG::solve
     // --- Calculate normalised residual norm
     solverPerf.initialResidual() = gSumMag(rA)/normFactor;
     solverPerf.finalResidual() = solverPerf.initialResidual();
+    solverPerf.normFactor() = normFactor;
 
     // --- Check convergence, solve if not converged
     if (!solverPerf.checkConvergence(tolerance_, relTol_))
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
index a3d9a8f..1412d8e 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
@@ -124,6 +124,7 @@ Foam::lduMatrix::solverPerformance Foam::smoothSolver::solve
             // Calculate residual magnitude
             solverPerf.initialResidual() = gSumMag(source - Apsi)/normFactor;
             solverPerf.finalResidual() = solverPerf.initialResidual();
+            solverPerf.normFactor() = normFactor;
         }
 
         if (lduMatrix::debug >= 2)
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
index d5086f5..017bdb7 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
@@ -1345,6 +1345,17 @@ Foam::lduMatrix::solverPerformance Foam::solve
 template<class Type>
 Foam::lduMatrix::solverPerformance Foam::solve
 (
+    fvMatrix<Type>& fvm,
+    const dictionary& solverControls,
+    Type& normFactor
+)
+{
+    return fvm.solve(solverControls, normFactor);
+}
+
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::solve
+(
     const tmp<fvMatrix<Type> >& tfvm,
     const dictionary& solverControls
 )
@@ -1357,6 +1368,21 @@ Foam::lduMatrix::solverPerformance Foam::solve
     return solverPerf;
 }
 
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::solve
+(
+    const tmp<fvMatrix<Type> >& tfvm,
+    const dictionary& solverControls,
+    Type& normFactor
+)
+{
+    lduMatrix::solverPerformance solverPerf =
+        const_cast<fvMatrix<Type>&>(tfvm()).solve(solverControls, normFactor);
+
+    tfvm.clear();
+
+    return solverPerf;
+}
 
 template<class Type>
 Foam::lduMatrix::solverPerformance Foam::solve(fvMatrix<Type>& fvm)
@@ -1365,6 +1391,12 @@ Foam::lduMatrix::solverPerformance Foam::solve(fvMatrix<Type>& fvm)
 }
 
 template<class Type>
+Foam::lduMatrix::solverPerformance Foam::solve(fvMatrix<Type>& fvm, Type& normFactor)
+{
+    return fvm.solve(normFactor);
+}
+
+template<class Type>
 Foam::lduMatrix::solverPerformance Foam::solve(const tmp<fvMatrix<Type> >& tfvm)
 {
     lduMatrix::solverPerformance solverPerf =
@@ -1375,6 +1407,16 @@ Foam::lduMatrix::solverPerformance Foam::solve(const tmp<fvMatrix<Type> >& tfvm)
     return solverPerf;
 }
 
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::solve(const tmp<fvMatrix<Type> >& tfvm, Type& normFactor)
+{
+    lduMatrix::solverPerformance solverPerf =
+        const_cast<fvMatrix<Type>&>(tfvm()).solve(normFactor);
+
+    tfvm.clear();
+
+    return solverPerf;
+}
 
 template<class Type>
 Foam::tmp<Foam::fvMatrix<Type> > Foam::correction
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
index 4763714..0946f09 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
@@ -142,6 +142,9 @@ class fvMatrix
         mutable GeometricField<Type, fvsPatchField, surfaceMesh>
             *faceFluxCorrectionPtr_;
 
+        // Normalization factors for residuals
+        Type normFactor_;
+
 protected:
 
     //- Declare friendship with the fvSolver class
@@ -240,9 +243,13 @@ public:
             //  Use the given solver controls
             lduMatrix::solverPerformance solve(const dictionary&);
 
+            lduMatrix::solverPerformance solve(const dictionary&, Type& normFactor);
+
             //- Solve returning the solution statistics.
             //  Solver controls read from fvSolution
             lduMatrix::solverPerformance solve();
+
+            lduMatrix::solverPerformance solve(Type& normFactor);
     };
 
 
@@ -323,6 +330,11 @@ public:
                 return faceFluxCorrectionPtr_;
             }
 
+            Type normFactor() const
+            {
+                return normFactor_;
+            }
+
 
         // Operations
 
@@ -385,14 +397,20 @@ public:
             //  Solver controls read from fvSolution
             autoPtr<fvSolver> solver();
 
+            autoPtr<fvSolver> solver(Type& normFactor);
+
             //- Solve returning the solution statistics.
             //  Use the given solver controls
             lduMatrix::solverPerformance solve(const dictionary&);
 
+            lduMatrix::solverPerformance solve(const dictionary&, Type& normFactor);
+
             //- Solve returning the solution statistics.
             //  Solver controls read from fvSolution
             lduMatrix::solverPerformance solve();
 
+            lduMatrix::solverPerformance solve(Type& normFactor);
+
             //- Return the matrix residual
             tmp<Field<Type> > residual() const;
 
@@ -559,6 +577,8 @@ lduMatrix::solverPerformance solve
 template<class Type>
 lduMatrix::solverPerformance solve(fvMatrix<Type>&);
 
+template<class Type>
+lduMatrix::solverPerformance solve(fvMatrix<Type>&, Type& normFactor);
 
 //- Solve returning the solution statistics given convergence tolerance,
 //  deleting temporary matrix after solution.
@@ -572,7 +592,6 @@ lduMatrix::solverPerformance solve(const tmp<fvMatrix<Type> >&);
 template<class Type>
 tmp<fvMatrix<Type> > correction(const fvMatrix<Type>&);
 
-
 //- Return the correction form of the given temporary matrix
 //  by subtracting the matrix multiplied by the current field
 template<class Type>
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
index 0b1f510..80df261 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
@@ -152,6 +152,9 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve
         solverPerfVec = max(solverPerfVec, solverPerf);
         solverPerfVec.solverName() = solverPerf.solverName();
 
+        // Save normalization factor for residuals
+        normFactor_.component(cmpt) = solverPerf.normFactor();
+
         psi.internalField().replace(cmpt, psiCmpt);
         diag() = saveDiag;
     }
@@ -163,6 +166,18 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve
     return solverPerfVec;
 }
 
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve
+(
+    const dictionary& solverControls,
+    Type& normFactor
+)
+{
+    lduMatrix::solverPerformance solverPerfVec = solve(solverControls);
+    normFactor = normFactor_;
+    return solverPerfVec;
+}
+
 
 template<class Type>
 Foam::autoPtr<typename Foam::fvMatrix<Type>::fvSolver>
@@ -181,6 +196,23 @@ Foam::fvMatrix<Type>::solver()
     );
 }
 
+template<class Type>
+Foam::autoPtr<typename Foam::fvMatrix<Type>::fvSolver>
+Foam::fvMatrix<Type>::solver(Type& normFactor)
+{
+    return solver
+    (
+        psi_.mesh().solverDict
+        (
+            psi_.select
+            (
+                psi_.mesh().data::template lookupOrDefault<bool>
+                ("finalIteration", false)
+            )
+        ),
+        normFactor
+    );
+}
 
 template<class Type>
 Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::fvSolver::solve()
@@ -198,6 +230,22 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::fvSolver::solve()
     );
 }
 
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::fvSolver::solve(Type& normFactor)
+{
+    return solve
+    (
+        fvMat_.psi_.mesh().solverDict
+        (
+            fvMat_.psi_.select
+            (
+                fvMat_.psi_.mesh().data::template lookupOrDefault<bool>
+                ("finalIteration", false)
+            )
+        ),
+        normFactor
+    );
+}
 
 template<class Type>
 Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve()
@@ -215,6 +263,22 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve()
     );
 }
 
+template<class Type>
+Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve(Type& normFactor)
+{
+    return solve
+    (
+        psi_.mesh().solverDict
+        (
+            psi_.select
+            (
+                psi_.mesh().data::template lookupOrDefault<bool>
+                ("finalIteration", false)
+            )
+        ),
+        normFactor
+    );
+}
 
 template<class Type>
 Foam::tmp<Foam::Field<Type> > Foam::fvMatrix<Type>::residual() const
diff --git a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
index c17cf14..91c556b 100644
--- a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
@@ -170,6 +170,9 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Foam::scalar>::solve
 
     psi.mesh().setSolverPerformance(psi.name(), solverPerf);
 
+    // Calculate normalization factor for residuals
+    normFactor_ = solverPerf.normFactor();
+
     return solverPerf;
 }
 
residual-patch-110805-hej.patch (13,100 bytes)   

jherb

2011-11-25 18:04

reporter   ~0000824

The links to the source have changed due to the new repository names to:
https://github.com/OpenFOAM/OpenFOAM-2.0.x/blob/master/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C#L179
https://github.com/OpenFOAM/OpenFOAM-2.0.x/blob/master/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C#L219

user4

2013-12-31 17:18

  ~0002728

Fixed through the data objectRegistry (src/OpenFOAM/lnInclude/data.H). See e.g. src/finiteVolume/lnInclude/simpleControl.C for use of residuals.

Issue History

Date Modified Username Field Change
2011-08-05 14:36 jherb New Issue
2011-08-05 14:36 jherb File Added: residual-patch-110805-hej.patch
2011-11-25 18:04 jherb Note Added: 0000824
2013-12-31 17:18 user4 Note Added: 0002728
2013-12-31 17:18 user4 Status new => resolved
2013-12-31 17:18 user4 Fixed in Version => Other
2013-12-31 17:18 user4 Resolution open => fixed
2013-12-31 17:18 user4 Assigned To => user4