diff --git a/applications/solvers/combustion/chemFoam/output.H b/applications/solvers/combustion/chemFoam/output.H
index 5620fe6820..2e84c59a48 100644
--- a/applications/solvers/combustion/chemFoam/output.H
+++ b/applications/solvers/combustion/chemFoam/output.H
@@ -1,6 +1,6 @@
     runTime.write();
 
-    Info<< "Sh = " << Sh
+    Info<< "Qdot = " << Qdot
         << ", T = " << thermo.T()[0]
         << ", p = " << thermo.p()[0]
         << ", " << Y[0].name() << " = " << Y[0][0]
@@ -8,4 +8,3 @@
 
     post<< runTime.value() << token::TAB << thermo.T()[0] << token::TAB
         << thermo.p()[0] << endl;
-
diff --git a/applications/solvers/combustion/chemFoam/solveChemistry.H b/applications/solvers/combustion/chemFoam/solveChemistry.H
index cb4258c31d..e32a1ab0bd 100644
--- a/applications/solvers/combustion/chemFoam/solveChemistry.H
+++ b/applications/solvers/combustion/chemFoam/solveChemistry.H
@@ -1,3 +1,3 @@
-    dtChem = chemistry.solve(runTime.deltaT().value());
-    scalar Sh = chemistry.Sh()()[0]/rho[0];
-    integratedHeat += Sh*runTime.deltaT().value();
+dtChem = chemistry.solve(runTime.deltaT().value());
+scalar Qdot = chemistry.Qdot()()[0]/rho[0];
+integratedHeat += Qdot*runTime.deltaT().value();
diff --git a/applications/solvers/combustion/fireFoam/YEEqn.H b/applications/solvers/combustion/fireFoam/YEEqn.H
index 1c9c2176c5..f72d574d95 100644
--- a/applications/solvers/combustion/fireFoam/YEEqn.H
+++ b/applications/solvers/combustion/fireFoam/YEEqn.H
@@ -11,7 +11,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 {
     radiation->correct();
     combustion->correct();
-    dQ = combustion->dQ();
+    Qdot = combustion->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
@@ -66,7 +66,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
         )
       - fvm::laplacian(turbulence->alphaEff(), he)
      ==
-        combustion->Sh()
+        Qdot
       + radiation->Sh(thermo)
       + parcels.Sh(he)
       + surfaceFilm.Sh()
diff --git a/applications/solvers/combustion/fireFoam/createFields.H b/applications/solvers/combustion/fireFoam/createFields.H
index 30e7d30377..85d26aed87 100644
--- a/applications/solvers/combustion/fireFoam/createFields.H
+++ b/applications/solvers/combustion/fireFoam/createFields.H
@@ -131,18 +131,18 @@ Switch solvePyrolysisRegion
     additionalControlsDict.lookupOrDefault<bool>("solvePyrolysisRegion", true)
 );
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 
diff --git a/applications/solvers/combustion/reactingFoam/EEqn.H b/applications/solvers/combustion/reactingFoam/EEqn.H
index 9267c9a9be..f92598bbed 100644
--- a/applications/solvers/combustion/reactingFoam/EEqn.H
+++ b/applications/solvers/combustion/reactingFoam/EEqn.H
@@ -17,7 +17,7 @@
         )
       - fvm::laplacian(turbulence->alphaEff(), he)
      ==
-        reaction->Sh()
+        Qdot
       + fvOptions(rho, he)
     );
 
