diff --git a/src/thermophysicalModels/radiation/radiationModel/P1/P1.C b/src/thermophysicalModels/radiation/radiationModel/P1/P1.C
index 1afa7a6f0..2cc76c07b 100644
--- a/src/thermophysicalModels/radiation/radiationModel/P1/P1.C
+++ b/src/thermophysicalModels/radiation/radiationModel/P1/P1.C
@@ -39,13 +39,7 @@ namespace Foam
     namespace radiation
     {
         defineTypeNameAndDebug(P1, 0);
-
-        addToRunTimeSelectionTable
-        (
-            radiationModel,
-            P1,
-            dictionary
-        );
+        addToRadiationRunTimeSelectionTables(P1);
     }
 }
 
@@ -122,6 +116,76 @@ Foam::radiation::P1::P1(const volScalarField& T)
 {}
 
 
+Foam::radiation::P1::P1(const dictionary& dict, const volScalarField& T)
+:
+    radiationModel(typeName, dict, T),
+    G_
+    (
+        IOobject
+        (
+            "G",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_
+    ),
+    Qr_
+    (
+        IOobject
+        (
+            "Qr",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("Qr", dimMass/pow3(dimTime), 0.0)
+    ),
+    a_
+    (
+        IOobject
+        (
+            "a",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("a", dimless/dimLength, 0.0)
+    ),
+    e_
+    (
+        IOobject
+        (
+            "e",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("a", dimless/dimLength, 0.0)
+    ),
+    E_
+    (
+        IOobject
+        (
+            "E",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("E", dimMass/dimLength/pow3(dimTime), 0.0)
+    )
+{}
+
+
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 Foam::radiation::P1::~P1()
@@ -150,7 +214,9 @@ void Foam::radiation::P1::calculate()
     a_ = absorptionEmission_->a();
     e_ = absorptionEmission_->e();
     E_ = absorptionEmission_->E();
-    const volScalarField sigmaEff = scatter_->sigmaEff();
+    const volScalarField sigmaEff(scatter_->sigmaEff());
+
+    const dimensionedScalar a0 ("a0", a_.dimensions(), ROOTVSMALL);
 
     // Construct diffusion
     const volScalarField gamma
@@ -163,7 +229,7 @@ void Foam::radiation::P1::calculate()
             IOobject::NO_READ,
             IOobject::NO_WRITE
         ),
-        1.0/(3.0*a_ + sigmaEff)
+        1.0/(3.0*a_ + sigmaEff + a0)
     );
 
     // Solve G transport equation
@@ -210,10 +276,8 @@ Foam::radiation::P1::Ru() const
 {
     const DimensionedField<scalar, volMesh>& G =
         G_.dimensionedInternalField();
-
     const DimensionedField<scalar, volMesh> E =
         absorptionEmission_->ECont()().dimensionedInternalField();
-
     const DimensionedField<scalar, volMesh> a =
         absorptionEmission_->aCont()().dimensionedInternalField();
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/P1/P1.H b/src/thermophysicalModels/radiation/radiationModel/P1/P1.H
index 5eeffb6ec..73321e92c 100644
--- a/src/thermophysicalModels/radiation/radiationModel/P1/P1.H
+++ b/src/thermophysicalModels/radiation/radiationModel/P1/P1.H
@@ -96,6 +96,8 @@ public:
         //- Construct from components
         P1(const volScalarField& T);
 
+        //- Construct from components
+        P1(const dictionary& dict, const volScalarField& T);
 
     // Destructor
     virtual ~P1();
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
index 880383e38..e856a3f5c 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
@@ -38,62 +38,14 @@ namespace Foam
     namespace radiation
     {
         defineTypeNameAndDebug(fvDOM, 0);
-
-        addToRunTimeSelectionTable
-        (
-            radiationModel,
-            fvDOM,
-            dictionary
-        );
+        addToRadiationRunTimeSelectionTables(fvDOM);
     }
 }
 
 
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
-:
-    radiationModel(typeName, T),
-    G_
-    (
-        IOobject
-        (
-            "G",
-            mesh().time().timeName(),
-            T.db(),
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh(),
-        dimensionedScalar("G", dimMass/pow3(dimTime), 0.0)
-    ),
-    a_
-    (
-        IOobject
-        (
-            "a",
-            mesh().time().timeName(),
-            T.db(),
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh(),
-        dimensionedScalar("a", dimless/dimLength, 0.0)
-    ),
-    nTheta_(readLabel(coeffs_.lookup("nTheta"))),
-    nPhi_(readLabel(coeffs_.lookup("nPhi"))),
-    nRay_(0),
-    nLambda_(absorptionEmission_->nBands()),
-    aLambda_(nLambda_),
-    blackBody_(nLambda_, T),
-    IRay_(0),
-    Qem_(nLambda_),
-    Qin_(nLambda_),
-    convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
-    maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
-    fvRayDiv_(nLambda_),
-    cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
-    omegaMax_(0)
+void Foam::radiation::fvDOM::initialise()
 {
     if (mesh().nSolutionD() == 3)    //3D
     {
@@ -208,7 +160,7 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
                 (
                     "aLambda_" + Foam::name(lambdaI) ,
                     mesh().time().timeName(),
-                    T.db(),
+                    T_.db(),
                     IOobject::NO_READ,
                     IOobject::NO_WRITE
                 ),
@@ -260,8 +212,7 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
             forAll(IRay_, rayId)
             {
                 const surfaceScalarField Ji(IRay_[rayId].dAve() & mesh().Sf());
-                volScalarField& iRayLambdaI =
-                    IRay_[rayId].ILambda(lambdaI);
+                volScalarField& iRayLambdaI = IRay_[rayId].ILambda(lambdaI);
 
                 fvRayDiv_[lambdaI].set
                 (
@@ -288,6 +239,194 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
     Info<< endl;
 }
 
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
+:
+    radiationModel(typeName, T),
+    G_
+    (
+        IOobject
+        (
+            "G",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("G", dimMass/pow3(dimTime), 0.0)
+    ),
+    Qr_
+    (
+        IOobject
+        (
+            "Qr",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::READ_IF_PRESENT,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("Qr", dimMass/pow3(dimTime), 0.0)
+    ),
+    a_
+    (
+        IOobject
+        (
+            "a",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("a", dimless/dimLength, 0.0)
+    ),
+    nTheta_(readLabel(coeffs_.lookup("nTheta"))),
+    nPhi_(readLabel(coeffs_.lookup("nPhi"))),
+    nRay_(0),
+    nLambda_(absorptionEmission_->nBands()),
+    aLambda_(nLambda_),
+    blackBody_(nLambda_, T),
+    IRay_(0),
+    Qem_(nLambda_),
+    Qin_(nLambda_),
+    convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
+    maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
+    fvRayDiv_(nLambda_),
+    cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
+    omegaMax_(0)
+{
+    initialise();
+}
+
+
+Foam::radiation::fvDOM::fvDOM(const word& type, const volScalarField& T)
+:
+    radiationModel(type, T),
+    G_
+    (
+        IOobject
+        (
+            "G",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("G", dimMass/pow3(dimTime), 0.0)
+    ),
+    Qr_
+    (
+        IOobject
+        (
+            "Qr",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::READ_IF_PRESENT,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("Qr", dimMass/pow3(dimTime), 0.0)
+    ),
+    a_
+    (
+        IOobject
+        (
+            "a",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("a", dimless/dimLength, 0.0)
+    ),
+    nTheta_(readLabel(coeffs_.lookup("nTheta"))),
+    nPhi_(readLabel(coeffs_.lookup("nPhi"))),
+    nRay_(0),
+    nLambda_(absorptionEmission_->nBands()),
+    aLambda_(nLambda_),
+    blackBody_(nLambda_, T),
+    IRay_(0),
+    Qem_(nLambda_),
+    Qin_(nLambda_),
+    convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
+    maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
+    fvRayDiv_(nLambda_),
+    cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
+    omegaMax_(0)
+{
+    initialise();
+}
+
+
+Foam::radiation::fvDOM::fvDOM
+(
+    const dictionary& dict,
+    const volScalarField& T
+)
+:
+    radiationModel(typeName, dict, T),
+    G_
+    (
+        IOobject
+        (
+            "G",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("G", dimMass/pow3(dimTime), 0.0)
+    ),
+    Qr_
+    (
+        IOobject
+        (
+            "Qr",
+            mesh_.time().timeName(),
+            T.db(),
+            IOobject::READ_IF_PRESENT,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_,
+        dimensionedScalar("Qr", dimMass/pow3(dimTime), 0.0)
+    ),
+    a_
+    (
+        IOobject
+        (
+            "a",
+            mesh().time().timeName(),
+            T.db(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh(),
+        dimensionedScalar("a", dimless/dimLength, 0.0)
+    ),
+    nTheta_(readLabel(coeffs_.lookup("nTheta"))),
+    nPhi_(readLabel(coeffs_.lookup("nPhi"))),
+    nRay_(0),
+    nLambda_(absorptionEmission_->nBands()),
+    aLambda_(nLambda_),
+    blackBody_(nLambda_, T),
+    IRay_(0),
+    Qem_(nLambda_),
+    Qin_(nLambda_),
+    convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
+    maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
+    fvRayDiv_(nLambda_),
+    cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
+    omegaMax_(0)
+{
+    initialise();
+}
+
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
@@ -346,53 +485,6 @@ Foam::radiation::fvDOM::Qem() const
     return tsumQem;
 }
 
-Foam::tmp<Foam::volScalarField::GeometricBoundaryField>
-Foam::radiation::fvDOM::Qr() const
-{
-    tmp<volScalarField::GeometricBoundaryField> tQr
-    (
-        new volScalarField::GeometricBoundaryField
-        (
-            mesh().boundary(),
-            mesh().V(),           // Dummy internal field,
-            calculatedFvPatchScalarField::typeName
-        )
-    );
-    volScalarField::GeometricBoundaryField& Qr = tQr();
-
-    Qr = 0;
-
-    forAll(Qem_, lambdaI)
-    {
-        Qr += Qin(lambdaI);
-        Qr -= Qem(lambdaI);
-    }
-
-    return tQr;
-}
-
-
-Foam::tmp<Foam::volScalarField::GeometricBoundaryField>
-Foam::radiation::fvDOM::Qr(scalar lambdaI) const
-{
-    tmp<volScalarField::GeometricBoundaryField> tQr
-    (
-        new volScalarField::GeometricBoundaryField
-        (
-            mesh().boundary(),
-            mesh().V(),           // Dummy internal field,
-            calculatedFvPatchScalarField::typeName
-        )
-    );
-    volScalarField::GeometricBoundaryField& Qr = tQr();
-
-    Qr = 0;
-    Qr += Qin(lambdaI);
-    Qr -= Qem(lambdaI);
-
-    return tQr;
-}
-
 
 bool Foam::radiation::fvDOM::read()
 {
@@ -562,12 +654,21 @@ void Foam::radiation::fvDOM::updateBlackBodyEmission()
 void Foam::radiation::fvDOM::updateG()
 {
     G_ = dimensionedScalar("zero",dimMass/pow3(dimTime), 0.0);
+    Qr_ = dimensionedScalar("zero",dimMass/pow3(dimTime), 0.0);
 
     forAll(IRay_, rayI)
     {
         IRay_[rayI].addIntensity();
         G_ += IRay_[rayI].I()*IRay_[rayI].omega();
     }
+
+    forAll(Qr_.boundaryField(), patchI)
+    {
+        forAll(aLambda_, lambdaI)
+        {
+            Qr_.boundaryField()[patchI] += Qin_[lambdaI][patchI] - Qem_[lambdaI][patchI];
+        }
+    }
 }
 
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
index f77b413a9..c5a531dc8 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
@@ -80,6 +80,9 @@ class fvDOM
         //- Incident radiation  [W/m2]
         volScalarField G_;
 
+        //- Total radiative heat flux [W/m2]
+        volScalarField Qr_;
+
         //- Total absorption coefficient [1/m]
         volScalarField a_;
 
@@ -128,6 +131,9 @@ class fvDOM
 
     // Private member functions
 
+        //- Initialise
+        void initialise();
+
         //- Disallow default bitwise copy construct
         fvDOM(const fvDOM&);
 
@@ -149,6 +155,12 @@ public:
         //- Construct from components
         fvDOM(const volScalarField& T);
 
+        //- Construct from components
+        fvDOM(const word& type, const volScalarField& T);
+
+        //- Construct from components
+        fvDOM(const dictionary& dict, const volScalarField& T);
+
 
     //- Destructor
     virtual ~fvDOM();
@@ -156,41 +168,31 @@ public:
 
     // Member functions
 
-        //- Total incident radiative heat flux field
-        tmp<volScalarField::GeometricBoundaryField> Qin() const;
+        // Edit
 
-        //- Total emitted radiative heat flux field
-        tmp<volScalarField::GeometricBoundaryField> Qem() const;
+            //- Solve radiation equation(s)
+            void calculate();
 
-        //- Total radiative heat flux as difference of Qem and Qin
-        tmp<volScalarField::GeometricBoundaryField> Qr() const;
+            //- Read radiation properties dictionary
+            virtual bool read();
 
-        //- Band-wise radiative heat flux as difference of Qem and Qin
-        tmp<volScalarField::GeometricBoundaryField> Qr(scalar lambda) const;
+            //- Update G and calculate total heat flux on boundary
+            void updateG();
 
-        //- Solve radiation equation(s)
-        void calculate();
+            //- Set the rayId and lambdaId from by decomposing an intensity
+            //  field name
+            void setRayIdLambdaId
+            (
+                const word& name,
+                label& rayId,
+                label& lambdaId
+            ) const;
 
-        //- Read radiation properties dictionary
-        bool read();
+            //- Source term component (for power of T^4)
+            virtual tmp<volScalarField> Rp() const;
 
-        //- Update G and calculate total heat flux on boundary
-        void updateG();
-
-        //- Set the rayId and lambdaId from by decomposing an intensity
-        //  field name
-        void setRayIdLambdaId
-        (
-            const word& name,
-            label& rayId,
-            label& lambdaId
-        ) const;
-
-        //- Source term component (for power of T^4)
-        virtual tmp<volScalarField> Rp() const;
-
-        //- Source term component (constant)
-        virtual tmp<DimensionedField<scalar, volMesh> > Ru() const;
+            //- Source term component (constant)
+            virtual tmp<DimensionedField<scalar, volMesh> > Ru() const;
 
 
         // Access
@@ -226,18 +228,27 @@ public:
             //- Const access to incident radiation field
             inline const volScalarField& G() const;
 
+            //- Const access to total radiative heat flux field
+            inline const volScalarField& Qr() const;
+
             //- Const access to band wise emitted radiative heat flux field
             inline const volScalarField::GeometricBoundaryField& Qem
             (
-                scalar lambda
+                const label lambda
             ) const;
 
             //- Const access to band wise incident radiative heat flux field
             inline const volScalarField::GeometricBoundaryField& Qin
             (
-                scalar lambda
+                const label lambda
             ) const;
 
+            //- Total incident radiative heat flux field
+            tmp<volScalarField::GeometricBoundaryField> Qin() const;
+
+            //- Total emitted radiative heat flux field
+            tmp<volScalarField::GeometricBoundaryField> Qem() const;
+
             //- Const access to black body
             inline const blackBodyEmission& blackBody() const;
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOMI.H b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOMI.H
index abf70e5ea..04f65dee7 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOMI.H
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOMI.H
@@ -29,6 +29,7 @@ Foam::radiation::fvDOM::IRay(const label rayI) const
     return  IRay_[rayI];
 }
 
+
 inline const Foam::volScalarField&
 Foam::radiation::fvDOM::IRayLambda
 (
@@ -85,34 +86,26 @@ inline const Foam::volScalarField& Foam::radiation::fvDOM::G() const
 }
 
 
+inline const Foam::volScalarField& Foam::radiation::fvDOM::Qr() const
+{
+    return Qr_;
+}
+
+
 inline const Foam::volScalarField::GeometricBoundaryField&
-Foam::radiation::fvDOM::Qin(scalar lambdaI) const
+Foam::radiation::fvDOM::Qin(const label lambdaI) const
 {
     return Qin_[lambdaI];
 }
 
 
-/*
-inline Foam::volScalarField::GeometricBoundaryField&
-Foam::radiation::fvDOM::Qin(scalar lambdaI)
-{
-    return Qin_[lambdaI];
-}
-*/
 inline const Foam::volScalarField::GeometricBoundaryField&
-Foam::radiation::fvDOM::Qem(scalar lambdaI) const
+Foam::radiation::fvDOM::Qem(const label lambdaI) const
 {
     return Qem_[lambdaI];
 }
 
 
-/*
-inline Foam::volScalarField::GeometricBoundaryField&
-Foam::radiation::fvDOM::Qem(scalar lambdaI)
-{
-    return Qem_[lambdaI];
-}
-*/
 inline const Foam::radiation::blackBodyEmission&
 Foam::radiation::fvDOM::blackBody() const
 {
diff --git a/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.C b/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.C
index 5fc61a5cc..9a8f53057 100644
--- a/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.C
+++ b/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.C
@@ -36,13 +36,7 @@ namespace Foam
     namespace radiation
     {
         defineTypeNameAndDebug(noRadiation, 0);
-
-        addToRunTimeSelectionTable
-        (
-            radiationModel,
-            noRadiation,
-            dictionary
-        );
+        addToRadiationRunTimeSelectionTables(noRadiation);
     }
 }
 
@@ -55,6 +49,16 @@ Foam::radiation::noRadiation::noRadiation(const volScalarField& T)
 {}
 
 
+Foam::radiation::noRadiation::noRadiation
+(
+    const dictionary& dict,
+    const volScalarField& T
+)
+:
+    radiationModel(T)
+{}
+
+
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 Foam::radiation::noRadiation::~noRadiation()
diff --git a/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.H b/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.H
index 778367031..77a2843a8 100644
--- a/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.H
+++ b/src/thermophysicalModels/radiation/radiationModel/noRadiation/noRadiation.H
@@ -73,6 +73,9 @@ public:
         //- Construct from components
         noRadiation(const volScalarField& T);
 
+        //- Construct from components
+        noRadiation(const dictionary& dict, const volScalarField& T);
+
 
     //- Destructor
     virtual ~noRadiation();
diff --git a/src/thermophysicalModels/radiation/radiationModel/radiationModel/newRadiationModel.C b/src/thermophysicalModels/radiation/radiationModel/radiationModel/newRadiationModel.C
index 012f5d134..95896c0a0 100644
--- a/src/thermophysicalModels/radiation/radiationModel/radiationModel/newRadiationModel.C
+++ b/src/thermophysicalModels/radiation/radiationModel/radiationModel/newRadiationModel.C
@@ -40,7 +40,7 @@ autoPtr<radiationModel> radiationModel::New
     const volScalarField& T
 )
 {
-    word radiationModelTypeName;
+    word modelType;
 
     // Note: no need to register/keep radiationProperties since models read
     // it themselves.
@@ -59,27 +59,56 @@ autoPtr<radiationModel> radiationModel::New
         );
 
         radiationPropertiesDict.lookup("radiationModel")
-            >> radiationModelTypeName;
+            >> modelType;
     }
 
-    Info<< "Selecting radiationModel " << radiationModelTypeName << endl;
+    Info<< "Selecting radiationModel " << modelType << endl;
+
+    TConstructorTable::iterator cstrIter =
+        TConstructorTablePtr_->find(modelType);
+
+    if (cstrIter == TConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "radiationModel::New(const volScalarField&)"
+        )   << "Unknown radiationModel type " << modelType
+            << nl << nl
+            << "Valid radiationModel types are:" << nl
+            << TConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    return autoPtr<radiationModel>(cstrIter()(T));
+}
+
+
+autoPtr<radiationModel> radiationModel::New
+(
+    const dictionary& dict,
+    const volScalarField& T
+)
+{
+    const word modelType(dict.lookup("radiationModel"));
+
+    Info<< "Selecting radiationModel " << modelType << endl;
 
     dictionaryConstructorTable::iterator cstrIter =
-        dictionaryConstructorTablePtr_->find(radiationModelTypeName);
+        dictionaryConstructorTablePtr_->find(modelType);
 
     if (cstrIter == dictionaryConstructorTablePtr_->end())
     {
         FatalErrorIn
         (
-            "radiationModel::New(const volScalarField&)"
-        )   << "Unknown radiationModel type " << radiationModelTypeName
-            << nl << nl
+            "radiationModel::New(const dictionary&, const volScalarField&)"
+        )   << "Unknown radiationModel type "
+            << modelType << nl << nl
             << "Valid radiationModel types are:" << nl
-            << dictionaryConstructorTablePtr_->sortedToc()
+            << dictionaryConstructorTablePtr_->toc()
             << exit(FatalError);
     }
 
-    return autoPtr<radiationModel>(cstrIter()(T));
+    return autoPtr<radiationModel>(cstrIter()(dict, T));
 }
 
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
index 3ddd47a72..8c9779f1c 100644
--- a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
+++ b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
@@ -35,11 +35,28 @@ namespace Foam
     namespace radiation
     {
         defineTypeNameAndDebug(radiationModel, 0);
+        defineRunTimeSelectionTable(radiationModel, T);
         defineRunTimeSelectionTable(radiationModel, dictionary);
     }
 }
 
 
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void Foam::radiation::radiationModel::initialise()
+{
+    if (radiation_)
+    {
+        solverFreq_ = max(1, lookupOrDefault<label>("solverFreq", 1));
+
+        absorptionEmission_.reset
+        (
+            absorptionEmissionModel::New(*this, mesh_).ptr()
+        );
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::radiation::radiationModel::radiationModel(const volScalarField& T)
@@ -89,10 +106,44 @@ Foam::radiation::radiationModel::radiationModel
     radiation_(lookup("radiation")),
     coeffs_(subDict(type + "Coeffs")),
     solverFreq_(readLabel(lookup("solverFreq"))),
-    absorptionEmission_(absorptionEmissionModel::New(*this, mesh_)),
+    absorptionEmission_(NULL),
     scatter_(scatterModel::New(*this, mesh_))
 {
-    solverFreq_ = max(1, solverFreq_);
+    initialise();
+}
+
+
+Foam::radiation::radiationModel::radiationModel
+(
+    const word& type,
+    const dictionary& dict,
+    const volScalarField& T
+)
+:
+    IOdictionary
+    (
+        IOobject
+        (
+            "radiationProperties",
+            T.time().constant(),
+            T.mesh(),
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        dict
+    ),
+    mesh_(T.mesh()),
+    time_(T.time()),
+    T_(T),
+    radiation_(lookup("radiation")),
+    coeffs_(subOrEmptyDict(type + "Coeffs")),
+    solverFreq_(1),
+    firstIter_(true),
+    absorptionEmission_(NULL),
+    scatter_(NULL)
+
+{
+    initialise();
 }
 
 
@@ -109,7 +160,10 @@ bool Foam::radiation::radiationModel::read()
     if (regIOobject::read())
     {
         lookup("radiation") >> radiation_;
-        coeffs_ = subDict(type() + "Coeffs");
+        coeffs_ = subOrEmptyDict(type() + "Coeffs");
+
+        solverFreq_ = lookupOrDefault<label>("solverFreq", 1);
+        solverFreq_ = max(1, solverFreq_);
 
         return true;
     }
@@ -127,9 +181,10 @@ void Foam::radiation::radiationModel::correct()
         return;
     }
 
-    if (time_.timeIndex() % solverFreq_ == 0)
+    if (firstIter_ || (time_.timeIndex() % solverFreq_ == 0))
     {
         calculate();
+        firstIter_ = false;
     }
 }
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.H b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.H
index 1052ccdbc..10ead0c03 100644
--- a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.H
+++ b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.H
@@ -35,6 +35,7 @@ Description
 
 SourceFiles
     radiationModel.C
+    radiationModelNew.C
 
 \*---------------------------------------------------------------------------*/
 
@@ -81,7 +82,7 @@ protected:
         //- Reference to the temperature field
         const volScalarField& T_;
 
-        //- Model specific dictionary input parameters
+        //- Radiation model on/off flag
         Switch radiation_;
 
         //- Radiation model dictionary
@@ -91,6 +92,9 @@ protected:
         //  radiation solver iteration
         label solverFreq_;
 
+        //- Flag to enable radiation model to be evaluated on first iteration
+        bool firstIter_;
+
 
         // References to the radiation sub-models
 
@@ -105,6 +109,12 @@ private:
 
     // Private Member Functions
 
+        //- Create IO object if dictionary is present
+        IOobject createIOobject(const fvMesh& mesh) const;
+
+        //- Initialise
+        void initialise();
+
         //- Disallow default bitwise copy construct
         radiationModel(const radiationModel&);
 
@@ -125,12 +135,24 @@ public:
          (
              autoPtr,
              radiationModel,
-             dictionary,
+             T,
              (
                  const volScalarField& T
              ),
              (T)
          );
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            radiationModel,
+            dictionary,
+            (
+                const dictionary& dict,
+                const volScalarField& T
+            ),
+            (dict, T)
+        );
 #endif
 
     // Constructors
@@ -141,11 +163,26 @@ public:
         //- Construct from components
         radiationModel(const word& type, const volScalarField& T);
 
+        //- Construct from components
+        radiationModel
+        (
+            const word& type,
+            const dictionary& dict,
+            const volScalarField& T
+        );
+
 
     // Selectors
 
-         //- Return a reference to the selected radiation model
-         static autoPtr<radiationModel> New(const volScalarField& T);
+        //- Return a reference to the selected radiation model
+        static autoPtr<radiationModel> New(const volScalarField& T);
+
+        //- Return a reference to the selected radiation model
+        static autoPtr<radiationModel> New
+        (
+            const dictionary& dict,
+            const volScalarField& T
+        );
 
 
     //- Destructor
@@ -168,6 +205,12 @@ public:
 
         // Access
 
+            //- Radiation model on/off flag
+            const Switch radiation() const
+            {
+                return radiation_;
+            }
+
             //- Const access to mesh
             const fvMesh& mesh() const
             {
@@ -188,6 +231,24 @@ public:
 };
 
 
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#define addToRadiationRunTimeSelectionTables(model)                           \
+                                                                              \
+    addToRunTimeSelectionTable                                                \
+    (                                                                         \
+        radiationModel,                                                       \
+        model,                                                                \
+        dictionary                                                            \
+    );                                                                        \
+                                                                              \
+    addToRunTimeSelectionTable                                                \
+    (                                                                         \
+        radiationModel,                                                       \
+        model,                                                                \
+        T                                                                     \
+    );
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace radiation