Diffusion number: Corrected in chtMultiRegionFoam and pyrolysisModels::reactingOneDim

Resolves bug-report https://bugs.openfoam.org/view.php?id=2512

applications/solvers/heatTransfer/chtMultiRegionFoam/solid/solidRegionDiffNo.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -24,7 +24,9 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "solidRegionDiffNo.H"
-#include "fvc.H"
+#include "surfaceInterpolate.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 Foam::scalar Foam::solidRegionDiffNo
 (
@@ -34,21 +36,16 @@ Foam::scalar Foam::solidRegionDiffNo
     const volScalarField& kappa
 )
 {
-    scalar DiNum = 0.0;
-    scalar meanDiNum = 0.0;
-
-    //- Take care: can have fluid domains with 0 cells so do not test for
-    //  zero internal faces.
     surfaceScalarField kapparhoCpbyDelta
     (
-        mesh.surfaceInterpolation::deltaCoeffs()
-      * fvc::interpolate(kappa)
-      / fvc::interpolate(Cprho)
+        sqr(mesh.surfaceInterpolation::deltaCoeffs())
+       *fvc::interpolate(kappa)
+       /fvc::interpolate(Cprho)
     );
 
-    DiNum = max(kapparhoCpbyDelta).value()*runTime.deltaT().value();
-
-    meanDiNum = (average(kapparhoCpbyDelta)).value()*runTime.deltaT().value();
+    const scalar DiNum = max(kapparhoCpbyDelta).value()*runTime.deltaTValue();
+    const scalar meanDiNum =
+        average(kapparhoCpbyDelta).value()*runTime.deltaTValue();
 
     Info<< "Region: " << mesh.name() << " Diffusion Number mean: " << meanDiNum
         << " max: " << DiNum << endl;

src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,13 +25,11 @@ License
 
 #include "reactingOneDim.H"
 #include "addToRunTimeSelectionTable.H"
-#include "surfaceInterpolate.H"
 #include "fvm.H"
 #include "fvcDiv.H"
 #include "fvcVolumeIntegrate.H"
-#include "fvMatrices.H"
-#include "absorptionEmissionModel.H"
 #include "fvcLaplacian.H"
+#include "absorptionEmissionModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -63,7 +61,6 @@ void reactingOneDim::readReactingOneDimControls()
     coeffs().lookup("QrHSource") >> QrHSource_;
     useChemistrySolvers_ =
         coeffs().lookupOrDefault<bool>("useChemistrySolvers", true);
-
 }
 
 
@@ -251,9 +248,7 @@ void reactingOneDim::solveContinuity()
 
     fvScalarMatrix rhoEqn
     (
-        fvm::ddt(rho_)
-        ==
-      - solidChemistry_->RRg()
+        fvm::ddt(rho_) == -solidChemistry_->RRg()
     );
 
     if (regionMesh().moving())
@@ -287,9 +282,7 @@ void reactingOneDim::solveSpeciesMass()
 
         fvScalarMatrix YiEqn
         (
-            fvm::ddt(rho_, Yi)
-         ==
-            solidChemistry_->RRs(i)
+            fvm::ddt(rho_, Yi) == solidChemistry_->RRs(i)
         );
 
         if (regionMesh().moving())
@@ -624,9 +617,9 @@ scalar reactingOneDim::solidRegionDiffNo() const
     {
         surfaceScalarField KrhoCpbyDelta
         (
-            regionMesh().surfaceInterpolation::deltaCoeffs()
-          * fvc::interpolate(kappa())
-          / fvc::interpolate(Cp()*rho_)
+            sqr(regionMesh().surfaceInterpolation::deltaCoeffs())
+           *fvc::interpolate(kappa())
+           /fvc::interpolate(Cp()*rho_)
         );
 
         DiNum = max(KrhoCpbyDelta.primitiveField())*time().deltaTValue();