Change Immersed Boundary dead cell handling in matrix: boost diagonal and setValues

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/basic/fixedValueIbFvPatchField/fixedValueIbFvPatchField.C

@@ -239,6 +239,16 @@ void fixedValueIbFvPatchField<Type>::evaluate
 }
 
 
+template<class Type>
+void Foam::fixedValueIbFvPatchField<Type>::manipulateMatrix
+(
+    fvMatrix<Type>& matrix
+)
+{
+    this->setDeadValues(matrix);
+}
+
+
 template<class Type>
 void fixedValueIbFvPatchField<Type>::write(Ostream& os) const
 {

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/basic/fixedValueIbFvPatchField/fixedValueIbFvPatchField.H

@@ -191,6 +191,9 @@ public:
                 const Pstream::commsTypes commsType = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvMatrix<Type>& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/basic/mixedIbFvPatchField/mixedIbFvPatchField.C

@@ -257,6 +257,16 @@ void mixedIbFvPatchField<Type>::evaluate
 }
 
 
+template<class Type>
+void Foam::mixedIbFvPatchField<Type>::manipulateMatrix
+(
+    fvMatrix<Type>& matrix
+)
+{
+    this->setDeadValues(matrix);
+}
+
+
 template<class Type>
 void mixedIbFvPatchField<Type>::write(Ostream& os) const
 {

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/basic/mixedIbFvPatchField/mixedIbFvPatchField.H

@@ -221,6 +221,9 @@ public:
                 const Pstream::commsTypes commsType = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvMatrix<Type>& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/derived/movingImmersedBoundaryVelocity/movingImmersedBoundaryVelocityFvPatchVectorField.C

@@ -220,6 +220,15 @@ void Foam::movingImmersedBoundaryVelocityFvPatchVectorField::evaluate
 }
 
 
+void Foam::movingImmersedBoundaryVelocityFvPatchVectorField::manipulateMatrix
+(
+    fvVectorMatrix& matrix
+)
+{
+    setDeadValues(matrix);
+}
+
+
 void Foam::movingImmersedBoundaryVelocityFvPatchVectorField::write
 (
     Ostream& os

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/derived/movingImmersedBoundaryVelocity/movingImmersedBoundaryVelocityFvPatchVectorField.H

@@ -163,6 +163,9 @@ public:
                 const Pstream::commsTypes commsType = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvVectorMatrix& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/immersedBoundaryFieldBase/immersedBoundaryFieldBase.C

@@ -63,6 +63,39 @@ void Foam::immersedBoundaryFieldBase<Type>::setDeadValues
 }
 
 
+template<class Type>
+void Foam::immersedBoundaryFieldBase<Type>::setDeadValues
+(
+    fvMatrix<Type>& matrix
+) const
+{
+    // Fix the value in dead cells
+    if (setDeadValue_)
+    {
+        const labelList& dc = ibPatch_.ibPolyPatch().deadCells();
+
+        // Boost the diagonal of dead cells by the volume ratio
+        // Volume ratio is set to SMALL; revert for diagonal
+        // This should also guarantee strong diagonal dominance.
+        // HJ, 19/Jun/2019
+
+        scalarField& diag = matrix.diag();
+
+        forAll (dc, dcI)
+        {
+            diag[dc[dcI]] *= GREAT;
+        }
+
+        // Set values
+        matrix.setValues
+        (
+            dc,
+            Field<Type>(dc.size(), deadValue_)
+        );
+    }
+}
+
+
 template<class Type>
 void Foam::immersedBoundaryFieldBase<Type>::writeDeadData(Ostream& os) const
 {

src/immersedBoundary/immersedBoundary/immersedBoundaryFvPatchFields/immersedBoundaryFieldBase/immersedBoundaryFieldBase.H

@@ -123,6 +123,9 @@ public:
             //- Set values in dead cells
             void setDeadValues(Field<Type>& psiI) const;
 
+            //- Set matrix constraints in dead cells
+            void setDeadValues(fvMatrix<Type>& matrix) const;
+
 
         // I/O
 

src/immersedBoundary/immersedBoundary/immersedBoundaryPolyPatch/immersedBoundaryPolyPatch.C

@@ -1262,8 +1262,8 @@ void Foam::immersedBoundaryPolyPatch::calcCorrectedGeometry() const
 
     forAll (dc, dcI)
     {
-        // Set dead volume to small
-        V[dc[dcI]] = SMALL;
+        // Scale dead volume to small
+        V[dc[dcI]] *= SMALL;
     }
 
     // Correct for all cut faces
@@ -1286,7 +1286,7 @@ void Foam::immersedBoundaryPolyPatch::calcCorrectedGeometry() const
 
     forAll (df, dfI)
     {
-        // Set dead area to small
+        // Scale dead area to small
         Sf[df[dfI]] *= SMALL;
     }
 

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryEpsilonWallFunction/immersedBoundaryEpsilonWallFunctionFvPatchScalarField.C

@@ -250,6 +250,17 @@ void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::evaluate
 }
 
 
+void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::manipulateMatrix
+(
+    fvScalarMatrix& matrix
+)
+{
+    setDeadValues(matrix);
+
+    epsilonWallFunctionFvPatchScalarField::manipulateMatrix(matrix);
+}
+
+
 void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::write
 (
     Ostream& os

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryEpsilonWallFunction/immersedBoundaryEpsilonWallFunctionFvPatchScalarField.H

@@ -170,6 +170,9 @@ public:
                 const Pstream::commsTypes  = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvScalarMatrix& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryKqRWallFunction/immersedBoundaryKqRWallFunctionFvPatchField.C

@@ -220,6 +220,16 @@ void immersedBoundaryKqRWallFunctionFvPatchField<Type>::evaluate
 }
 
 
