diff --git a/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C b/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
index 7e5da847ca..108f68fe8f 100644
--- a/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
+++ b/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2017-2022 OpenCFD Ltd.
+    Copyright (C) 2017-2023 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -32,6 +32,7 @@ License
 #include "volFields.H"
 #include "mappedPatchBase.H"
 #include "basicThermo.H"
+#include "IOField.H"
 #include "mappedPatchFieldBase.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -41,6 +42,82 @@ namespace Foam
 namespace compressible
 {
 
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+volScalarField&
+turbulentTemperatureRadCoupledMixedFvPatchScalarField::getOrCreateField
+(
+    const word& fieldName
+) const
+{
+    const fvMesh& mesh = patch().boundaryMesh().mesh();
+
+    auto* ptr = mesh.getObjectPtr<volScalarField>(fieldName);
+
+    if (!ptr)
+    {
+        ptr = new volScalarField
+        (
+            IOobject
+            (
+                fieldName,
+                mesh.time().timeName(),
+                mesh,
+                IOobject::NO_READ,
+                IOobject::AUTO_WRITE
+            ),
+            mesh,
+            dimensionedScalar(dimless, Zero)
+        );
+        mesh.objectRegistry::store(ptr);
+    }
+
+    return *ptr;
+}
+
+
+void turbulentTemperatureRadCoupledMixedFvPatchScalarField::storeHTCFields
+(
+    const word& prefix,
+    const scalarField& shtc,
+    const scalarField& shtcPatch
+)
+const
+{
+    volScalarField& htc =
+        getOrCreateField(IOobject::scopedName(prefix, "htc"));
+    htc.boundaryFieldRef()[patch().index()] = shtc;
+
+    volScalarField& htcPatch =
+        getOrCreateField(IOobject::scopedName(prefix, "htcPatch"));
+    htcPatch.boundaryFieldRef()[patch().index()] = shtcPatch;
+}
+
+
+void turbulentTemperatureRadCoupledMixedFvPatchScalarField::writeFileHeader
+(
+    Ostream& os
+)
+{
+    writeCommented(os, "Time");
+    writeTabbed(os, "Q_[W]");
+    writeTabbed(os, "q_[W/m^2]");
+    writeTabbed(os, "HTCavg_[W/m^2/K]");
+    writeTabbed(os, "patchHTCavg_[W/m^2/K]");
+    writeTabbed(os, "TpMin_[K]");
+    writeTabbed(os, "TpMax_[K]");
+    writeTabbed(os, "TpAvg_[K]");
+    writeTabbed(os, "TpNbrMin_[K]");
+    writeTabbed(os, "TpNbrMax_[K]");
+    writeTabbed(os, "TpNbrAvg_[K]");
+
+    os  << endl;
+
+    writtenHeader_ = true;
+    updateHeader_ = false;
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 turbulentTemperatureRadCoupledMixedFvPatchScalarField::
@@ -57,10 +134,21 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
         mappedPatchFieldBase<scalar>::mapper(p, iF),
         *this
     ),
+    functionObjects::writeFile
+    (
+        db(),
+        "turbulentTemperatureRadCoupledMixed",
+        "undefined",
+        false
+    ),
     TnbrName_("undefined-Tnbr"),
     qrNbrName_("undefined-qrNbr"),
     qrName_("undefined-qr"),
-    thermalInertia_(false)
+    logInterval_(-1),
+    executionIndex_(0),
+    thermalInertia_(false),
+    verbose_(false),
+    prefix_(word::null)
 {
     this->refValue() = Zero;
     this->refGrad() = Zero;
@@ -86,6 +174,7 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
         *this,
         psf
     ),
+    functionObjects::writeFile(psf),
     TnbrName_(psf.TnbrName_),
     qrNbrName_(psf.qrNbrName_),
     qrName_(psf.qrName_),
@@ -93,7 +182,11 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
     thicknessLayer_(psf.thicknessLayer_.clone(p.patch())),
     kappaLayers_(psf.kappaLayers_),
     kappaLayer_(psf.kappaLayer_.clone(p.patch())),
-    thermalInertia_(psf.thermalInertia_)
+    logInterval_(psf.logInterval_),
+    executionIndex_(psf.executionIndex_),
+    thermalInertia_(psf.thermalInertia_),
+    verbose_(psf.verbose_),
+    prefix_(psf.prefix_)
 {}
 
 
@@ -113,10 +206,21 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
         *this,
         dict
     ),