diff --git a/applications/solvers/combustion/reactingFoam/YEqn.H b/applications/solvers/combustion/reactingFoam/YEqn.H
index 7c3f0f117b..4d209a3d81 100644
--- a/applications/solvers/combustion/reactingFoam/YEqn.H
+++ b/applications/solvers/combustion/reactingFoam/YEqn.H
@@ -11,7 +11,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 
 {
     reaction->correct();
-    dQ = reaction->dQ();
+    Qdot = reaction->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
diff --git a/applications/solvers/combustion/reactingFoam/createFields.H b/applications/solvers/combustion/reactingFoam/createFields.H
index e4475e5d1b..6fcf83624c 100644
--- a/applications/solvers/combustion/reactingFoam/createFields.H
+++ b/applications/solvers/combustion/reactingFoam/createFields.H
@@ -117,18 +117,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/combustion/reactingFoam/rhoReactingBuoyantFoam/createFields.H b/applications/solvers/combustion/reactingFoam/rhoReactingBuoyantFoam/createFields.H
index a4b328f367..675aef3f5a 100644
--- a/applications/solvers/combustion/reactingFoam/rhoReactingBuoyantFoam/createFields.H
+++ b/applications/solvers/combustion/reactingFoam/rhoReactingBuoyantFoam/createFields.H
@@ -117,18 +117,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/combustion/reactingFoam/rhoReactingFoam/createFields.H b/applications/solvers/combustion/reactingFoam/rhoReactingFoam/createFields.H
index fbd45dd20c..687dffcc13 100644
--- a/applications/solvers/combustion/reactingFoam/rhoReactingFoam/createFields.H
+++ b/applications/solvers/combustion/reactingFoam/rhoReactingFoam/createFields.H
@@ -96,18 +96,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/combustion/reactingFoam/setRDeltaT.H b/applications/solvers/combustion/reactingFoam/setRDeltaT.H
index 5c24294ef6..28300dc134 100644
--- a/applications/solvers/combustion/reactingFoam/setRDeltaT.H
+++ b/applications/solvers/combustion/reactingFoam/setRDeltaT.H
@@ -77,7 +77,7 @@ License
     {
         volScalarField::Internal rDeltaTT
         (
-            mag(reaction->Sh())/(alphaTemp*rho*thermo.Cp()*T)
+            mag(Qdot)/(alphaTemp*rho*thermo.Cp()*T)
         );
 
         Info<< "    Temperature = "
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/EEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/EEqn.H
index 8d32932607..1a1f014da0 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/EEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/EEqn.H
@@ -18,7 +18,7 @@
       - fvm::laplacian(turbulence->alphaEff(), he)
      ==
         rho*(U&g)
-      + combustion->Sh()
+      + Qdot
       + coalParcels.Sh(he)
       + limestoneParcels.Sh(he)
       + radiation->Sh(thermo)
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H
index 4edbce9b42..49e3f70d89 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H
@@ -12,7 +12,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 
 {
     combustion->correct();
-    dQ = combustion->dQ();
+    Qdot = combustion->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/createFields.H b/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
index d3daf7265c..4e4f45e721 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
@@ -131,18 +131,18 @@ volScalarField dpdt
 Info<< "Creating field kinetic energy K\n" << endl;
 volScalarField K("K", 0.5*magSqr(U));
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/setRDeltaT.H b/applications/solvers/lagrangian/coalChemistryFoam/setRDeltaT.H
index 7e2ed6074d..c31ad040a1 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/setRDeltaT.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/setRDeltaT.H
@@ -80,7 +80,7 @@ License
             (
                (coalParcels.hsTrans() + limestoneParcels.hsTrans())
               /(mesh.V()*runTime.deltaT())
-             + combustion->Sh()()
+             + Qdot
             )
            /(
                alphaTemp
diff --git a/applications/solvers/lagrangian/reactingParcelFilmFoam/EEqn.H b/applications/solvers/lagrangian/reactingParcelFilmFoam/EEqn.H
index 4c2d95da9d..325009f067 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmFoam/EEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFilmFoam/EEqn.H
@@ -21,7 +21,7 @@
       + parcels.Sh(he)
       + surfaceFilm.Sh()
       + radiation->Sh(thermo)
-      + combustion->Sh()
+      + Qdot
       + fvOptions(rho, he)
     );
 
diff --git a/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H
index 555b5bc0bc..eb24899013 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H
@@ -12,7 +12,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 
 {
     combustion->correct();
-    dQ = combustion->dQ();
+    Qdot = combustion->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
diff --git a/applications/solvers/lagrangian/reactingParcelFilmFoam/createFields.H b/applications/solvers/lagrangian/reactingParcelFilmFoam/createFields.H
index 8ec3c94169..0159cdba15 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmFoam/createFields.H
+++ b/applications/solvers/lagrangian/reactingParcelFilmFoam/createFields.H
@@ -134,18 +134,18 @@ Switch solvePrimaryRegion
     additionalControlsDict.lookup("solvePrimaryRegion")
 );
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/EEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/EEqn.H
index 194e03c59d..2ba1ef693f 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/EEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/EEqn.H
@@ -20,7 +20,7 @@
         rho*(U&g)
       + parcels.Sh(he)
       + radiation->Sh(thermo)
-      + combustion->Sh()
+      + Qdot
       + fvOptions(rho, he)
     );
 
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H
index db6628936c..60a27fec85 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H
@@ -11,7 +11,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 
 {
     combustion->correct();
-    dQ = combustion->dQ();
+    Qdot = combustion->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/createFields.H b/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
index e5ed566eba..2da47d4ced 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
@@ -121,18 +121,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/setRDeltaT.H b/applications/solvers/lagrangian/reactingParcelFoam/setRDeltaT.H
index 63398fb001..28bc211e05 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/setRDeltaT.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/setRDeltaT.H
@@ -79,7 +79,7 @@ License
             mag
             (
                parcels.hsTrans()/(mesh.V()*runTime.deltaT())
-             + combustion->Sh()()
+             + Qdot
             )
            /(
                alphaTemp
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/EEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/EEqn.H
index 8239d44621..9bc1254428 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/EEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/EEqn.H
@@ -14,7 +14,7 @@
         rho*(U&g)
       + parcels.Sh(he)
       + radiation->Sh(thermo)
-      + combustion->Sh()
+      + Qdot
       + fvOptions(rho, he)
     );
 
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H
index 7c7cf9a4c0..1510f8e391 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H
@@ -11,7 +11,7 @@ tmp<fv::convectionScheme<scalar>> mvConvection
 
 {
     combustion->correct();
-    dQ = combustion->dQ();
+    Qdot = combustion->Qdot();
     volScalarField Yt(0.0*Y[0]);
 
     forAll(Y, i)
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createFields.H b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createFields.H
index 16f297c995..328fc67224 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createFields.H
@@ -103,18 +103,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
-        IOobject::NO_READ,
+        IOobject::READ_IF_PRESENT,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/lagrangian/sprayFoam/createFields.H b/applications/solvers/lagrangian/sprayFoam/createFields.H
index 91b7bc2b6d..e2124b3685 100644
--- a/applications/solvers/lagrangian/sprayFoam/createFields.H
+++ b/applications/solvers/lagrangian/sprayFoam/createFields.H
@@ -117,18 +117,18 @@ forAll(Y, i)
 }
 fields.add(thermo.he());
 
-volScalarField dQ
+volScalarField Qdot
 (
     IOobject
     (
-        "dQ",
+        "Qdot",
         runTime.timeName(),
         mesh,
         IOobject::NO_READ,
         IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+    dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
 );
 
 #include "createMRF.H"
diff --git a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C
index 1fe2b8da1e..ea8f5c10db 100644
--- a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C
@@ -135,7 +135,7 @@ int main(int argc, char *argv[])
 
         if (runTime.write())
         {
-            combustion->dQ()().write();
+            combustion->Qdot()().write();
         }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
diff --git a/applications/solvers/lagrangian/sprayFoam/sprayEngineFoam/sprayEngineFoam.C b/applications/solvers/lagrangian/sprayFoam/sprayEngineFoam/sprayEngineFoam.C
index 82cf4f6cab..a000aacb9b 100644
--- a/applications/solvers/lagrangian/sprayFoam/sprayEngineFoam/sprayEngineFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/sprayEngineFoam/sprayEngineFoam.C
@@ -110,7 +110,7 @@ int main(int argc, char *argv[])
 
         if (runTime.write())
         {
-            combustion->dQ()().write();
+            combustion->Qdot()().write();
         }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
diff --git a/applications/solvers/lagrangian/sprayFoam/sprayFoam.C b/applications/solvers/lagrangian/sprayFoam/sprayFoam.C
index 30dbe63830..187477ef32 100644
--- a/applications/solvers/lagrangian/sprayFoam/sprayFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/sprayFoam.C
@@ -102,7 +102,7 @@ int main(int argc, char *argv[])
 
             if (runTime.write())
             {
-                combustion->dQ()().write();
+                combustion->Qdot()().write();
             }
         }
         else
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
index b5dd4ab37a..cabb9e2ebf 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/AnisothermalPhaseModel/AnisothermalPhaseModel.C
@@ -135,7 +135,7 @@ Foam::AnisothermalPhaseModel<BasePhaseModel>::heEqn()
             he
         )
      ==
-        this->Sh()
+        this->Qdot()
     );
 
     // Add the appropriate pressure-work term
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.C
index bd4dda226f..e92cae4c9f 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.C
@@ -65,7 +65,7 @@ Foam::InertPhaseModel<BasePhaseModel>::R
 
 template<class BasePhaseModel>
 Foam::tmp<Foam::volScalarField>
-Foam::InertPhaseModel<BasePhaseModel>::Sh() const
+Foam::InertPhaseModel<BasePhaseModel>::Qdot() const
 {
     return tmp<volScalarField>
     (
@@ -73,7 +73,7 @@ Foam::InertPhaseModel<BasePhaseModel>::Sh() const
         (
             IOobject
             (
-                IOobject::groupName("Sh", this->name()),
+                IOobject::groupName("Qdot", this->name()),
                 this->mesh().time().timeName(),
                 this->mesh()
             ),
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.H
index 26c9661ab7..c152da488a 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.H
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/InertPhaseModel/InertPhaseModel.H
@@ -78,8 +78,8 @@ public:
                 volScalarField& Yi
             ) const;
 
-            //- Return the reaction heat source
-            virtual tmp<volScalarField> Sh() const;
+            //- Return the heat release rate
+            virtual tmp<volScalarField> Qdot() const;
 };
 
 
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.C
index b15f2898b5..aef23273f2 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.C
@@ -105,9 +105,9 @@ Foam::ReactingPhaseModel<BasePhaseModel, ReactionType>::R
 
 template<class BasePhaseModel, class ReactionType>
 Foam::tmp<Foam::volScalarField>
-Foam::ReactingPhaseModel<BasePhaseModel, ReactionType>::Sh() const
+Foam::ReactingPhaseModel<BasePhaseModel, ReactionType>::Qdot() const
 {
-    return reaction_->Sh();
+    return reaction_->Qdot();
 }
 
 
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.H
index e05c71b44a..fdda42ca09 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.H
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/ReactingPhaseModel/ReactingPhaseModel.H
@@ -87,8 +87,8 @@ public:
             volScalarField& Yi
         ) const;
 
-        //- Return the reacton heat source
-        virtual tmp<volScalarField> Sh() const;
+        //- Return heat release rate
+        virtual tmp<volScalarField> Qdot() const;
 };
 
 
