diff --git a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrix.H b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrix.H
index bae2bd10d5901cc78b91e6c5b9d7156074bd509d..c0673eb34f7878e715df9c913702dea7bf39f43d 100644
--- a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrix.H
+++ b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrix.H
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -134,6 +134,9 @@ public:
             //- Minimum number of iterations in the solver
             label minIter_;
 
+            //- Level of verbosity in the solver output statements
+            label log_;
+
             //- Final convergence tolerance
             Type tolerance_;
 
diff --git a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrixSolver.C b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrixSolver.C
index 26e0dab299c2f149b6b1018cbc55d5735b245ac2..54a5923603c457f2dfdebc7549112719a4ccdded 100644
--- a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrixSolver.C
+++ b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/LduMatrixSolver.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -131,6 +131,7 @@ Foam::LduMatrix<Type, DType, LUType>::solver::solver
 
     maxIter_(defaultMaxIter_),
     minIter_(0),
+    log_(1),
     tolerance_(1e-6*pTraits<Type>::one),
     relTol_(Zero)
 {
@@ -145,6 +146,7 @@ void Foam::LduMatrix<Type, DType, LUType>::solver::readControls()
 {
     readControl(controlDict_, maxIter_, "maxIter");
     readControl(controlDict_, minIter_, "minIter");
+    readControl(controlDict_, log_, "log");
     readControl(controlDict_, tolerance_, "tolerance");
     readControl(controlDict_, relTol_, "relTol");
 }
diff --git a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C
index e291b0289ab718a256635e3997f550790fc9862c..2c53aedbee84fe524e2c4939e137d8e133a64035 100644
--- a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C
+++ b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -62,10 +63,11 @@ template<class Type>
 bool Foam::SolverPerformance<Type>::checkConvergence
 (
     const Type& Tolerance,
-    const Type& RelTolerance
+    const Type& RelTolerance,
+    const label log
 )
 {
-    if (debug >= 2)
+    if ((log >= 2) || (debug >= 2))
     {
         Info<< solverName_
             << ":  Iteration " << nIterations_
diff --git a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.H b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.H
index 9d835599ce525408836b4784b3ad575b7f21f507..476e842cb3c08942ebcbde2f165e5d3a0557e13b 100644
--- a/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.H
+++ b/src/OpenFOAM/matrices/LduMatrix/LduMatrix/SolverPerformance.H
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2012-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -218,7 +219,8 @@ public:
         bool checkConvergence
         (
             const Type& tolerance,
-            const Type& relTolerance
+            const Type& relTolerance,
+            const label log
         );
 
         //- Singularity test
diff --git a/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCCCG/PBiCCCG.C b/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCCCG/PBiCCCG.C
index f20a64a02600c794f615c847f2f7dda85c7beba6..c219c276a60b9d6a5303e9e7624d580757aa3395 100644
--- a/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCCCG/PBiCCCG.C
+++ b/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCCCG/PBiCCCG.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -98,7 +99,7 @@ Foam::PBiCCCG<Type, DType, LUType>::solve
     // --- Calculate normalisation factor
     Type normFactor = this->normFactor(psi, wA, pA);
 
-    if (LduMatrix<Type, DType, LUType>::debug >= 2)
+    if ((this->log_ >= 2) || (LduMatrix<Type, DType, LUType>::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -111,7 +112,12 @@ Foam::PBiCCCG<Type, DType, LUType>::solve
     if
     (
         this->minIter_ > 0
-     || !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+     || !solverPerf.checkConvergence
+        (
+            this->tolerance_,
+            this->relTol_,
+            this->log_
+        )
     )
     {
         // --- Select and construct the preconditioner
@@ -192,7 +198,12 @@ Foam::PBiCCCG<Type, DType, LUType>::solve
         (
             (
                 nIter++ < this->maxIter_
-            && !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+            && !solverPerf.checkConvergence
+                (
+                    this->tolerance_,
+                    this->relTol_,
+                    this->log_
+                )
             )
          || nIter < this->minIter_
         );
diff --git a/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCICG/PBiCICG.C b/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCICG/PBiCICG.C
index 5796597fac9e3e549528d20dab642c750b5d84b7..6dee6ad57ee8dc4f3392d7bee810209f09b472db 100644
--- a/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCICG/PBiCICG.C
+++ b/src/OpenFOAM/matrices/LduMatrix/Solvers/PBiCICG/PBiCICG.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -95,7 +96,7 @@ Foam::PBiCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
     // --- Calculate normalisation factor
     Type normFactor = this->normFactor(psi, wA, pA);
 
-    if (LduMatrix<Type, DType, LUType>::debug >= 2)
+    if ((this->log_ >= 2) || (LduMatrix<Type, DType, LUType>::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -105,7 +106,15 @@ Foam::PBiCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
     solverPerf.finalResidual() = solverPerf.initialResidual();
 
     // --- Check convergence, solve if not converged
-    if (!solverPerf.checkConvergence(this->tolerance_, this->relTol_))
+    if
+    (
+        !solverPerf.checkConvergence
+        (
+            this->tolerance_,
+            this->relTol_,
+            this->log_
+        )
+    )
     {
         // --- Select and construct the preconditioner
         autoPtr<typename LduMatrix<Type, DType, LUType>::preconditioner>
@@ -192,7 +201,12 @@ Foam::PBiCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
         } while
         (
             nIter++ < this->maxIter_
-        && !(solverPerf.checkConvergence(this->tolerance_, this->relTol_))
+        && !solverPerf.checkConvergence
+            (
+                this->tolerance_,
+                this->relTol_,
+                this->log_
+            )
         );
     }
 
diff --git a/src/OpenFOAM/matrices/LduMatrix/Solvers/PCICG/PCICG.C b/src/OpenFOAM/matrices/LduMatrix/Solvers/PCICG/PCICG.C
index 7adf70b392f88515b18dc04d975f887abd44e4de..6dec5d953d32225bf4644879aec3ea88d2da3b77 100644
--- a/src/OpenFOAM/matrices/LduMatrix/Solvers/PCICG/PCICG.C
+++ b/src/OpenFOAM/matrices/LduMatrix/Solvers/PCICG/PCICG.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -86,7 +87,7 @@ Foam::PCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
     // --- Calculate normalisation factor
     Type normFactor = this->normFactor(psi, wA, pA);
 
-    if (LduMatrix<Type, DType, LUType>::debug >= 2)
+    if ((this->log_ >= 2) || (LduMatrix<Type, DType, LUType>::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -99,7 +100,12 @@ Foam::PCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
     if
     (
         this->minIter_ > 0
-     || !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+     || !solverPerf.checkConvergence
+        (
+            this->tolerance_,
+            this->relTol_,
+            this->log_
+        )
     )
     {
         // --- Select and construct the preconditioner
@@ -184,7 +190,12 @@ Foam::PCICG<Type, DType, LUType>::solve(Field<Type>& psi) const
         (
             (
                 nIter++ < this->maxIter_
-            && !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+            && !solverPerf.checkConvergence
+                (
+                    this->tolerance_,
+                    this->relTol_,
+                    this->log_
+                )
             )
          || nIter < this->minIter_
         );
diff --git a/src/OpenFOAM/matrices/LduMatrix/Solvers/SmoothSolver/SmoothSolver.C b/src/OpenFOAM/matrices/LduMatrix/Solvers/SmoothSolver/SmoothSolver.C
index 9c5a5c80fd7e895f7ca0a7bf4420a30d8d577d45..ed10c2ae1211b1350579ff069cf6318e32b30ffc 100644
--- a/src/OpenFOAM/matrices/LduMatrix/Solvers/SmoothSolver/SmoothSolver.C
+++ b/src/OpenFOAM/matrices/LduMatrix/Solvers/SmoothSolver/SmoothSolver.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -110,7 +111,7 @@ Foam::SmoothSolver<Type, DType, LUType>::solve(Field<Type>& psi) const
             solverPerf.finalResidual() = solverPerf.initialResidual();
         }
 
-        if (LduMatrix<Type, DType, LUType>::debug >= 2)
+        if ((this->log_ >= 2) || (LduMatrix<Type, DType, LUType>::debug >= 2))
         {
             Info<< "   Normalisation factor = " << normFactor << endl;
         }
@@ -120,7 +121,12 @@ Foam::SmoothSolver<Type, DType, LUType>::solve(Field<Type>& psi) const
         if
         (
             this->minIter_ > 0
-         || !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+         || !solverPerf.checkConvergence
+            (
+                this->tolerance_,
+                this->relTol_,
+                this->log_
+            )
         )
         {
             autoPtr<typename LduMatrix<Type, DType, LUType>::smoother>
@@ -150,7 +156,12 @@ Foam::SmoothSolver<Type, DType, LUType>::solve(Field<Type>& psi) const
             (
                 (
                     (nIter += nSweeps_) < this->maxIter_
-                && !solverPerf.checkConvergence(this->tolerance_, this->relTol_)
+                && !solverPerf.checkConvergence
+                    (
+                        this->tolerance_,
+                        this->relTol_,
+                        this->log_
+                    )
                 )
              || nIter < this->minIter_
             );
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
index 644de9caa0df36b69d8bbfae87b259800a3ae7d8..f7642bad246a93e90e73d4e0ab4f46e1d8eab844 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2016-2019 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -120,6 +120,9 @@ public:
             //- Minimum number of iterations in the solver
             label minIter_;
 
+            //- Level of verbosity in the solver output statements
+            label log_;
+
             //- Final convergence tolerance
             scalar tolerance_;
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
index b750b74a65b1b736d5076f2c9079268abab6f76f..8a3f5ad6a1ed4953bf76f3aa2a554bfeabcd9e81 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -164,6 +164,7 @@ void Foam::lduMatrix::solver::readControls()
 {
     maxIter_ = controlDict_.getOrDefault<label>("maxIter", defaultMaxIter_);
     minIter_ = controlDict_.getOrDefault<label>("minIter", 0);
+    log_ = controlDict_.getOrDefault<label>("log", 1);
     tolerance_ = controlDict_.getOrDefault<scalar>("tolerance", 1e-6);
     relTol_ = controlDict_.getOrDefault<scalar>("relTol", 0);
 }
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
index 946ae1d49f62f71abe8a71ec656ac4a343f3ddd5..bd79b8c07fd3666e802d6fd6bd068ec21c9695f0 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
@@ -6,6 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -208,8 +209,7 @@ Foam::GAMGSolver::GAMGSolver
         }
     }
 
-
-    if (debug & 2)
+    if ((log_ >= 2) || (debug & 2))
     {
         for
         (
@@ -368,7 +368,7 @@ void Foam::GAMGSolver::readControls()
     controlDict_.readIfPresent("scaleCorrection", scaleCorrection_);
     controlDict_.readIfPresent("directSolveCoarsest", directSolveCoarsest_);
 
-    if (debug)
+    if ((log_ >= 2) || debug)
     {
         Info<< "GAMGSolver settings :"
             << " cacheAgglomeration:" << cacheAgglomeration_
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
index a8d1dbf175e33414451659f7c1d203690979b4e1..58076216e33e31c7bf2cf1cccb946e5087db59c6 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2016-2019 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -59,7 +59,7 @@ Foam::solverPerformance Foam::GAMGSolver::solve
     solveScalar normFactor =
         this->normFactor(psi, tsource(), Apsi, finestCorrection);
 
-    if (debug >= 2)
+    if ((log_ >= 2) || (debug >= 2))
     {
         Pout<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -87,7 +87,7 @@ Foam::solverPerformance Foam::GAMGSolver::solve
     if
     (
         minIter_ > 0
-     || !solverPerf.checkConvergence(tolerance_, relTol_)
+     || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
     )
     {
         // Create coarse grid correction fields
@@ -144,7 +144,7 @@ Foam::solverPerformance Foam::GAMGSolver::solve
                 matrix().mesh().comm()
             )/normFactor;
 
-            if (debug >= 2)
+            if ((log_ >= 2) || (debug >= 2))
             {
                 solverPerf.print(Info.masterStream(matrix().mesh().comm()));
             }
@@ -152,7 +152,7 @@ Foam::solverPerformance Foam::GAMGSolver::solve
         (
             (
               ++solverPerf.nIterations() < maxIter_
-            && !solverPerf.checkConvergence(tolerance_, relTol_)
+            && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
             )
          || solverPerf.nIterations() < minIter_
         );
@@ -193,7 +193,7 @@ void Foam::GAMGSolver::Vcycle
     // Restrict finest grid residual for the next level up.
     agglomeration_.restrictField(coarseSources[0], finestResidual, 0, true);
 
-    if (debug >= 2 && nPreSweeps_)
+    if (nPreSweeps_ && ((log_ >= 2) || (debug >= 2)))
     {
         Pout<< "Pre-smoothing scaling factors: ";
     }
@@ -274,7 +274,7 @@ void Foam::GAMGSolver::Vcycle
         }
     }
 
-    if (debug >= 2 && nPreSweeps_)
+    if (nPreSweeps_ && ((log_ >= 2) || (debug >= 2)))
     {
         Pout<< endl;
     }
@@ -290,7 +290,7 @@ void Foam::GAMGSolver::Vcycle
         );
     }
 
-    if (debug >= 2)
+    if ((log_ >= 2) || (debug >= 2))
     {
         Pout<< "Post-smoothing scaling factors: ";
     }
@@ -703,7 +703,7 @@ void Foam::GAMGSolver::solveCoarsestLevel
             )
         );
 
-        if (debug)
+        if ((log_ >= 2) || debug)
         {
             coarseSolverPerf.print(Info.masterStream(coarseComm));
         }
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
index cba20eba933b24c45f0bf88ceccc11efd181fa2f..d033d292fc336e8d890e90046c92ed205238a243 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -111,7 +111,7 @@ Foam::solverPerformance Foam::PBiCG::solve
     // --- Calculate normalisation factor
     const solveScalar normFactor = this->normFactor(psi, tsource(), wA, pA);
 
-    if (lduMatrix::debug >= 2)
+    if ((log_ >= 2) || (lduMatrix::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -126,7 +126,7 @@ Foam::solverPerformance Foam::PBiCG::solve
     if
     (
         minIter_ > 0
-     || !solverPerf.checkConvergence(tolerance_, relTol_)
+     || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
     )
     {
         solveScalarField pT(nCells, 0);
@@ -217,7 +217,7 @@ Foam::solverPerformance Foam::PBiCG::solve
         (
             (
               ++solverPerf.nIterations() < maxIter_
-            && !solverPerf.checkConvergence(tolerance_, relTol_)
+            && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
             )
          || solverPerf.nIterations() < minIter_
         );
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C
index 086fba5c07b52f00b5517cfcef6996a0cf368f16..ccf535206fedd9b582dceb0aa81f6f1476d732aa 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCGStab/PBiCGStab.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2016-2017 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -110,7 +110,7 @@ Foam::solverPerformance Foam::PBiCGStab::scalarSolve
     // --- Calculate normalisation factor
     const solveScalar normFactor = this->normFactor(psi, source, yA, pA);
 
-    if (lduMatrix::debug >= 2)
+    if ((log_ >= 2) || (lduMatrix::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -125,7 +125,7 @@ Foam::solverPerformance Foam::PBiCGStab::scalarSolve
     if
     (
         minIter_ > 0
-     || !solverPerf.checkConvergence(tolerance_, relTol_)
+     || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
     )
     {
         solveScalarField AyA(nCells);
@@ -219,7 +219,7 @@ Foam::solverPerformance Foam::PBiCGStab::scalarSolve
             if
             (
                 solverPerf.nIterations() >= minIter_
-             && solverPerf.checkConvergence(tolerance_, relTol_)
+             && solverPerf.checkConvergence(tolerance_, relTol_, log_)
             )
             {
                 for (label cell=0; cell<nCells; cell++)
@@ -258,7 +258,7 @@ Foam::solverPerformance Foam::PBiCGStab::scalarSolve
         (
             (
               ++solverPerf.nIterations() < maxIter_
-            && !solverPerf.checkConvergence(tolerance_, relTol_)
+            && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
             )
          || solverPerf.nIterations() < minIter_
         );
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
index 2b34a6aad378ef4286b8407114d622e72c2d8a56..ce33493e17cc807243cd896073c60c9e561a5743 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2019 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -110,7 +110,7 @@ Foam::solverPerformance Foam::PCG::scalarSolve
     // --- Calculate normalisation factor
     solveScalar normFactor = this->normFactor(psi, source, wA, pA);
 
-    if (lduMatrix::debug >= 2)
+    if ((log_ >= 2) || (lduMatrix::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -125,7 +125,7 @@ Foam::solverPerformance Foam::PCG::scalarSolve
     if
     (
         minIter_ > 0
-     || !solverPerf.checkConvergence(tolerance_, relTol_)
+     || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
     )
     {
         // --- Select and construct the preconditioner
@@ -193,7 +193,7 @@ Foam::solverPerformance Foam::PCG::scalarSolve
         (
             (
               ++solverPerf.nIterations() < maxIter_
-            && !solverPerf.checkConvergence(tolerance_, relTol_)
+            && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
             )
          || solverPerf.nIterations() < minIter_
         );
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PPCG/PPCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PPCG/PPCG.C
index d231a2daac6ff7d1b754b677aba1726e5a81d953..6ad660f76f1406b498f62f0afea668c1e7cd3bce 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PPCG/PPCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PPCG/PPCG.C
@@ -107,7 +107,7 @@ Foam::solverPerformance Foam::PPCG::scalarSolveCG
     solveScalarField p(nCells);
     const solveScalar normFactor = this->normFactor(psi, source, w, p);
 
-    if (lduMatrix::debug >= 2)
+    if ((log_ >= 2) || (lduMatrix::debug >= 2))
     {
         Info<< "   Normalisation factor = " << normFactor << endl;
     }
@@ -190,7 +190,7 @@ Foam::solverPerformance Foam::PPCG::scalarSolveCG
         if
         (
             (minIter_ <= 0 || solverPerf.nIterations() >= minIter_)
-         && solverPerf.checkConvergence(tolerance_, relTol_)
+         && solverPerf.checkConvergence(tolerance_, relTol_, log_)
         )
         {
             break;
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
index 78e6a7c62d8cf6235d31102b30ca73c9313c103e..1c25d08ccba87939071be3f4ab529a985f6c2dd2 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2014 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -149,7 +149,7 @@ Foam::solverPerformance Foam::smoothSolver::solve
             solverPerf.finalResidual() = solverPerf.initialResidual();
         }
 
-        if (lduMatrix::debug >= 2)
+        if ((log_ >= 2) || (lduMatrix::debug >= 2))
         {
             Info.masterStream(matrix().mesh().comm())
                 << "   Normalisation factor = " << normFactor << endl;
@@ -160,7 +160,7 @@ Foam::solverPerformance Foam::smoothSolver::solve
         if
         (
             minIter_ > 0
-         || !solverPerf.checkConvergence(tolerance_, relTol_)
+         || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
         )
         {
             addProfiling(solve, "lduMatrix::smoother." + fieldName_);
@@ -203,7 +203,7 @@ Foam::solverPerformance Foam::smoothSolver::solve
             (
                 (
                     (solverPerf.nIterations() += nSweeps_) < maxIter_
-                && !solverPerf.checkConvergence(tolerance_, relTol_)
+                && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
                 )
              || solverPerf.nIterations() < minIter_
             );
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
index 07c552b52deb58a857065cc2f6d4e0b771a679d8..111761404ed579b792dbd68ac92edb579d0bc852 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -222,7 +222,14 @@ Foam::SolverPerformance<Type> Foam::fvMatrix<Type>::solveSegregated
             solverControls
         )->solve(psiCmpt, sourceCmpt, cmpt);
 
-        if (SolverPerformance<Type>::debug)
+        const label log =
+            solverControls.getOrDefault<label>
+            (
+                "log",
+                SolverPerformance<Type>::debug
+            );
+
+        if (log)
         {
             solverPerf.print(Info.masterStream(this->mesh().comm()));
         }
@@ -289,7 +296,14 @@ Foam::SolverPerformance<Type> Foam::fvMatrix<Type>::solveCoupled
         coupledMatrixSolver->solve(psi)
     );
 
-    if (SolverPerformance<Type>::debug)
+    const label log =
+        solverControls.getOrDefault<label>
+        (
+            "log",
+            SolverPerformance<Type>::debug
+        );
+
+    if (log)
     {
         solverPerf.print(Info.masterStream(this->mesh().comm()));
     }
diff --git a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
index 99f4baf92c28b045c4ea24424eeb5b245b3aa39b..029b0d004e7b8ff56b57a069b89c36882eb360ad 100644
--- a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C
@@ -137,7 +137,10 @@ Foam::solverPerformance Foam::fvMatrix<Foam::scalar>::fvSolver::solve
         totalSource
     );
 
-    if (solverPerformance::debug)
+    const label log =
+        solverControls.getOrDefault<label>("log", solverPerformance::debug);
+
+    if (log)
     {
         solverPerf.print(Info.masterStream(fvMat_.mesh().comm()));
     }
@@ -264,7 +267,10 @@ Foam::solverPerformance Foam::fvMatrix<Foam::scalar>::solveSegregated
         }
     }
 
-    if (solverPerformance::debug)
+    const label log =
+        solverControls.getOrDefault<label>("log", solverPerformance::debug);
+
+    if (log)
     {
         solverPerf.print(Info.masterStream(mesh().comm()));
     }