+    functionObjects::writeFile
+    (
+        db(),
+        "turbulentTemperatureRadCoupledMixed",
+        patch().name(),
+        false
+    ),
     TnbrName_(dict.getOrDefault<word>("Tnbr", "T")),
     qrNbrName_(dict.getOrDefault<word>("qrNbr", "none")),
     qrName_(dict.getOrDefault<word>("qr", "none")),
-    thermalInertia_(dict.getOrDefault<Switch>("thermalInertia", false))
+    logInterval_(dict.getOrDefault<scalar>("logInterval", -1)),
+    executionIndex_(0),
+    thermalInertia_(dict.getOrDefault<Switch>("thermalInertia", false)),
+    verbose_(dict.getOrDefault<bool>("verbose", false)),
+    prefix_(dict.getOrDefault<word>("prefix", "multiWorld"))
 {
     if (!isA<mappedPatchBase>(this->patch().patch()))
     {
@@ -218,6 +322,8 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
     {
         source() = 0.0;
     }
+
+    writeFile::read(dict);
 }
 
 
@@ -236,6 +342,7 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
         *this,
         psf
     ),
+    functionObjects::writeFile(psf),
     TnbrName_(psf.TnbrName_),
     qrNbrName_(psf.qrNbrName_),
     qrName_(psf.qrName_),
@@ -243,7 +350,11 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
     thicknessLayer_(psf.thicknessLayer_.clone(patch().patch())),
     kappaLayers_(psf.kappaLayers_),
     kappaLayer_(psf.kappaLayer_.clone(patch().patch())),
-    thermalInertia_(psf.thermalInertia_)
+    logInterval_(psf.logInterval_),
+    executionIndex_(psf.executionIndex_),
+    thermalInertia_(psf.thermalInertia_),
+    verbose_(psf.verbose_),
+    prefix_(psf.prefix_)
 {}
 
 
@@ -261,6 +372,7 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
         *this,
         psf
     ),
+    functionObjects::writeFile(psf),
     TnbrName_(psf.TnbrName_),
     qrNbrName_(psf.qrNbrName_),
     qrName_(psf.qrName_),
@@ -268,7 +380,11 @@ turbulentTemperatureRadCoupledMixedFvPatchScalarField
     thicknessLayer_(psf.thicknessLayer_.clone(patch().patch())),
     kappaLayers_(psf.kappaLayers_),
     kappaLayer_(psf.kappaLayer_.clone(patch().patch())),
-    thermalInertia_(psf.thermalInertia_)
+    logInterval_(psf.logInterval_),
+    executionIndex_(psf.executionIndex_),
+    thermalInertia_(psf.thermalInertia_),
+    verbose_(psf.verbose_),
+    prefix_(psf.prefix_)
 {}
 
 
@@ -358,6 +474,7 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
 
 
     scalarField TcNbr;
+    scalarField TpNbr;  // only if verbose_
     scalarField KDeltaNbr;
 
     if (mpp.sameWorld())
@@ -375,6 +492,10 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
 
         // Swap to obtain full local values of neighbour K*delta
         TcNbr = nbrField.patchInternalField();
+        if (verbose_)
+        {
+            TpNbr = nbrField;
+        }
         KDeltaNbr = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
     }
     else
@@ -382,9 +503,17 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
         // Different world so use my region,patch. Distribution below will
         // do the reordering.
         TcNbr = patchInternalField();
+        if (verbose_)
+        {
+            TpNbr = Tp;
+        }
         KDeltaNbr = KDelta;
     }
     distribute(this->internalField().name() + "_value", TcNbr);
+    if (verbose_)
+    {
+        distribute(this->internalField().name() + "_patchValue", TpNbr);
+    }
     distribute(this->internalField().name() + "_weights", KDeltaNbr);
 
     scalarField KDeltaC(this->size(), GREAT);
@@ -545,22 +674,111 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
 
     mixedFvPatchScalarField::updateCoeffs();
 