diff --git a/src/combustionModels/PaSR/PaSR.C b/src/combustionModels/PaSR/PaSR.C
index 5fd8fdf104..fa3b7feaed 100644
--- a/src/combustionModels/PaSR/PaSR.C
+++ b/src/combustionModels/PaSR/PaSR.C
@@ -116,29 +116,14 @@ Foam::combustionModels::PaSR<Type>::R(volScalarField& Y) const
 
 template<class Type>
 Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::PaSR<Type>::dQ() const
+Foam::combustionModels::PaSR<Type>::Qdot() const
 {
     return tmp<volScalarField>
     (
         new volScalarField
         (
             IOobject::groupName("PaSR:dQ", this->phaseName_),
-            kappa_*laminar<Type>::dQ()
-        )
-    );
-}
-
-
-template<class Type>
-Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::PaSR<Type>::Sh() const
-{
-    return tmp<volScalarField>
-    (
-        new volScalarField
-        (
-            IOobject::groupName("PaSR:Sh", this->phaseName_),
-            kappa_*laminar<Type>::Sh()
+            kappa_*laminar<Type>::Qdot()
         )
     );
 }
diff --git a/src/combustionModels/PaSR/PaSR.H b/src/combustionModels/PaSR/PaSR.H
index c855ea0dc0..88e347134a 100644
--- a/src/combustionModels/PaSR/PaSR.H
+++ b/src/combustionModels/PaSR/PaSR.H
@@ -107,11 +107,8 @@ public:
         //- Fuel consumption rate matrix.
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const;
-
-        //-  Return source for enthalpy equation [kg/m/s3]
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/combustionModels/combustionModel/combustionModel.C b/src/combustionModels/combustionModel/combustionModel.C
index b27a5f2320..3094e4bdff 100644
--- a/src/combustionModels/combustionModel/combustionModel.C
+++ b/src/combustionModels/combustionModel/combustionModel.C
@@ -81,7 +81,6 @@ Foam::combustionModel::~combustionModel()
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-
 bool Foam::combustionModel::read()
 {
     if (regIOobject::read())
@@ -97,25 +96,4 @@ bool Foam::combustionModel::read()
 }
 
 