+template<class Type>
+void immersedBoundaryKqRWallFunctionFvPatchField<Type>::manipulateMatrix
+(
+    fvMatrix<Type>& matrix
+)
+{
+    this->setDeadValues(matrix);
+}
+
+
 template<class Type>
 void immersedBoundaryKqRWallFunctionFvPatchField<Type>::write(Ostream& os) const
 {

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryKqRWallFunction/immersedBoundaryKqRWallFunctionFvPatchField.H

@@ -167,6 +167,9 @@ public:
                 const Pstream::commsTypes commsType = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvMatrix<Type>& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryOmegaWallFunctions/immersedBoundaryOmegaWallFunctionFvPatchScalarField.C

@@ -66,7 +66,7 @@ immersedBoundaryOmegaWallFunctionFvPatchScalarField
         p,
         Switch(dict.lookup("setDeadValue")),
         readScalar(dict.lookup("deadValue"))
-    )    
+    )
 {
     this->readPatchType(dict);
 
@@ -249,6 +249,17 @@ void immersedBoundaryOmegaWallFunctionFvPatchScalarField::evaluate
 }
 
 
+void immersedBoundaryOmegaWallFunctionFvPatchScalarField::manipulateMatrix
+(
+    fvScalarMatrix& matrix
+)
+{
+    setDeadValues(matrix);
+
+    omegaWallFunctionFvPatchScalarField::manipulateMatrix(matrix);
+}
+
+
 void immersedBoundaryOmegaWallFunctionFvPatchScalarField::write
 (
     Ostream& os

src/immersedBoundary/immersedBoundaryTurbulence/compressible/wallFunctions/immersedBoundaryOmegaWallFunctions/immersedBoundaryOmegaWallFunctionFvPatchScalarField.H

@@ -184,6 +184,9 @@ public:
                 const Pstream::commsTypes  = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvScalarMatrix& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryEpsilonWallFunction/immersedBoundaryEpsilonWallFunctionFvPatchScalarField.C

@@ -222,10 +222,6 @@ void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::updateCoeffs()
     {
         epsilonWallFunctionFvPatchScalarField::updateCoeffs();
     }
-
-    const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
-    const scalarField& y = rasModel.y()[patch().index()];
-    Info<< "Patch y: (" << min(y) << " " << max(y) << ")" << endl;
 }
 
 
@@ -254,6 +250,17 @@ void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::evaluate
 }
 
 
+void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::manipulateMatrix
+(
+    fvScalarMatrix& matrix
+)
+{
+    setDeadValues(matrix);
+
+    epsilonWallFunctionFvPatchScalarField::manipulateMatrix(matrix);
+}
+
+
 void immersedBoundaryEpsilonWallFunctionFvPatchScalarField::write
 (
     Ostream& os

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryEpsilonWallFunction/immersedBoundaryEpsilonWallFunctionFvPatchScalarField.H

@@ -170,6 +170,9 @@ public:
                 const Pstream::commsTypes  = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvScalarMatrix& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryKqRWallFunction/immersedBoundaryKqRWallFunctionFvPatchField.C

@@ -220,6 +220,16 @@ void immersedBoundaryKqRWallFunctionFvPatchField<Type>::evaluate
 }
 
 
+template<class Type>
+void immersedBoundaryKqRWallFunctionFvPatchField<Type>::manipulateMatrix
+(
+    fvMatrix<Type>& matrix
+)
+{
+    this->setDeadValues(matrix);
+}
+
+
 template<class Type>
 void immersedBoundaryKqRWallFunctionFvPatchField<Type>::write(Ostream& os) const
 {

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryKqRWallFunction/immersedBoundaryKqRWallFunctionFvPatchField.H

@@ -167,6 +167,9 @@ public:
                 const Pstream::commsTypes commsType = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvMatrix<Type>& matrix);
+
 
         // I-O
 

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryOmegaWallFunctions/immersedBoundaryOmegaWallFunctionFvPatchScalarField.C

@@ -66,7 +66,7 @@ immersedBoundaryOmegaWallFunctionFvPatchScalarField
         p,
         Switch(dict.lookup("setDeadValue")),
         readScalar(dict.lookup("deadValue"))
-    )    
+    )
 {
     this->readPatchType(dict);
 
@@ -249,6 +249,17 @@ void immersedBoundaryOmegaWallFunctionFvPatchScalarField::evaluate
 }
 
 
+void immersedBoundaryOmegaWallFunctionFvPatchScalarField::manipulateMatrix
+(
+    fvScalarMatrix& matrix
+)
+{
+    setDeadValues(matrix);
+
+    omegaWallFunctionFvPatchScalarField::manipulateMatrix(matrix);
+}
+
+
 void immersedBoundaryOmegaWallFunctionFvPatchScalarField::write
 (
     Ostream& os

src/immersedBoundary/immersedBoundaryTurbulence/incompressible/wallFunctions/immersedBoundaryOmegaWallFunctions/immersedBoundaryOmegaWallFunctionFvPatchScalarField.H

@@ -184,6 +184,9 @@ public:
                 const Pstream::commsTypes  = Pstream::blocking
             );
 
+            //- Manipulate a matrix
+            virtual void manipulateMatrix(fvScalarMatrix& matrix);
+
 
         // I-O