-    if (debug)
-    {
-        scalar Q = gSum(kappaTp*patch().magSf()*snGrad());
 
-        Info<< patch().boundaryMesh().mesh().name() << ':'
+    if (verbose_)
+    {
+        // Calculate heat-transfer rate and heat flux
+        const scalar Q = gSum(kappaTp*patch().magSf()*snGrad());
+        const scalar magSf = gSum(patch().magSf());
+        const scalar q = Q/max(magSf, SMALL);
+
+
+        // Calculate heat-transfer coeff based on the first definition
+        // [W/m^2] = [W/m/K K * 1/m]
+        const scalarField qField
+        (
+            kappaTp*snGrad()
+        );
+        const scalarField deltaT(TcNbr - Tc);
+        scalarField htc(deltaT.size(), Zero);
+
+        forAll(deltaT, i)
+        {
+            if (mag(deltaT[i]) > SMALL)
+            {
+                htc[i] = qField[i]/deltaT[i];
+            }
+        }
+        const scalar aveHtc = gSum(htc*patch().magSf())/max(magSf, SMALL);
+
+
+        // Calculate heat-transfer coeff based on the second definition
+        const scalarField deltaTPatch(TpNbr - Tp);
+        scalarField htcPatch(deltaTPatch.size(), Zero);
+
+        forAll(deltaTPatch, i)
+        {
+            if (mag(deltaTPatch[i]) > SMALL)
+            {
+                htcPatch[i] = qField[i]/deltaTPatch[i];
+            }
+        }
+        const scalar aveHtcPatch =
+            gSum(htcPatch*patch().magSf())/max(magSf, SMALL);
+
+        // Calculate various averages of temperature
+        const scalarMinMax TpMinMax = gMinMax(Tp);
+        const scalar TpAvg = gAverage(Tp);
+        const scalarMinMax TpNbrMinMax = gMinMax(TpNbr);
+        const scalar TpNbrAvg = gAverage(TpNbr);
+
+
+        Info<< nl
+            << patch().boundaryMesh().mesh().name() << ':'
             << patch().name() << ':'
             << this->internalField().name() << " <- "
             << mpp.sampleRegion() << ':'
             << mpp.samplePatch() << ':'
-            << this->internalField().name() << " :"
-            << " heat transfer rate:" << Q
-            << " walltemperature "
-            << " min:" << gMin(Tp)
-            << " max:" << gMax(Tp)
-            << " avg:" << gAverage(Tp)
-            << endl;
+            << this->internalField().name() << " :" << nl
+            << " Heat transfer rate [W]:" << Q << nl
+            << " Area [m^2]:" << magSf << nl
+            << " Heat flux [W/m^2]:" << q << nl
+            << " Area-averaged heat-transfer coefficient [W/m^2/K]:"
+            << aveHtc << nl
+            << " Area-averaged patch heat-transfer coefficient [W/m^2/K]:"
+            << aveHtcPatch << nl
+            << " Wall temperature [K]"
+            << " min:" << TpMinMax.min()
+            << " max:" << TpMinMax.max()
+            << " avg:" << TpAvg << nl
+            << " Neighbour wall temperature [K]"
+            << " min:" << TpNbrMinMax.min()
+            << " max:" << TpNbrMinMax.max()
+            << " avg:" << TpNbrAvg
+            << nl << endl;
+
+
+        // Handle data for file output
+        if (canResetFile())
+        {
+            resetFile(patch().name());
+        }
+
+        if (canWriteHeader())
+        {
+            writeFileHeader(file());
+        }
+
+        if (canWriteToFile() && writeFile())
+        {
+            file()
+                << db().time().timeOutputValue() << token::TAB
+                << Q << token::TAB
+                << q << token::TAB
+                << aveHtc << token::TAB
+                << aveHtcPatch << token::TAB
+                << TpMinMax.min() << token::TAB
+                << TpMinMax.max() << token::TAB
+                << TpAvg << token::TAB
+                << TpNbrMinMax.min() << token::TAB
+                << TpNbrMinMax.max() << token::TAB
+                << TpNbrAvg << token::TAB
+                << endl;
+        }
+
+
+        // Store htc fields as patch fields of a volScalarField
+        storeHTCFields(prefix_, htc, htcPatch);
     }
 
     // Restore tag
@@ -713,6 +931,39 @@ deltaH() const
 }
 
 