-Foam::tmp<Foam::volScalarField> Foam::combustionModel::Sh() const
-{
-    return tmp<Foam::volScalarField>
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                IOobject::groupName("Sh", phaseName_),
-                mesh_.time().timeName(),
-                mesh_,
-                IOobject::NO_READ,
-                IOobject::NO_WRITE
-            ),
-            mesh_,
-            dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
-        )
-    );
-}
-
-
 // ************************************************************************* //
diff --git a/src/combustionModels/combustionModel/combustionModel.H b/src/combustionModels/combustionModel/combustionModel.H
index d0927b7d01..9845dad4f5 100644
--- a/src/combustionModels/combustionModel/combustionModel.H
+++ b/src/combustionModels/combustionModel/combustionModel.H
@@ -137,11 +137,8 @@ public:
         //- Fuel consumption rate matrix, i.e. source term for fuel equation
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const = 0;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const = 0;
-
-        //-  Return source for enthalpy equation [kg/m/s3]
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const = 0;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/combustionModels/laminar/laminar.C b/src/combustionModels/laminar/laminar.C
index dcf9c6450c..679a086780 100644
--- a/src/combustionModels/laminar/laminar.C
+++ b/src/combustionModels/laminar/laminar.C
@@ -136,15 +136,15 @@ Foam::combustionModels::laminar<Type>::R(volScalarField& Y) const
 
 template<class Type>
 Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::laminar<Type>::dQ() const
+Foam::combustionModels::laminar<Type>::Qdot() const
 {
-    tmp<volScalarField> tdQ
+    tmp<volScalarField> tQdot
     (
         new volScalarField
         (
             IOobject
             (
-                IOobject::groupName(typeName + ":dQ", this->phaseName_),
+                IOobject::groupName(typeName + ":Qdot", this->phaseName_),
                 this->mesh().time().timeName(),
                 this->mesh(),
                 IOobject::NO_READ,
@@ -152,47 +152,16 @@ Foam::combustionModels::laminar<Type>::dQ() const
                 false
             ),
             this->mesh(),
-            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+            dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
         )
     );
 
     if (this->active())
     {
-        tdQ.ref() = this->chemistryPtr_->dQ();
+        tQdot.ref() = this->chemistryPtr_->Qdot();
     }
 
-    return tdQ;
-}
-
-
-template<class Type>
-Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::laminar<Type>::Sh() const
-{
-    tmp<volScalarField> tSh
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                IOobject::groupName(typeName + ":Sh", this->phaseName_),
-                this->mesh().time().timeName(),
-                this->mesh(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            this->mesh(),
-            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
-        )
-    );
-
-    if (this->active())
-    {
-        tSh.ref() = this->chemistryPtr_->Sh();
-    }
-
-    return tSh;
+    return tQdot;
 }
 
 
