From d5b8f549b9da1e13954ad1d293db61a7f75a3ae0 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Wed, 30 May 2012 15:23:33 +0100
Subject: [PATCH] fvMatrix: Further improvement to relaxation procedure

---
 .../fvMatrices/fvMatrix/fvMatrix.C            | 26 +++++++------------
 1 file changed, 10 insertions(+), 16 deletions(-)

diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
index 87a6647ca9..57fac754e9 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
@@ -570,24 +570,11 @@ void Foam::fvMatrix<Type>::relax(const scalar alpha)
             }
             else
             {
-                // For non-coupled boundaries subtract the diagonal
-                // contribution off-diagonal sum which avoids having to remove
-                // it from the diagonal later.
-                // Also add the source contribution from the relaxation
+                // For non-coupled boundaries add the maximum magnitude diagonal
+                // contribution to ensure stability
                 forAll(pa, face)
                 {
-                    // Type iCoeff0 = iCoeffs[face];
-                    // iCoeffs[face] = cmptMag(iCoeffs[face]);
-                    // sumOff[pa[face]] -= cmptMin(iCoeffs[face]);
-                    // iCoeffs[face] /= alpha;
-                    D[pa[face]] +=
-                        cmptMag(cmptMin(iCoeffs[face]))
-                      - cmptMin(iCoeffs[face]);
-                    sumOff[pa[face]] +=
-                        cmptMag(cmptMin(iCoeffs[face]))
-                      - cmptMin(iCoeffs[face]);
-                    // S[pa[face]] +=
-                    // cmptMultiply(iCoeffs[face] - iCoeff0, psi_[pa[face]]);
+                    D[pa[face]] += cmptMax(cmptMag(iCoeffs[face]));
                 }
             }
         }
@@ -653,6 +640,13 @@ void Foam::fvMatrix<Type>::relax(const scalar alpha)
                     D[pa[face]] -= component(iCoeffs[face], 0);
                 }
             }
+            else
+            {
+                forAll(pa, face)
+                {
+                    D[pa[face]] -= cmptMin(iCoeffs[face]);
+                }
+            }
         }
     }