+bool turbulentTemperatureRadCoupledMixedFvPatchScalarField::writeFile()
+{
+    if (!verbose_ || (logInterval_ <= 0))
+    {
+        return false;
+    }
+
+    const auto& time = patch().boundaryMesh().mesh().time();
+
+    const scalar t = time.timeOutputValue();
+    const scalar ts = time.startTime().value();
+    const scalar deltaT = time.deltaTValue();
+
+    const label executionIndex = label
+    (
+        (
+            (t - ts)
+          + 0.5*deltaT
+        )
+        /logInterval_
+    );
+
+    bool write = false;
+    if (executionIndex > executionIndex_)
+    {
+        executionIndex_ = executionIndex;
+        write = true;
+    }
+
+    return write;
+}
+
+
 void turbulentTemperatureRadCoupledMixedFvPatchScalarField::write
 (
     Ostream& os
@@ -720,17 +971,17 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::write
 {
     mixedFvPatchField<scalar>::write(os);
 
-    //os.writeEntry("Tnbr", TnbrName_);
     os.writeEntryIfDifferent<word>("Tnbr", "T", TnbrName_);
-
-    //os.writeEntry("qrNbr", qrNbrName_);
     os.writeEntryIfDifferent<word>("qrNbr", "none", qrNbrName_);
-    //os.writeEntry("qr", qrName_);
     os.writeEntryIfDifferent<word>("qr", "none", qrName_);
+    os.writeEntry<scalar>("logInterval", logInterval_);
+
     if (thermalInertia_)
     {
         os.writeEntry("thermalInertia", thermalInertia_);
     }
+    os.writeEntryIfDifferent<bool>("verbose", false, verbose_);
+    os.writeEntryIfDifferent<word>("prefix", "multiWorld", prefix_);
 
     if (thicknessLayer_)
     {
@@ -743,6 +994,12 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::write
         kappaLayers_.writeEntry("kappaLayers", os);
     }
 
+    // Write writeFile entries
+    os.writeEntry<label>("writePrecision", writePrecision_);
+    os.writeEntry<bool>("updateHeader", updateHeader_);
+    os.writeEntry<bool>("writeToFile", writeToFile_);
+    os.writeEntry<bool>("useUserTime", useUserTime_);
+
     temperatureCoupledBase::write(os);
     mappedPatchFieldBase<scalar>::write(os);
 }
diff --git a/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.H b/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.H
index 0dbabd4d7b..44a721daa2 100644
--- a/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.H
+++ b/src/thermoTools/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.H
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2017-2021 OpenCFD Ltd.
+    Copyright (C) 2017-2023 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -25,41 +25,59 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::compressible::
-        turbulentTemperatureRadCoupledMixedFvPatchScalarField
+    Foam::compressible::turbulentTemperatureRadCoupledMixedFvPatchScalarField
 
 Description
-    Mixed boundary condition for temperature and radiation heat transfer
-    to be used for in multiregion cases. Optional thin thermal layer
-    resistances can be specified through thicknessLayers and kappaLayers
-    entries.
+    Mixed boundary condition for temperature and radiation heat transfer,
+    suitable for multiregion cases. The optional specification of
+    thin thermal layer resistances can be made using the entries
+    \c thicknessLayers and \c kappaLayers.
 
-    The thermal conductivity \c kappa can either be retrieved from various
-    possible sources, as detailed in the class temperatureCoupledBase.
+    The thermal conductivity \c kappa can be obtained from a variety
+    of sources, as explained in the \c temperatureCoupledBase class.
 
 Usage
-    \table
-        Property     | Description             | Required    | Default value
-        Tnbr         | name of the field    | no | T
-        qrNbr      | name of the radiative flux in the nbr region | no | none
-        qr         | name of the radiative flux in this region | no | none
-        thicknessLayers | list of thicknesses per layer [m] | no |
-        kappaLayers  | list of thermal conductivites per layer [W/m/K] | no |
-        thicknessLayer | single thickness of layer [m] | no |
-        kappaLayer  | corresponding thermal conductivity [W/m/K] | no |
-        kappaMethod  | inherited from temperatureCoupledBase | inherited |
-        kappa        | inherited from temperatureCoupledBase | inherited |
-        thermalInertia | Add thermal inertia to wall node | no | false
-    \endtable
-
     Example of the boundary condition specification:
     \verbatim
     <patchName>
+    {
+        // Mandatory entries
+        type                compressible::turbulentTemperatureRadCoupledMixed;
+
+        // Optional entries
+        Tnbr                <word>;
+        qrNbr               <word>;
+        qr                  <word>;
+        logInterval         <scalar>;
+        thermalInertia      <bool>;
+        verbose             <bool>;
+        prefix              <word>;
+
+        // Conditional entries
+
+            // Option-1
+            thicknessLayers <scalarList>;
+            kappaLayers     <scalarList>;
+
+            // Option-2
+            thicknessLayer  <PatchFunction1<scalar>>;
+            kappaLayer      <PatchFunction1<scalar>>;
+
+        // Inherited entries
+        kappaMethod         <word>;
+        kappa               <word>;
+        ...
+    }
+    \endverbatim
+
+    For example:
+    \verbatim
+    <patchName>
     {
         type            compressible::turbulentTemperatureRadCoupledMixed;
         Tnbr            T;
         qrNbr           qr; // or none. Name of qr field on neighbour region
-        qr              qr; // or none. Name of qr field on local region
+        qr              qr; // or none. Name of qr field on this region
         thicknessLayers (0.1 0.2 0.3 0.4);
         kappaLayers     (1 2 3 4);
         thermalInertia  false/true;
@@ -69,23 +87,53 @@ Usage
     }
     \endverbatim
 
-    Needs to be on underlying mapped(Wall)FvPatch.
+    where the entries mean:
+    \table
+      Property | Description                            | Type  | Reqd | Deflt
+      type   | Type name: compressible::turbulentTemperatureRadCoupledMixed <!--
+               -->                                      | word  | yes  | -
+      Tnbr     | Name of the nbr temperature field      | word  | no   | T
+      qrNbr    | Name of the radiative flux in the nbr region | word | no | none
+      qr       | Name of the radiative flux in this region | word | no | none
+      logInterval | Log-file output frequency [s]          | scalar | no  | -1
+      thermalInertia | Flag to add thermal inertia to wall node | bool | no <!--
+               -->                                                     | false
+      verbose | Flag to enable verbose output with extra fields | bool | no <!--
+              -->                                                      | false
+      prefix  | Name of output-field prefix (in verbose mode)         <!--
+              -->                                       | word | no | multiWorld
+      thicknessLayers | List of thicknesses per layer [m] | scalarList <!--
+              -->                                             | choice | -
+      kappaLayers  | List of thermal conductivites per layer [W/m/K]   <!--
+              -->                                | scalarList | choice | -
+      thicknessLayer | Single thickness of layer [m]                   <!--
+              -->                   | PatchFunction1\<scalar\>| choice | -
+      kappaLayer     | Corresponding thermal conductivity [W/m/K]      <!--
+              -->                   | PatchFunction1\<scalar\>| choice | -
+    \endtable
 
-See also
-    Foam::temperatureCoupledBase
+    The inherited entries are elaborated in:
+      - \link mixedFvPatchFields.H \endlink
+      - \link temperatureCoupledBase.H \endlink
+      - \link mappedPatchFieldBase.H \endlink
+      - \link writeFile.H \endlink
+
+Note
+  - Needs to be on underlying \c mapped(Wall)FvPatch.
 
 SourceFiles
     turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef turbulentTemperatureRadCoupledMixedFvPatchScalarField_H
-#define turbulentTemperatureRadCoupledMixedFvPatchScalarField_H
+#ifndef Foam_turbulentTemperatureRadCoupledMixedFvPatchScalarField_H
+#define Foam_turbulentTemperatureRadCoupledMixedFvPatchScalarField_H
 
 #include "mixedFvPatchFields.H"
 #include "temperatureCoupledBase.H"
 #include "scalarList.H"
 #include "mappedPatchFieldBase.H"
+#include "writeFile.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -102,18 +150,18 @@ class turbulentTemperatureRadCoupledMixedFvPatchScalarField
 :
     public mixedFvPatchScalarField,
     public temperatureCoupledBase,
-    public mappedPatchFieldBase<scalar>
+    public mappedPatchFieldBase<scalar>,
+    public functionObjects::writeFile
 {
-
     // Private Data
 
-        //- Name of field on the neighbour region
+        //- Name of temperature field on the neighbour region
         const word TnbrName_;
 
         //- Name of the radiative heat flux in the neighbour region
         const word qrNbrName_;
 
-        //- Name of the radiative heat flux in local region
+        //- Name of the radiative heat flux in this region
         const word qrName_;
 
         //- Thickness of layers (either scalarList or a single PatchFunction1)
@@ -124,29 +172,61 @@ class turbulentTemperatureRadCoupledMixedFvPatchScalarField
         scalarList kappaLayers_;
         autoPtr<PatchFunction1<scalar>> kappaLayer_;
 
-        //- Thermal inertia term
-        Switch thermalInertia_;
+        //- Log-file output frequency
+        scalar logInterval_;
+
+        //- Index of previous execution
+        label executionIndex_;
+
+        //- Flag to enable thermal inertia term
+        const Switch thermalInertia_;
+
+        //- Flag to enable verbose output with extra fields
+        const bool verbose_;
+
+        //- Name of output-field prefix (in verbose mode)
+        const word prefix_;
 
 
      // Private Functions
 
-            //- Return local alphaSfDelta
-            tmp<scalarField> alphaSfDelta() const;
+        //- Return local alphaSfDelta
+        tmp<scalarField> alphaSfDelta() const;
 
-            //- Return delta enthalpy between regions
-            tmp<scalarField> deltaH() const;
+        //- Return delta enthalpy between regions
+        tmp<scalarField> deltaH() const;
 
-            //- Return the sum of deltaCoeff*alpha from nbr and local
-            tmp<scalarField> beta() const;
+        //- Return the sum of deltaCoeff*alpha from nbr and local
+        tmp<scalarField> beta() const;
+
+        //- Calculate coefficients for assembly matrix
+        tmp<Field<scalar>> coeffs
+        (
+            fvMatrix<scalar>& matrix,
+            const Field<scalar>&,
+            const label
+        ) const;
+
+        //- Return requested field from the object registry
+        //- or create+register the field to the object registry
+        volScalarField& getOrCreateField(const word& fieldName) const;
+
+        //- Store heat-transfer coeff fields to facilitate postprocessing
+        void storeHTCFields
+        (
+            const word& prefix,
+            const scalarField& htc,
+            const scalarField& htcPatch
+        ) const;
+
+        //- Flag to allow writing to file
+        bool writeFile();
 
 
-            //- Calculate coefficients for assembly matrix
-            tmp<Field<scalar>> coeffs
-            (
-                fvMatrix<scalar>& matrix,
-                const Field<scalar>&,
-                const label
-            ) const;
+        // functionObjects::writeFile
+
+            //- Write file header information
+            void writeFileHeader(Ostream& os);
 
 
 public:
@@ -173,8 +253,8 @@ public:
         );
 
         //- Construct by mapping given
-        //  turbulentTemperatureCoupledBaffleMixedFvPatchScalarField onto a
-        //  new patch
+        //- turbulentTemperatureCoupledBaffleMixedFvPatchScalarField onto a
+        //- new patch
         turbulentTemperatureRadCoupledMixedFvPatchScalarField
         (
             const
@@ -228,22 +308,24 @@ public:
 
     // Member Functions
 
-        // Mapping functions
+    // Mapping
 
-            //- Map (and resize as needed) from self given a mapping object
-            virtual void autoMap
-            (
-                const fvPatchFieldMapper&
-            );
+        //- Map (and resize as needed) from self given a mapping object
+        virtual void autoMap
+        (
+            const fvPatchFieldMapper&
+        );
 
-            //- Reverse map the given fvPatchField onto this fvPatchField
-            virtual void rmap
-            (
-                const fvPatchField<scalar>&,
-                const labelList&
-            );
+        //- Reverse map the given fvPatchField onto this fvPatchField
+        virtual void rmap
+        (
+            const fvPatchField<scalar>&,
+            const labelList&
+        );
 
 
+    // Evaluation
+
         //- Given patch temperature calculate corresponding K field. Override
         //- temperatureCoupledBase::kappa to includes effect of any
         //- explicit kappaThickness
@@ -261,6 +343,8 @@ public:
         );
 
 
+    // I-O
+
         //- Write
         virtual void write(Ostream& os) const;
 };