diff --git a/src/combustionModels/laminar/laminar.H b/src/combustionModels/laminar/laminar.H
index 4b1fad879f..41b00799cf 100644
--- a/src/combustionModels/laminar/laminar.H
+++ b/src/combustionModels/laminar/laminar.H
@@ -105,11 +105,8 @@ public:
         //- Fuel consumption rate matrix.
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const;
-
-        //-  Return source for enthalpy equation [kg/m/s3]
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/combustionModels/noCombustion/noCombustion.C b/src/combustionModels/noCombustion/noCombustion.C
index 7caeac0fe4..190459f5ae 100644
--- a/src/combustionModels/noCombustion/noCombustion.C
+++ b/src/combustionModels/noCombustion/noCombustion.C
@@ -75,15 +75,15 @@ Foam::combustionModels::noCombustion<CombThermoType>::R
 
 template<class CombThermoType>
 Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::noCombustion<CombThermoType>::dQ() const
+Foam::combustionModels::noCombustion<CombThermoType>::Qdot() const
 {
-    tmp<volScalarField> tdQ
+    tmp<volScalarField> tQdot
     (
         new volScalarField
         (
             IOobject
             (
-                IOobject::groupName("dQ", this->phaseName_),
+                IOobject::groupName("Qdot", this->phaseName_),
                 this->mesh().time().timeName(),
                 this->mesh(),
                 IOobject::NO_READ,
@@ -91,37 +91,11 @@ Foam::combustionModels::noCombustion<CombThermoType>::dQ() const
                 false
             ),
             this->mesh(),
-            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+            dimensionedScalar("Qdot", dimEnergy/dimVolume/dimTime, 0.0)
         )
     );
 
-    return tdQ;
-}
-
-
-template<class CombThermoType>
-Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::noCombustion<CombThermoType>::Sh() const
-{
-    tmp<volScalarField> tSh
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                IOobject::groupName("Sh", this->phaseName_),
-                this->mesh().time().timeName(),
-                this->mesh(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            this->mesh(),
-            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
-        )
-    );
-
-    return tSh;
+    return tQdot;
 }
 
 
diff --git a/src/combustionModels/noCombustion/noCombustion.H b/src/combustionModels/noCombustion/noCombustion.H
index 1ea6d410aa..27783065cb 100644
--- a/src/combustionModels/noCombustion/noCombustion.H
+++ b/src/combustionModels/noCombustion/noCombustion.H
@@ -89,11 +89,8 @@ public:
         //- Fuel consumption rate matrix
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const;
-
-        //-  Return source for enthalpy equation [kg/m/s3]
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/combustionModels/singleStepCombustion/singleStepCombustion.C b/src/combustionModels/singleStepCombustion/singleStepCombustion.C
index c2a243747b..668a3c0726 100644
--- a/src/combustionModels/singleStepCombustion/singleStepCombustion.C
+++ b/src/combustionModels/singleStepCombustion/singleStepCombustion.C
@@ -127,7 +127,7 @@ tmp<fvScalarMatrix> singleStepCombustion<CombThermoType, ThermoType>::R
 
 template<class CombThermoType, class ThermoType>
 tmp<volScalarField>
-singleStepCombustion<CombThermoType, ThermoType>::Sh() const
+singleStepCombustion<CombThermoType, ThermoType>::Qdot() const
 {
     const label fuelI = singleMixturePtr_->fuelIndex();
     volScalarField& YFuel =
@@ -137,37 +137,6 @@ singleStepCombustion<CombThermoType, ThermoType>::Sh() const
 }
 
 
-template<class CombThermoType, class ThermoType>
-tmp<volScalarField>
-singleStepCombustion<CombThermoType, ThermoType>::dQ() const
-{
-    tmp<volScalarField> tdQ
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                IOobject::groupName("dQ", this->phaseName_),
-                this->mesh_.time().timeName(),
-                this->mesh_,
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            this->mesh_,
-            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
-        )
-    );
-
-    if (this->active())
-    {
-        volScalarField& dQ = tdQ.ref();
-        dQ.ref() = this->mesh().V()*Sh()();
-    }
-    return tdQ;
-}
-
-
 template<class CombThermoType, class ThermoType>
 bool singleStepCombustion<CombThermoType, ThermoType>::read()
 {
diff --git a/src/combustionModels/singleStepCombustion/singleStepCombustion.H b/src/combustionModels/singleStepCombustion/singleStepCombustion.H
index 181965bd52..3c5d8fb7e0 100644
--- a/src/combustionModels/singleStepCombustion/singleStepCombustion.H
+++ b/src/combustionModels/singleStepCombustion/singleStepCombustion.H
@@ -99,11 +99,8 @@ public:
         //- Fuel consumption rate matrix
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const;
-
-        //- Sensible enthalpy source term
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/combustionModels/zoneCombustion/zoneCombustion.C b/src/combustionModels/zoneCombustion/zoneCombustion.C
index 3cec4e1197..3284c2aca5 100644
--- a/src/combustionModels/zoneCombustion/zoneCombustion.C
+++ b/src/combustionModels/zoneCombustion/zoneCombustion.C
@@ -165,17 +165,9 @@ Foam::combustionModels::zoneCombustion<Type>::R(volScalarField& Y) const
 
 template<class Type>
 Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::zoneCombustion<Type>::dQ() const
+Foam::combustionModels::zoneCombustion<Type>::Qdot() const
 {
-    return filter(combustionModelPtr_->dQ());
-}
-
-
-template<class Type>
-Foam::tmp<Foam::volScalarField>
-Foam::combustionModels::zoneCombustion<Type>::Sh() const
-{
-    return filter(combustionModelPtr_->Sh());
+    return filter(combustionModelPtr_->Qdot());
 }
 
 
diff --git a/src/combustionModels/zoneCombustion/zoneCombustion.H b/src/combustionModels/zoneCombustion/zoneCombustion.H
index ca014f98dd..351b40cab0 100644
--- a/src/combustionModels/zoneCombustion/zoneCombustion.H
+++ b/src/combustionModels/zoneCombustion/zoneCombustion.H
@@ -114,11 +114,8 @@ public:
         //- Fuel consumption rate matrix.
         virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        //- Heat release rate calculated from fuel consumption rate matrix
-        virtual tmp<volScalarField> dQ() const;
-
-        //-  Return source for enthalpy equation [kg/m/s3]
-        virtual tmp<volScalarField> Sh() const;
+        //- Heat release rate [kg/m/s3]
+        virtual tmp<volScalarField> Qdot() const;
 
         //- Update properties from given dictionary
         virtual bool read();
diff --git a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C
index bfd963079d..be8584b5da 100644
--- a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C
+++ b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.C
@@ -329,7 +329,7 @@ void reactingOneDim::solveEnergy()
       + fvc::laplacian(alpha, h_)
       - fvc::laplacian(kappa(), T())
      ==
-        chemistrySh_
+        chemistryQdot_
       - fvm::Sp(solidChemistry_->RRg(), h_)
     );
 
@@ -371,7 +371,7 @@ void reactingOneDim::calculateMassTransfer()
 
     if (infoOutput_)
     {
-        totalHeatRR_ = fvc::domainIntegrate(chemistrySh_);
+        totalHeatRR_ = fvc::domainIntegrate(chemistryQdot_);
 
         addedGasMass_ +=
             fvc::domainIntegrate(solidChemistry_->RRg())*time_.deltaT();
@@ -440,11 +440,11 @@ reactingOneDim::reactingOneDim
         dimensionedScalar("zero", dimEnergy/dimTime, 0.0)
     ),
 
-    chemistrySh_
+    chemistryQdot_
     (
         IOobject
         (
-            "chemistrySh",
+            "chemistryQdot",
             time().timeName(),
             regionMesh(),
             IOobject::NO_READ,
@@ -540,11 +540,11 @@ reactingOneDim::reactingOneDim
         dimensionedScalar("zero", dimEnergy/dimTime, 0.0)
     ),
 
-    chemistrySh_
+    chemistryQdot_
     (
         IOobject
         (
-            "chemistrySh",
+            "chemistryQdot",
             time().timeName(),
             regionMesh(),
             IOobject::NO_READ,
@@ -705,7 +705,7 @@ void reactingOneDim::evolveRegion()
 
     solveContinuity();
 
-    chemistrySh_ = solidChemistry_->Sh()();
+    chemistryQdot_ = solidChemistry_->Qdot()();
 
     updateFields();
 
diff --git a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H
index 0826d6fc23..a5b6ecfbce 100644
--- a/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H
+++ b/src/regionModels/pyrolysisModels/reactingOneDim/reactingOneDim.H
@@ -117,8 +117,8 @@ protected:
             //- Sensible enthalpy gas flux [J/m2/s]
             volScalarField phiHsGas_;
 
-            //- Heat release [J/s/m3]
-            volScalarField chemistrySh_;
+            //- Heat release rate [J/s/m3]
+            volScalarField chemistryQdot_;
 
 
         // Source term fields
diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H b/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H
index 17b64f542b..c85cf37644 100644
--- a/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H
+++ b/src/thermophysicalModels/chemistryModel/chemistryModel/basicChemistryModel/basicChemistryModel.H
@@ -174,11 +174,8 @@ public:
                 //- Return the chemical time scale
                 virtual tmp<volScalarField> tc() const = 0;
 
-                //- Return source for enthalpy equation [kg/m/s3]
-                virtual tmp<volScalarField> Sh() const = 0;
-
-                //- Return the heat release, i.e. enthalpy/sec [m2/s3]
-                virtual tmp<volScalarField> dQ() const = 0;
+                //- Return the heat release rate [kg/m/s3]
+                virtual tmp<volScalarField> Qdot() const = 0;
 };
 
 
diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C
index a779d0c5b2..661da67c3b 100644
--- a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C
+++ b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.C
@@ -545,15 +545,15 @@ Foam::chemistryModel<CompType, ThermoType>::tc() const
 
 template<class CompType, class ThermoType>
 Foam::tmp<Foam::volScalarField>
-Foam::chemistryModel<CompType, ThermoType>::Sh() const
+Foam::chemistryModel<CompType, ThermoType>::Qdot() const
 {
-    tmp<volScalarField> tSh
+    tmp<volScalarField> tQdot
     (
         new volScalarField
         (
             IOobject
             (
-                "Sh",
+                "Qdot",
                 this->mesh_.time().timeName(),
                 this->mesh_,
                 IOobject::NO_READ,
@@ -561,57 +561,25 @@ Foam::chemistryModel<CompType, ThermoType>::Sh() const
                 false
             ),
             this->mesh_,
-            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
+            dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
         )
     );
 
     if (this->chemistry_)
     {
-        scalarField& Sh = tSh.ref();
+        scalarField& Qdot = tQdot.ref();
 
         forAll(Y_, i)
         {
-            forAll(Sh, celli)
+            forAll(Qdot, celli)
             {
                 const scalar hi = specieThermo_[i].Hc();
-                Sh[celli] -= hi*RR_[i][celli];
+                Qdot[celli] -= hi*RR_[i][celli];
             }
         }
     }
 
-    return tSh;
-}
-
-
-template<class CompType, class ThermoType>
-Foam::tmp<Foam::volScalarField>
-Foam::chemistryModel<CompType, ThermoType>::dQ() const
-{
-    tmp<volScalarField> tdQ
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                "dQ",
-                this->mesh_.time().timeName(),
-                this->mesh_,
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            this->mesh_,
-            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
-        )
-    );
-
-    if (this->chemistry_)
-    {
-        volScalarField& dQ = tdQ.ref();
-        dQ.ref() = this->mesh_.V()*Sh()();
-    }
-
-    return tdQ;
+    return tQdot;
 }
 
 
diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H
index 98b2ec2a91..552aed094f 100644
--- a/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H
+++ b/src/thermophysicalModels/chemistryModel/chemistryModel/chemistryModel/chemistryModel.H
@@ -231,11 +231,8 @@ public:
             //- Return the chemical time scale
             virtual tmp<volScalarField> tc() const;
 
-            //- Return source for enthalpy equation [kg/m/s3]
-            virtual tmp<volScalarField> Sh() const;
-
-            //- Return the heat release, i.e. enthalpy/sec [kg/m2/s3]
-            virtual tmp<volScalarField> dQ() const;
+            //- Return the heat release rate [kg/m/s3]
+            virtual tmp<volScalarField> Qdot() const;
 
 
         // ODE functions (overriding abstract functions in ODE.H)
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
index 69c68ad184..bcda522186 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
@@ -293,32 +293,32 @@ Foam::radiation::greyMeanAbsorptionEmission::ECont(const label bandI) const
         )
     );
 
-    if (mesh_.foundObject<volScalarField>("dQ"))
+    if (mesh_.foundObject<volScalarField>("Qdot"))
     {
-        const volScalarField& dQ =
-            mesh_.lookupObject<volScalarField>("dQ");
+        const volScalarField& Qdot =
+            mesh_.lookupObject<volScalarField>("Qdot");
 
-        if (dQ.dimensions() == dimEnergy/dimTime)
+        if (Qdot.dimensions() == dimEnergy/dimTime)
         {
-            E.ref().primitiveFieldRef() = EhrrCoeff_*dQ/mesh_.V();
+            E.ref().primitiveFieldRef() = EhrrCoeff_*Qdot/mesh_.V();
         }
-        else if (dQ.dimensions() == dimEnergy/dimTime/dimVolume)
+        else if (Qdot.dimensions() == dimEnergy/dimTime/dimVolume)
         {
-            E.ref().primitiveFieldRef() = EhrrCoeff_*dQ;
+            E.ref().primitiveFieldRef() = EhrrCoeff_*Qdot;
         }
         else
         {
             if (debug)
             {
                 WarningInFunction
-                    << "Incompatible dimensions for dQ field" << endl;
+                    << "Incompatible dimensions for Qdot field" << endl;
             }
         }
     }
     else
     {
         WarningInFunction
-          << "dQ field not found in mesh" << endl;
+          << "Qdot field not found in mesh" << endl;
     }
 
     return E;
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
index d7899338d0..6e31e837ab 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
@@ -249,31 +249,32 @@ Foam::radiation::wideBandAbsorptionEmission::ECont(const label bandI) const
         )
     );
 
-    if (mesh().foundObject<volScalarField>("dQ"))
+    if (mesh().foundObject<volScalarField>("Qdot"))
     {
-        const volScalarField& dQ = mesh().lookupObject<volScalarField>("dQ");
+        const volScalarField& Qdot =
+            mesh().lookupObject<volScalarField>("Qdot");
 
-        if (dQ.dimensions() == dimEnergy/dimTime)
+        if (Qdot.dimensions() == dimEnergy/dimTime)
         {
             E.ref().primitiveFieldRef() =
                 iEhrrCoeffs_[bandI]
-               *dQ.primitiveField()
+               *Qdot.primitiveField()
                *(iBands_[bandI][1] - iBands_[bandI][0])
                /totalWaveLength_
                /mesh_.V();
         }
-        else if (dQ.dimensions() == dimEnergy/dimTime/dimVolume)
+        else if (Qdot.dimensions() == dimEnergy/dimTime/dimVolume)
         {
             E.ref().primitiveFieldRef() =
                 iEhrrCoeffs_[bandI]
-               *dQ.primitiveField()
+               *Qdot.primitiveField()
                *(iBands_[bandI][1] - iBands_[bandI][0])
                /totalWaveLength_;
         }
         else
         {
             WarningInFunction
-                << "Incompatible dimensions for dQ field" << endl;
+                << "Incompatible dimensions for Qdot field" << endl;
         }
     }
 
diff --git a/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.C b/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.C
index c7c07a760f..47664c4fc9 100644
--- a/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.C
+++ b/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.C
@@ -114,15 +114,15 @@ Foam::solidChemistryModel<CompType, SolidThermo>::tc() const
 
 template<class CompType, class SolidThermo>
 Foam::tmp<Foam::volScalarField>
-Foam::solidChemistryModel<CompType, SolidThermo>::Sh() const
+Foam::solidChemistryModel<CompType, SolidThermo>::Qdot() const
 {
-    tmp<volScalarField> tSh
+    tmp<volScalarField> tQdot
     (
         new volScalarField
         (
             IOobject
             (
-                "Sh",
+                "Qdot",
                 this->mesh_.time().timeName(),
                 this->mesh_,
                 IOobject::NO_READ,
@@ -130,57 +130,25 @@ Foam::solidChemistryModel<CompType, SolidThermo>::Sh() const
                 false
             ),
             this->mesh_,
-            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
+            dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
         )
     );
 
     if (this->chemistry_)
     {
-        scalarField& Sh = tSh.ref();
+        scalarField& Qdot = tQdot.ref();
 
         forAll(Ys_, i)
         {
-            forAll(Sh, celli)
+            forAll(Qdot, celli)
             {
                 scalar hf = solidThermo_[i].Hc();
-                Sh[celli] -= hf*RRs_[i][celli];
+                Qdot[celli] -= hf*RRs_[i][celli];
             }
         }
     }
 
-    return tSh;
-}
-
-
-template<class CompType, class SolidThermo>
-Foam::tmp<Foam::volScalarField>
-Foam::solidChemistryModel<CompType, SolidThermo>::dQ() const
-{
-    tmp<volScalarField> tdQ
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                "dQ",
-                this->mesh_.time().timeName(),
-                this->mesh_,
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            this->mesh_,
-            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
-        )
-    );
-
-    if (this->chemistry_)
-    {
-        volScalarField& dQ = tdQ.ref();
-        dQ.ref() = this->mesh_.V()*Sh()();
-    }
-
-    return tdQ;
+    return tQdot;
 }
 
 
diff --git a/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.H b/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.H
index f3692027ce..0a04c4f4f8 100644
--- a/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.H
+++ b/src/thermophysicalModels/solidChemistryModel/solidChemistryModel/solidChemistryModel.H
@@ -197,11 +197,8 @@ public:
             //- Return the chemical time scale
             virtual tmp<volScalarField> tc() const;
 
-            //- Return source for enthalpy equation [kg/m/s3]
-            virtual tmp<volScalarField> Sh() const;
-
-            //- Return the heat release, i.e. enthalpy/sec [m2/s3]
-            virtual tmp<volScalarField> dQ() const;
+            //- Return the heat release rate [kg/m/s3]
+            virtual tmp<volScalarField> Qdot() const;
 
 
         // ODE functions (overriding abstract functions in ODE.H)