diff --git a/src/lagrangian/coalCombustion/Make/options b/src/lagrangian/coalCombustion/Make/options
index 04b133ba05..87e1d6a710 100644
--- a/src/lagrangian/coalCombustion/Make/options
+++ b/src/lagrangian/coalCombustion/Make/options
@@ -18,6 +18,7 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
@@ -43,6 +44,7 @@ LIB_LIBS = \
     -lincompressibleTransportModels \
     -lregionModels \
     -lsurfaceFilmModels \
+    -ldynamicMesh \
     -ldynamicFvMesh \
     -lfiniteVolume \
     -lmeshTools
diff --git a/src/lagrangian/intermediate/Make/files b/src/lagrangian/intermediate/Make/files
index 9cadf55117..3abb93961e 100644
--- a/src/lagrangian/intermediate/Make/files
+++ b/src/lagrangian/intermediate/Make/files
@@ -1,3 +1,5 @@
+submodels/CloudFunctionObjects/CloudToVTK/vtkTools.C
+
 PARCELS=parcels
 BASEPARCELS=$(PARCELS)/baseClasses
 DERIVEDPARCELS=$(PARCELS)/derived
diff --git a/src/lagrangian/intermediate/Make/options b/src/lagrangian/intermediate/Make/options
index a5a3b3c8eb..dc893c1148 100644
--- a/src/lagrangian/intermediate/Make/options
+++ b/src/lagrangian/intermediate/Make/options
@@ -14,6 +14,7 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
@@ -35,6 +36,7 @@ LIB_LIBS = \
     -lincompressibleTransportModels \
     -lregionModels \
     -lsurfaceFilmModels \
+    -ldynamicMesh \
     -ldynamicFvMesh \
     -lsampling \
     -lfiniteVolume \
diff --git a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.C b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.C
index 8a0001f842..975ab88366 100644
--- a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.C
+++ b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -28,11 +28,7 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::cloudSolution::cloudSolution
-(
-    const fvMesh& mesh,
-    const dictionary& dict
-)
+Foam::cloudSolution::cloudSolution(const fvMesh& mesh, const dictionary& dict)
 :
     mesh_(mesh),
     dict_(dict),
@@ -42,6 +38,7 @@ Foam::cloudSolution::cloudSolution
     maxCo_(0.3),
     iter_(1),
     trackTime_(0.0),
+    deltaTMax_(GREAT),
     coupled_(false),
     cellValueSourceCorrection_(false),
     maxTrackTime_(0.0),
@@ -52,13 +49,25 @@ Foam::cloudSolution::cloudSolution
     {
         read();
     }
+    else
+    {
+        // see if existing source terms should be reset
+        const dictionary sourceTerms(dict_.subOrEmptyDict("sourceTerms"));
+        sourceTerms.readIfPresent("resetOnStartup", resetSourcesOnStartup_);
+
+        if (resetSourcesOnStartup_)
+        {
+            Info<< "Cloud source terms will be reset" << endl;
+        }
+        else
+        {
+            Info<< "Cloud source terms will be held constant" << endl;
+        }
+    }
 }
 
 
-Foam::cloudSolution::cloudSolution
-(
-    const cloudSolution& cs
-)
+Foam::cloudSolution::cloudSolution(const cloudSolution& cs)
 :
     mesh_(cs.mesh_),
     dict_(cs.dict_),
@@ -68,6 +77,7 @@ Foam::cloudSolution::cloudSolution
     maxCo_(cs.maxCo_),
     iter_(cs.iter_),
     trackTime_(cs.trackTime_),
+    deltaTMax_(cs.deltaTMax_),
     coupled_(cs.coupled_),
     cellValueSourceCorrection_(cs.cellValueSourceCorrection_),
     maxTrackTime_(cs.maxTrackTime_),
@@ -76,10 +86,7 @@ Foam::cloudSolution::cloudSolution
 {}
 
 
-Foam::cloudSolution::cloudSolution
-(
-    const fvMesh& mesh
-)
+Foam::cloudSolution::cloudSolution(const fvMesh& mesh)
 :
     mesh_(mesh),
     dict_(dictionary::null),
@@ -89,6 +96,7 @@ Foam::cloudSolution::cloudSolution
     maxCo_(GREAT),
     iter_(0),
     trackTime_(0.0),
+    deltaTMax_(GREAT),
     coupled_(false),
     cellValueSourceCorrection_(false),
     maxTrackTime_(0.0),
@@ -111,6 +119,7 @@ void Foam::cloudSolution::read()
     dict_.lookup("coupled") >> coupled_;
     dict_.lookup("cellValueSourceCorrection") >> cellValueSourceCorrection_;
     dict_.readIfPresent("maxCo", maxCo_);
+    dict_.readIfPresent("deltaTMax", deltaTMax_);
 
     if (steadyState())
     {
@@ -229,4 +238,23 @@ bool Foam::cloudSolution::output() const
 }
 
 
+Foam::scalar Foam::cloudSolution::deltaTMax(const scalar trackTime) const
+{
+    if (transient_)
+    {
+        return min(deltaTMax_, maxCo_*trackTime);
+    }
+    else
+    {
+        return min(deltaTMax_, trackTime);
+    }
+}
+
+
+Foam::scalar Foam::cloudSolution::deltaLMax(const scalar lRef) const
+{
+    return maxCo_*lRef;
+}
+
+
 // ************************************************************************* //
diff --git a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.H b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.H
index e155abc375..01e5c88025 100644
--- a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.H
+++ b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolution.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -80,6 +80,9 @@ class cloudSolution
         //- Particle track time
         scalar trackTime_;
 
+        //- Maximum integration time step (optional)
+        scalar deltaTMax_;
+
 
         // Run-time options
 
@@ -171,6 +174,9 @@ public:
             //- Return the particle track time
             inline scalar trackTime() const;
 
+            //- Return the maximum integation time step
+            inline scalar deltaTMax() const;
+
             //- Return const access to the coupled flag
             inline const Switch coupled() const;
 
@@ -207,6 +213,12 @@ public:
 
         //- Returns true if writing this step
         bool output() const;
+
+        //- Return the maximum integration time
+        scalar deltaTMax(const scalar trackTime) const;
+
+        //- Return the maximum integration length
+        scalar deltaLMax(const scalar lRef) const;
 };
 
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolutionI.H b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolutionI.H
index 6280b97837..27b0b38fde 100644
--- a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolutionI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/cloudSolution/cloudSolutionI.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -101,6 +101,12 @@ inline Foam::scalar Foam::cloudSolution::trackTime() const
 }
 
 
+inline Foam::scalar Foam::cloudSolution::deltaTMax() const
+{
+    return deltaTMax_;
+}
+
+
 inline Foam::Switch& Foam::cloudSolution::coupled()
 {
     return coupled_;
diff --git a/src/lagrangian/intermediate/clouds/Templates/MPPICCloud/MPPICCloud.C b/src/lagrangian/intermediate/clouds/Templates/MPPICCloud/MPPICCloud.C
index 8c7c2cbb55..49d60b1c6a 100644
--- a/src/lagrangian/intermediate/clouds/Templates/MPPICCloud/MPPICCloud.C
+++ b/src/lagrangian/intermediate/clouds/Templates/MPPICCloud/MPPICCloud.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013-2014 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -80,25 +80,25 @@ Foam::MPPICCloud<CloudType>::MPPICCloud
     dampingModel_(NULL),
     isotropyModel_(NULL)
 {
-    if (this->solution().steadyState())
-    {
-        FatalErrorIn
-        (
-            "Foam::MPPICCloud<CloudType>::MPPICCloud"
-            "("
-                "const word&, "
-                "const volScalarField&, "
-                "const volVectorField&, "
-                "const volScalarField&, "
-                "const dimensionedVector&, "
-                "bool"
-            ")"
-        )   << "MPPIC modelling not available for steady state calculations"
-            << exit(FatalError);
-    }
-
     if (this->solution().active())
     {
+        if (this->solution().steadyState())
+        {
+            FatalErrorIn
+            (
+                "Foam::MPPICCloud<CloudType>::MPPICCloud"
+                "("
+                    "const word&, "
+                    "const volScalarField&, "
+                    "const volVectorField&, "
+                    "const volScalarField&, "
+                    "const dimensionedVector&, "
+                    "bool"
+                ")"
+            )   << "MPPIC modelling not available for steady state calculations"
+                << exit(FatalError);
+        }
+
         setModels();
 
         if (readFields)
diff --git a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.C b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.C
index 69666fedde..24c9cc1b9d 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.C
+++ b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -51,10 +51,7 @@ void Foam::ThermoCloud<CloudType>::setModels()
         ).ptr()
     );
 
-    if (this->solution().coupled())
-    {
-        this->subModelProperties().lookup("radiation") >> radiation_;
-    }
+    this->subModelProperties().lookup("radiation") >> radiation_;
 
     if (radiation_)
     {
diff --git a/src/lagrangian/intermediate/parcels/Templates/CollidingParcel/CollidingParcel.C b/src/lagrangian/intermediate/parcels/Templates/CollidingParcel/CollidingParcel.C
index 36f146d77a..b765d4d285 100644
--- a/src/lagrangian/intermediate/parcels/Templates/CollidingParcel/CollidingParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/CollidingParcel/CollidingParcel.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -69,6 +69,8 @@ bool Foam::CollidingParcel<ParcelType>::move
     typename TrackData::cloudType::parcelType& p =
         static_cast<typename TrackData::cloudType::parcelType&>(*this);
 
+    td.keepParticle = true;
+
     switch (td.part())
     {
         case TrackData::tpVelocityHalfStep:
diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
index f34f4b2f28..6a96921be3 100644
--- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,6 +27,7 @@ License
 #include "forceSuSp.H"
 #include "IntegrationScheme.H"
 #include "meshTools.H"
+#include "cloudSolution.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -271,15 +272,11 @@ bool Foam::KinematicParcel<ParcelType>::move
 
     const polyMesh& mesh = td.cloud().pMesh();
     const polyBoundaryMesh& pbMesh = mesh.boundaryMesh();
+    const cloudSolution& solution = td.cloud().solution();
     const scalarField& cellLengthScale = td.cloud().cellLengthScale();
-    const scalar maxCo = td.cloud().solution().maxCo();
 
     scalar tEnd = (1.0 - p.stepFraction())*trackTime;
-    scalar dtMax = trackTime;
-    if (td.cloud().solution().transient())
-    {
-        dtMax *= maxCo;
-    }
+    scalar dtMax = solution.deltaTMax(trackTime);
 
     bool tracking = true;
     label nTrackingStalled = 0;
@@ -301,7 +298,8 @@ bool Foam::KinematicParcel<ParcelType>::move
         if (p.active() && tracking && (magU > ROOTVSMALL))
         {
             const scalar d = dt*magU;
-            const scalar dCorr = min(d, maxCo*cellLengthScale[cellI]);
+            const scalar deltaLMax = solution.deltaLMax(cellLengthScale[cellI]);
+            const scalar dCorr = min(d, deltaLMax);
             dt *=
                 dCorr/d
                *p.trackToFace(p.position() + dCorr*U_/magU, td);
@@ -309,7 +307,7 @@ bool Foam::KinematicParcel<ParcelType>::move
 
         tEnd -= dt;
 
-        scalar newStepFraction = 1.0 - tEnd/trackTime;
+        const scalar newStepFraction = 1.0 - tEnd/trackTime;
 
         if (tracking)
         {
@@ -335,7 +333,7 @@ bool Foam::KinematicParcel<ParcelType>::move
         p.stepFraction() = newStepFraction;
 
         bool calcParcel = true;
-        if (!tracking && td.cloud().solution().steadyState())
+        if (!tracking && solution.steadyState())
         {
             calcParcel = false;
         }
@@ -346,7 +344,7 @@ bool Foam::KinematicParcel<ParcelType>::move
             // Update cell based properties
             p.setCellValues(td, dt, cellI);
 
-            if (td.cloud().solution().cellValueSourceCorrection())
+            if (solution.cellValueSourceCorrection())
             {
                 p.cellValueSourceCorrection(td, dt, cellI);
             }
@@ -466,6 +464,8 @@ void Foam::KinematicParcel<ParcelType>::hitPatch
 )
 {
     td.keepParticle = false;
+
+    td.cloud().patchInteraction().addToEscapedParcels(nParticle_*mass());
 }
 
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
index 289127b18b..2f050e26ae 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
@@ -121,7 +121,13 @@ Foam::scalar Foam::ReactingMultiphaseParcel<ParcelType>::updateMassFractions
 
     YMix[GAS] = massGas/massNew;
     YMix[LIQ] = massLiquid/massNew;
-    YMix[SLD] = 1.0 - YMix[GAS] - YMix[LIQ];
+    YMix[SLD] = massSolid/massNew;
+
+    scalar Ytotal = sum(YMix);
+
+    YMix[GAS] /= Ytotal;
+    YMix[LIQ] /= Ytotal;
+    YMix[SLD] /= Ytotal;
 
     return massNew;
 }
diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
index 4198420d38..ae5f386411 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
@@ -80,6 +80,9 @@ void Foam::ThermoParcel<ParcelType>::cellValueSourceCorrection
     this->Uc_ += td.cloud().UTrans()[cellI]/this->massCell(cellI);
 
     const scalar CpMean = td.CpInterp().psi()[cellI];
+
+    tetIndices tetIs = this->currentTetIndices();
+    Tc_ = td.TInterp().interpolate(this->position(), tetIs);
     Tc_ += td.cloud().hsTrans()[cellI]/(CpMean*this->massCell(cellI));
 
     if (Tc_ < td.cloud().constProps().TMin())
diff --git a/src/lagrangian/intermediate/parcels/include/makeParcelCloudFunctionObjects.H b/src/lagrangian/intermediate/parcels/include/makeParcelCloudFunctionObjects.H
index 783878f480..a181ae19b1 100644
--- a/src/lagrangian/intermediate/parcels/include/makeParcelCloudFunctionObjects.H
+++ b/src/lagrangian/intermediate/parcels/include/makeParcelCloudFunctionObjects.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -28,6 +28,7 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+#include "CloudToVTK.H"
 #include "FacePostProcessing.H"
 #include "ParticleCollector.H"
 #include "ParticleErosion.H"
@@ -42,6 +43,7 @@ License
                                                                               \
     makeCloudFunctionObject(CloudType);                                       \
                                                                               \
+    makeCloudFunctionObjectType(CloudToVTK, CloudType);                       \
     makeCloudFunctionObjectType(FacePostProcessing, CloudType);               \
     makeCloudFunctionObjectType(ParticleCollector, CloudType);                \
     makeCloudFunctionObjectType(ParticleErosion, CloudType);                  \
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.C
new file mode 100644
index 0000000000..240154edf3
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.C
@@ -0,0 +1,164 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CloudToVTK.H"
+#include "vtkTools.H"
+#include "floatScalar.H"
+
+// * * * * * * * * * * * * * Protectd Member Functions * * * * * * * * * * * //
+
+template<class CloudType>
+void Foam::CloudToVTK<CloudType>::writeData
+(
+    std::ostream& vtkOs,
+    const bool binary,
+    const List<floatScalar>& data
+) const
+{
+    const label procI = Pstream::myProcNo();
+
+    List<List<floatScalar> > allProcData(Pstream::nProcs());
+    allProcData[procI] = data;
+    Pstream::gatherList(allProcData);
+    List<floatScalar> allData =
+        ListListOps::combine<List<floatScalar> >
+        (
+            allProcData,
+            accessOp<List<floatScalar> >()
+        );
+
+    vtkTools::write(vtkOs, binary, allData);
+}
+
+
+template<class CloudType>
+template<class Type>
+void Foam::CloudToVTK<CloudType>::writeFieldData
+(
+    std::ostream& vtkOs,
+    const bool binary,
+    const List<floatScalar>& data,
+    const word& title,
+    const label nParcels
+) const
+{
+    vtkOs
+        << title << ' ' << pTraits<Type>::nComponents << ' '
+        << nParcels << " float" << std::endl;
+    writeData(vtkOs, binary, data);
+}
+
+
+template<class CloudType>
+void Foam::CloudToVTK<CloudType>::write()
+{
+    label nParcels = this->owner().size();
+    DynamicList<floatScalar> position(3*nParcels);
+    DynamicList<floatScalar> U(3*nParcels);
+    DynamicList<floatScalar> d(nParcels);
+    DynamicList<floatScalar> age(nParcels);
+    DynamicList<floatScalar> rho(nParcels);
+
+    forAllConstIter(typename CloudType, this->owner(), iter)
+    {
+        vtkTools::insert(iter().position(), position);
+        vtkTools::insert(iter().U(), U);
+        vtkTools::insert(iter().d(), d);
+        vtkTools::insert(iter().age(), age);
+        vtkTools::insert(iter().rho(), rho);
+    }
+
+    reduce(nParcels, sumOp<label>());
+
+
+binary_ = false;
+    if (Pstream::master())
+    {
+        // Create directory if does not exist
+        mkDir(this->outputTimeDir());
+
+        // Open new file at start up
+
+        const fileName fName = this->outputTimeDir()/(type() + ".vtk");
+        this->setModelProperty("file", fName);
+
+        OFstream os(fName, binary_ ? IOstream::BINARY : IOstream::ASCII);
+        std::ostream& vtkOs = os.stdStream();
+
+
+        vtkTools::writeHeader(vtkOs, binary_, this->modelName().c_str());
+        vtkOs
+            << "DATASET POLYDATA" << std::endl
+            << "POINTS " << nParcels << " float" << std::endl;
+
+        writeData(vtkOs, binary_, position);
+
+        vtkOs
+            << "POINT_DATA " << nParcels << std::endl
+            << "FIELD attributes " << 4
+            << std::endl;
+
+        writeFieldData<vector>(vtkOs, binary_, U, "U", nParcels);
+        writeFieldData<scalar>(vtkOs, binary_, d, "d", nParcels);
+        writeFieldData<scalar>(vtkOs, binary_, age, "age", nParcels);
+        writeFieldData<scalar>(vtkOs, binary_, rho, "rho", nParcels);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::CloudToVTK<CloudType>::CloudToVTK
+(
+    const dictionary& dict,
+    CloudType& owner,
+    const word& modelName
+)
+:
+    CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
+    binary_(dict.lookupOrDefault<bool>("binary", true))
+{}
+
+
+template<class CloudType>
+Foam::CloudToVTK<CloudType>::CloudToVTK
+(
+    const CloudToVTK<CloudType>& c
+)
+:
+    CloudFunctionObject<CloudType>(c),
+    binary_(c.binary_)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::CloudToVTK<CloudType>::~CloudToVTK()
+{}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.H b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.H
new file mode 100644
index 0000000000..13a68f1879
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/CloudToVTK.H
@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::CloudToVTK
+
+Description
+    Write cloud data in VTK format
+
+SourceFiles
+    CloudToVTK.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef CloudToVTK_H
+#define CloudToVTK_H
+
+#include "CloudFunctionObject.H"
+#include "volFields.H"
+#include "OFstream.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class CloudToVTK Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CloudType>
+class CloudToVTK
+:
+    public CloudFunctionObject<CloudType>
+{
+
+protected:
+
+    // Protected data
+
+        //- Ascii/binary output flag
+        bool binary_;
+
+
+    // Protected Member Functions
+
+        //- Write post-processing info
+        virtual void write();
+
+        //- Helper function to write VTK data
+        void writeData
+        (
+            std::ostream& vtkOs,
+            const bool binary,
+            const List<floatScalar>& data
+        ) const;
+
+        //- Helper function to write VTK field data
+        template<class Type>
+        void writeFieldData
+        (
+            std::ostream& vtkOs,
+            const bool binary,
+            const List<floatScalar>& data,
+            const word& title,
+            const label nParcels
+        ) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("cloudToVTK");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        CloudToVTK
+        (
+            const dictionary& dict,
+            CloudType& owner,
+            const word& modelName
+        );
+
+        //- Construct copy
+        CloudToVTK(const CloudToVTK<CloudType>& c);
+
+        //- Construct and return a clone
+        virtual autoPtr<CloudFunctionObject<CloudType> > clone() const
+        {
+            return autoPtr<CloudFunctionObject<CloudType> >
+            (
+                new CloudToVTK<CloudType>(*this)
+            );
+        }
+
+
+    //- Destructor
+    virtual ~CloudToVTK();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "CloudToVTK.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.C
new file mode 100644
index 0000000000..71de625526
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.C
@@ -0,0 +1,287 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "vtkTools.H"
+
+#if defined(__mips)
+#include <standards.h>
+#include <sys/endian.h>
+#endif
+
+// MacOSX
+#ifdef __DARWIN_BYTE_ORDER
+#if __DARWIN_BYTE_ORDER==__DARWIN_BIG_ENDIAN
+#undef LITTLE_ENDIAN
+#else
+#undef BIG_ENDIAN
+#endif
+#endif
+
+#if defined(LITTLE_ENDIAN) \
+ || defined(_LITTLE_ENDIAN) \
+ || defined(__LITTLE_ENDIAN)
+#   define LITTLEENDIAN 1
+#elif defined(BIG_ENDIAN) || defined(_BIG_ENDIAN) || defined(__BIG_ENDIAN)
+#   undef LITTLEENDIAN
+#else
+#   error "Cannot find LITTLE_ENDIAN or BIG_ENDIAN symbol defined."
+#   error "Please add to compilation options"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+void Foam::vtkTools::swapWord(label& word32)
+{
+    char* mem = reinterpret_cast<char*>(&word32);
+
+    char a = mem[0];
+    mem[0] = mem[3];
+    mem[3] = a;
+
+    a = mem[1];
+    mem[1] = mem[2];
+    mem[2] = a;
+}
+
+
+void Foam::vtkTools::swapWords(const label nWords, label* words32)
+{
+    for (label i = 0; i < nWords; i++)
+    {
+        swapWord(words32[i]);
+    }
+}
+
+
+void Foam::vtkTools::write
+(
+    std::ostream& os,
+    const bool binary,
+    List<floatScalar>& fField
+)
+{
+    if (binary)
+    {
+#       ifdef LITTLEENDIAN
+        swapWords(fField.size(), reinterpret_cast<label*>(fField.begin()));
+#       endif
+        os.write
+        (
+            reinterpret_cast<char*>(fField.begin()),
+            fField.size()*sizeof(float)
+        );
+
+        os  << std::endl;
+    }
+    else
+    {
+        forAll(fField, i)
+        {
+            os  << fField[i];
+
+            if (i > 0 && (i % 10) == 0)
+            {
+                os  << std::endl;
+            }
+            else
+            {
+                os  << ' ';
+            }
+        }
+        os  << std::endl;
+    }
+}
+
+
+void Foam::vtkTools::write
+(
+    std::ostream& os,
+    const bool binary,
+    DynamicList<floatScalar>& fField
+)
+{
+    List<floatScalar>& fld = fField.shrink();
+
+    write(os, binary, fld);
+}
+
+
+void Foam::vtkTools::write
+(
+    std::ostream& os,
+    const bool binary,
+    labelList& elems
+)
+{
+    if (binary)
+    {
+#       ifdef LITTLEENDIAN
+        swapWords(elems.size(), reinterpret_cast<label*>(elems.begin()));
+#       endif
+        os.write
+        (
+            reinterpret_cast<char*>(elems.begin()),
+            elems.size()*sizeof(label)
+        );
+
+        os  << std::endl;
+    }
+    else
+    {
+        forAll(elems, i)
+        {
+            os  << elems[i];
+
+            if (i > 0 && (i % 10) == 0)
+            {
+                os  << std::endl;
+            }
+            else
+            {
+                os  << ' ';
+            }
+        }
+        os  << std::endl;
+    }
+}
+
+
+void Foam::vtkTools::write
+(
+    std::ostream& os,
+    const bool binary,
+    DynamicList<label>& elems
+)
+{
+    labelList& fld = elems.shrink();
+
+    write(os, binary, fld);
+}
+
+
+void Foam::vtkTools::writeHeader
+(
+    std::ostream& os,
+    const bool binary,
+    const std::string& title
+)
+{
+    os  << "# vtk DataFile Version 2.0" << std::endl
+        << title << std::endl;
+
+    if (binary)
+    {
+        os  << "BINARY" << std::endl;
+    }
+    else
+    {
+        os  << "ASCII" << std::endl;
+    }
+}
+
+
+void Foam::vtkTools::writeCellDataHeader
+(
+    std::ostream& os,
+    const label nCells,
+    const label nFields
+)
+{
+    os  << "CELL_DATA " << nCells << std::endl
+        << "FIELD attributes " << nFields << std::endl;
+}
+
+
+void Foam::vtkTools::writePointDataHeader
+(
+    std::ostream& os,
+    const label nPoints,
+    const label nFields
+)
+{
+    os  << "POINT_DATA  " << nPoints << std::endl
+        << "FIELD attributes " << nFields << std::endl;
+}
+
+
+void Foam::vtkTools::insert(const scalar src, DynamicList<floatScalar>& dest)
+{
+    dest.append(float(src));
+}
+
+
+void Foam::vtkTools::insert(const vector& src, DynamicList<floatScalar>& dest)
+{
+    for (direction cmpt = 0; cmpt < vector::nComponents; ++cmpt)
+    {
+        dest.append(float(src[cmpt]));
+    }
+}
+
+
+void Foam::vtkTools::insert
+(
+    const sphericalTensor& src,
+    DynamicList<floatScalar>& dest
+)
+{
+    for (direction cmpt = 0; cmpt < sphericalTensor::nComponents; ++cmpt)
+    {
+        dest.append(float(src[cmpt]));
+    }
+}
+
+
+void Foam::vtkTools::insert
+(
+    const symmTensor& src,
+    DynamicList<floatScalar>& dest
+)
+{
+    dest.append(float(src.xx()));
+    dest.append(float(src.yy()));
+    dest.append(float(src.zz()));
+    dest.append(float(src.xy()));
+    dest.append(float(src.yz()));
+    dest.append(float(src.xz()));
+}
+
+
+void Foam::vtkTools::insert(const tensor& src, DynamicList<floatScalar>& dest)
+{
+    for (direction cmpt = 0; cmpt < tensor::nComponents; ++cmpt)
+    {
+        dest.append(float(src[cmpt]));
+    }
+}
+
+
+void Foam::vtkTools::insert(const labelList& src, DynamicList<label>& dest)
+{
+    dest.append(src);
+}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.H b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.H
new file mode 100644
index 0000000000..12bdda1f38
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkTools.H
@@ -0,0 +1,126 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::vtkTools
+
+Description
+    Various functions for collecting and writing binary data.
+
+SourceFiles
+    vtkTools.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef vtkTools_H
+#define vtkTools_H
+
+#include "floatScalar.H"
+#include "DynamicList.H"
+#include "volFieldsFwd.H"
+#include "pointFieldsFwd.H"
+#include "volPointInterpolation.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class vtkTools Declaration
+\*---------------------------------------------------------------------------*/
+
+class vtkTools
+{
+    // Private Member Functions
+
+        // Swap halves of word.
+
+            static void swapWord(label& word32);
+            static void swapWords(const label nWords, label* words32);
+
+
+public:
+
+    // Write ascii or binary. If binary optionally in-place swaps argument
+
+        static void write(std::ostream&, const bool, List<floatScalar>&);
+        static void write(std::ostream&, const bool, DynamicList<floatScalar>&);
+        static void write(std::ostream&, const bool, labelList&);
+        static void write(std::ostream&, const bool, DynamicList<label>&);
+
+
+    // Write header
+
+        static void writeHeader
+        (
+            std::ostream&,
+            const bool isBinary,
+            const std::string& title
+        );
+        static void writeCellDataHeader
+        (
+            std::ostream&,
+            const label nCells,
+            const label nFields
+        );
+        static void writePointDataHeader
+        (
+            std::ostream&,
+            const label nPoints,
+            const label nFields
+        );
+
+
+    // Convert to VTK and store
+
+        static void insert(const scalar, DynamicList<floatScalar>&);
+        static void insert(const point&, DynamicList<floatScalar>&);
+        static void insert(const sphericalTensor&, DynamicList<floatScalar>&);
+        static void insert(const symmTensor&, DynamicList<floatScalar>&);
+        static void insert(const tensor&, DynamicList<floatScalar>&);
+
+
+    //- Append elements to DynamicList
+    static void insert(const labelList&, DynamicList<label>&);
+
+    template<class Type>
+    static void insert(const List<Type>&, DynamicList<floatScalar>&);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "vtkToolsTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkToolsTemplates.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkToolsTemplates.C
new file mode 100644
index 0000000000..53ea67e398
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/CloudToVTK/vtkToolsTemplates.C
@@ -0,0 +1,44 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "vtkTools.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type>
+void Foam::vtkTools::insert
+(
+    const List<Type>& source,
+    DynamicList<floatScalar>& dest
+)
+{
+    forAll(source, i)
+    {
+        insert(source[i], dest);
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/ParticleCollector/ParticleCollector.C b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/ParticleCollector/ParticleCollector.C
index bf347c9f30..97f66d6bb3 100644
--- a/src/lagrangian/intermediate/submodels/CloudFunctionObjects/ParticleCollector/ParticleCollector.C
+++ b/src/lagrangian/intermediate/submodels/CloudFunctionObjects/ParticleCollector/ParticleCollector.C
@@ -457,7 +457,15 @@ void Foam::ParticleCollector<CloudType>::write()
     {
         if (Pstream::master())
         {
-            autoPtr<surfaceWriter> writer(surfaceWriter::New(surfaceFormat_));
+            autoPtr<surfaceWriter> writer
+            (
+                surfaceWriter::New
+                (
+                    surfaceFormat_,
+                    this->coeffDict().subOrEmptyDict("formatOptions").
+                        subOrEmptyDict(surfaceFormat_)
+                )
+            );
 
             writer->write
             (
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H
index f5e4b31ff9..01336c8518 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H
@@ -41,7 +41,7 @@ Description
 namespace Foam
 {
 /*---------------------------------------------------------------------------*\
-                    Class PairSpringSliderDashpot Declaration
+                   Class PairSpringSliderDashpot Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class CloudType>
@@ -71,7 +71,7 @@ class PairSpringSliderDashpot
         //- Cohesion energy density [J/m^3]
         scalar cohesionEnergyDensity_;
 
-        // Switch cohesion on and off
+        //- Switch cohesion on and off
         bool cohesion_;
 
         //- The number of steps over which to resolve the minimum
@@ -104,7 +104,7 @@ class PairSpringSliderDashpot
     // Private Member Functions
 
         //- Find the appropriate properties for determining the minimum
-        //- Allowable timestep
+        //  allowable timestep
         void findMinMaxProperties
         (
             scalar& RMin,
@@ -136,8 +136,8 @@ public:
             return volumeFactor_;
         }
 
-        // Return the area of overlap between two spheres of radii rA and rB,
-        // centres separated by a distance rAB.  Assumes rAB < (rA + rB).
+        //- Return the area of overlap between two spheres of radii rA and rB,
+        //  centres separated by a distance rAB.  Assumes rAB < (rA + rB).
         inline scalar overlapArea(scalar rA, scalar rB, scalar rAB) const
         {
             // From:
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.H b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.H
index e48e949f57..699c9cda2c 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.H
@@ -68,7 +68,7 @@ class WallLocalSpringSliderDashpot
         //- Cohesion energy density [J/m^3]
         scalarList cohesionEnergyDensity_;
 
-        // Switch cohesion on and off
+        //- Switch cohesion on and off
         boolList cohesion_;
 
         //- Mapping the patch index to the model data
@@ -107,7 +107,7 @@ class WallLocalSpringSliderDashpot
     // Private Member Functions
 
         //- Find the appropriate properties for determining the minimum
-        //- Allowable timestep
+        //  allowable timestep
         void findMinMaxProperties
         (
             scalar& rMin,
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H
index 88dca9b66c..0087d01197 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H
@@ -101,7 +101,7 @@ class WallSpringSliderDashpot
     // Private Member Functions
 
         //- Find the appropriate properties for determining the minimum
-        //- Allowable timestep
+        //  allowable timestep
         void findMinMaxProperties
         (
             scalar& rMin,
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.H
index 9edd88736e..ece119404d 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.H
@@ -27,10 +27,12 @@ Class
 Description
     Injection positions specified by a particle number density within a cell set
 
-    - User specifies
-      - Number density of particles in cell set (effective)
-      - Total mass to inject
-      - Initial parcel velocity
+    User specifies
+    - Number density of particles in cell set (effective)
+    - Total mass to inject
+    - Initial parcel velocity
+
+    Additional
     - Parcel diameters obtained by PDF model
     - All parcels introduced at SOI
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.C
index d71e1b070e..29c0c30812 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -203,8 +203,6 @@ Foam::label Foam::ConeInjection<CloudType>::parcelsToInject
 
         const label nToInject = targetParcels - nInjected_;
 
-        nInjected_ += nToInject;
-
         return positionAxis_.size()*nToInject;
     }
     else
@@ -264,7 +262,7 @@ void Foam::ConeInjection<CloudType>::setProperties
 {
     cachedRandom& rnd = this->owner().rndGen();
 
-    // set particle velocity
+    // Set particle velocity
     const label i = parcelI % positionAxis_.size();
 
     scalar t = time - this->SOI_;
@@ -283,8 +281,11 @@ void Foam::ConeInjection<CloudType>::setProperties
 
     parcel.U() = Umag_.value(t)*dirVec;
 
-    // set particle diameter
+    // Set particle diameter
     parcel.d() = sizeDistribution_().sample();
+
+    // Increment number of particles injected
+    nInjected_++;
 }
 
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.H
index f45fd4a4f0..56519243e4 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeInjection/ConeInjection.H
@@ -26,12 +26,15 @@ Class
 
 Description
     Multi-point cone injection model
-    - User specifies
-      - time of start of injection
-      - list of injector positions and directions (along injection axes)
-      - number of parcels to inject per injector
-      - parcel velocities
-      - inner and outer half-cone angles
+
+    User specifies
+    - time of start of injection
+    - list of injector positions and directions (along injection axes)
+    - number of parcels to inject per injector
+    - parcel velocities
+    - inner and outer half-cone angles
+
+    Additional
     - Parcel diameters obtained by distribution model
 
 SourceFiles
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C
index 97c3c819a0..4bde7b454f 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -170,8 +170,11 @@ Foam::ConeNozzleInjection<CloudType>::ConeNozzleInjection
                 "CloudType&, "
                 "const word&"
             ")"
-        )<< "innerNozzleDiameter >= outerNozzleDiameter" << nl
-         << exit(FatalError);
+        )
+            << "Inner diameter must be less than the outer diameter:" << nl
+            << "    innerDiameter: " << innerDiameter_ << nl
+            << "    outerDiameter: " << outerDiameter_
+            << exit(FatalError);
     }
 
     duration_ = owner.db().time().userTimeToTime(duration_);
@@ -265,7 +268,7 @@ void Foam::ConeNozzleInjection<CloudType>::updateMesh()
         }
         default:
         {
-            // do nothing
+            // Do nothing
         }
     }
 }
@@ -344,9 +347,10 @@ void Foam::ConeNozzleInjection<CloudType>::setPositionAndCell
         }
         case imDisc:
         {
-            scalar frac = rndGen.sample01<scalar>();
+            scalar frac = rndGen.globalSample01<scalar>();
             scalar dr = outerDiameter_ - innerDiameter_;
             scalar r = 0.5*(innerDiameter_ + frac*dr);
+
             position = position_ + r*normal_;
 
             this->findCellAtPosition
@@ -354,8 +358,7 @@ void Foam::ConeNozzleInjection<CloudType>::setPositionAndCell
                 cellOwner,
                 tetFaceI,
                 tetPtI,
-                position,
-                false
+                position
             );
             break;
         }
@@ -391,7 +394,7 @@ void Foam::ConeNozzleInjection<CloudType>::setProperties
 {
     cachedRandom& rndGen = this->owner().rndGen();
 
-    // set particle velocity
+    // Set particle velocity
     const scalar deg2Rad = mathematical::pi/180.0;
 
     scalar t = time - this->SOI_;
@@ -441,7 +444,7 @@ void Foam::ConeNozzleInjection<CloudType>::setProperties
         }
     }
 
-    // set particle diameter
+    // Set particle diameter
     parcel.d() = sizeDistribution_->sample();
 }
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.H
index bd67443ce1..15713700ee 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.H
@@ -27,22 +27,21 @@ Class
 Description
     Cone injection
 
-    - User specifies
-      - time of start of injection
-      - injector position
-      - direction (along injection axis)
-      - parcel flow rate
-      - inner and outer half-cone angles
+    User specifies
+    - time of start of injection
+    - injector position
+    - direction (along injection axis)
+    - parcel flow rate
+    - inner and outer half-cone angles
 
+    Additional
     - Parcel diameters obtained by size distribution model
-
     - Parcel velocity is calculated as:
-
-        - Constant velocity
-            U = <specified by user>
-        - Pressure driven velocity
-            U = sqrt(2*(Pinj - Pamb)/rho)
-        - Flow rate and discharge
+       - Constant velocity
+           U = \<specified by user\>
+       - Pressure driven velocity
+           U = sqrt(2*(Pinj - Pamb)/rho)
+       - Flow rate and discharge
             U = V_dot/(A*Cd)
 
 SourceFiles
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/FieldActivatedInjection/FieldActivatedInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/FieldActivatedInjection/FieldActivatedInjection.H
index 704345d8f5..8e7c7d682d 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/FieldActivatedInjection/FieldActivatedInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/FieldActivatedInjection/FieldActivatedInjection.H
@@ -25,18 +25,19 @@ Class
     Foam::FieldActivatedInjection
 
 Description
-    Injection at specified positions, with the conditions:
+    Conditional injection at specified positions.
 
-    - for injection to be allowed
+    For injection to be allowed, the equality must be satisfied:
 
         factor*referenceField[cellI] >= thresholdField[cellI]
 
-        where:
-          - referenceField is the field used to supply the look-up values
-          - thresholdField supplies the values beyond which the injection is
+        Where:
+        - referenceField is the field used to supply the look-up values
+        - thresholdField supplies the values beyond which the injection is
             permitted
 
-    - limited to a user-supllied number of injections per injector location
+    Note:
+    - Limited to a user-supplied number of injections per injector location
 
 SourceFiles
     FieldActivatedInjection.C
@@ -56,7 +57,7 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                    Class FieldActivatedInjection Declaration
+                   Class FieldActivatedInjection Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class CloudType>
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
index 5be9e15ebf..43dc0bc5e5 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
@@ -285,7 +285,8 @@ Foam::InjectionModel<CloudType>::InjectionModel(CloudType& owner)
     nParticleFixed_(0.0),
     time0_(0.0),
     timeStep0_(this->template getModelProperty<scalar>("timeStep0")),
-    delayedVolume_(0.0)
+    delayedVolume_(0.0),
+    injectorID_(-1)
 {}
 
 
@@ -313,7 +314,8 @@ Foam::InjectionModel<CloudType>::InjectionModel
     nParticleFixed_(0.0),
     time0_(owner.db().time().value()),
     timeStep0_(this->template getModelProperty<scalar>("timeStep0")),
-    delayedVolume_(0.0)
+    delayedVolume_(0.0),
+    injectorID_(this->coeffDict().lookupOrDefault("injectorID", -1))
 {
     // Provide some info
     // - also serves to initialise mesh dimensions - needed for parallel runs
@@ -321,18 +323,25 @@ Foam::InjectionModel<CloudType>::InjectionModel
     Info<< "    Constructing " << owner.mesh().nGeometricD() << "-D injection"
         << endl;
 
+    if (injectorID_ != -1)
+    {
+        Info<< "    injector ID: " << injectorID_ << endl;
+    }
+
     if (owner.solution().transient())
     {
         this->coeffDict().lookup("massTotal") >> massTotal_;
         this->coeffDict().lookup("SOI") >> SOI_;
-        SOI_ = owner.db().time().userTimeToTime(SOI_);
     }
     else
     {
         massFlowRate_.reset(this->coeffDict());
         massTotal_ = massFlowRate_.value(owner.db().time().value());
+        this->coeffDict().readIfPresent("SOI", SOI_);
     }
 
+    SOI_ = owner.db().time().userTimeToTime(SOI_);
+
     const word parcelBasisType = this->coeffDict().lookup("parcelBasisType");
 
     if (parcelBasisType == "mass")
@@ -387,7 +396,8 @@ Foam::InjectionModel<CloudType>::InjectionModel
     nParticleFixed_(im.nParticleFixed_),
     time0_(im.time0_),
     timeStep0_(im.timeStep0_),
-    delayedVolume_(im.delayedVolume_)
+    delayedVolume_(im.delayedVolume_),
+    injectorID_(im.injectorID_)
 {}
 
 
@@ -440,10 +450,10 @@ void Foam::InjectionModel<CloudType>::inject(TrackData& td)
     scalar massAdded = 0.0;
     label newParcels = 0;
     scalar newVolumeFraction = 0.0;
+    scalar delayedVolume = 0;
 
     if (prepareForNextTimeStep(time, newParcels, newVolumeFraction))
     {
-        scalar delayedVolume = 0;
 
         const scalar trackTime = this->owner().solution().trackTime();
         const polyMesh& mesh = this->owner().mesh();
@@ -529,7 +539,8 @@ void Foam::InjectionModel<CloudType>::inject(TrackData& td)
 
                         if (pPtr->move(td, dt))
                         {
-                            td.cloud().addParticle(pPtr);
+                            pPtr->typeId() = injectorID_;
+                            cloud.addParticle(pPtr);
                         }
                         else
                         {
@@ -544,10 +555,10 @@ void Foam::InjectionModel<CloudType>::inject(TrackData& td)
                 }
             }
         }
-
-        delayedVolume_ = delayedVolume;
     }
 
+    delayedVolume_ = returnReduce(delayedVolume, sumOp<scalar>());
+
     postInjectCheck(parcelsAdded, massAdded);
 }
 
@@ -565,6 +576,13 @@ void Foam::InjectionModel<CloudType>::injectSteadyState
         return;
     }
 
+    const scalar time = this->owner().db().time().value();
+
+    if (time < SOI_)
+    {
+        return;
+    }
+
     const polyMesh& mesh = this->owner().mesh();
     typename TrackData::cloudType& cloud = td.cloud();
 
@@ -634,8 +652,10 @@ void Foam::InjectionModel<CloudType>::injectSteadyState
                     pPtr->rho()
                 );
 
+            pPtr->typeId() = injectorID_;
+
             // Add the new parcel
-            td.cloud().addParticle(pPtr);
+            cloud.addParticle(pPtr);
 
             massAdded += pPtr->nParticle()*pPtr->mass();
             parcelsAdded++;
@@ -649,9 +669,9 @@ void Foam::InjectionModel<CloudType>::injectSteadyState
 template<class CloudType>
 void Foam::InjectionModel<CloudType>::info(Ostream& os)
 {
-    os  << "    " << this->modelName() << ":" << nl
-        << "        number of parcels added     = " << parcelsAddedTotal_ << nl
-        << "        mass introduced             = " << massInjected_ << nl;
+    os  << "    Injector " << this->modelName() << ":" << nl
+        << "      - parcels added               = " << parcelsAddedTotal_ << nl
+        << "      - mass introduced             = " << massInjected_ << nl;
 
     if (this->outputTime())
     {
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
index b570a6bf4a..fc0e285608 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
@@ -138,6 +138,9 @@ protected:
             //  less than 1 particle per parcel
             scalar delayedVolume_;
 
+            //- Optional injector ID
+            label injectorID_;
+
 
     // Protected Member Functions
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ManualInjection/ManualInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ManualInjection/ManualInjection.H
index bbadab5409..c9a91e8be6 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ManualInjection/ManualInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ManualInjection/ManualInjection.H
@@ -26,10 +26,13 @@ Class
 
 Description
     Manual injection
-    - User specifies
-      - Total mass to inject
-      - Parcel positions in file \c positionsFile
-      - Initial parcel velocity
+
+    User specifies
+    - Total mass to inject
+    - Parcel positions in file \c positionsFile
+    - Initial parcel velocity
+
+    Additional
     - Parcel diameters obtained by distribution model
     - All parcels introduced at SOI
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
index c858a27d01..2483e6b3a1 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
@@ -121,9 +121,8 @@ Foam::label Foam::PatchInjection<CloudType>::parcelsToInject
     if ((time0 >= 0.0) && (time0 < duration_))
     {
         scalar nParcels = (time1 - time0)*parcelsPerSecond_;
-
         cachedRandom& rnd = this->owner().rndGen();
-
+        scalar rndPos = rnd.globalPosition(scalar(0), scalar(1));
         label nParcelsToInject = floor(nParcels);
 
         // Inject an additional parcel with a probability based on the
@@ -131,10 +130,7 @@ Foam::label Foam::PatchInjection<CloudType>::parcelsToInject
         if
         (
             nParcelsToInject > 0
-         && (
-               nParcels - scalar(nParcelsToInject)
-             > rnd.globalPosition(scalar(0), scalar(1))
-            )
+         && (nParcels - scalar(nParcelsToInject) > rndPos)
         )
         {
             ++nParcelsToInject;
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H
index c782ad860d..e32b291851 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H
@@ -26,12 +26,15 @@ Class
 
 Description
     Patch injection
-    - User specifies
-      - Total mass to inject
-      - Name of patch
-      - Injection duration
-      - Initial parcel velocity
-      - Injection volume flow rate
+
+    User specifies
+    - Total mass to inject
+    - Name of patch
+    - Injection duration
+    - Initial parcel velocity
+    - Injection volume flow rate
+
+    Additional
     - Parcel diameters obtained by distribution model
     - Parcels injected randomly across the patch
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Drag/PlessisMasliyahDrag/PlessisMasliyahDragForce.C b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Drag/PlessisMasliyahDrag/PlessisMasliyahDragForce.C
index 970dee90f9..0dd49b7f0a 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Drag/PlessisMasliyahDrag/PlessisMasliyahDragForce.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Drag/PlessisMasliyahDrag/PlessisMasliyahDragForce.C
@@ -104,7 +104,7 @@ Foam::forceSuSp Foam::PlessisMasliyahDragForce<CloudType>::calcCoupled
 {
     scalar alphac(alphac_[p.cell()]);
 
-    scalar cbrtAlphap(pow(1.0 - alphac, 1.0/3.0));
+    scalar cbrtAlphap(cbrt(1.0 - alphac));
 
     scalar A =
         26.8*pow3(alphac)
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
index 7a2a3edf49..fba5e28ef3 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
@@ -2,8 +2,8 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.C
index 5ccaf9bc02..d219177107 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -123,7 +123,9 @@ Foam::PatchInteractionModel<CloudType>::PatchInteractionModel
 )
 :
     CloudSubModelBase<CloudType>(owner, dict, typeName, type),
-    UName_(this->coeffDict().lookupOrDefault("UName", word("U")))
+    UName_(this->coeffDict().lookupOrDefault("UName", word("U"))),
+    escapedParcels_(0),
+    escapedMass_(0.0)
 {}
 
 
@@ -134,7 +136,9 @@ Foam::PatchInteractionModel<CloudType>::PatchInteractionModel
 )
 :
     CloudSubModelBase<CloudType>(pim),
-    UName_(pim.UName_)
+    UName_(pim.UName_),
+    escapedParcels_(pim.escapedParcels_),
+    escapedMass_(pim.escapedMass_)
 {}
 
 
@@ -154,10 +158,42 @@ const Foam::word& Foam::PatchInteractionModel<CloudType>::UName() const
 }
 
 
+template<class CloudType>
+void Foam::PatchInteractionModel<CloudType>::addToEscapedParcels
+(
+    const scalar mass
+)
+{
+    escapedMass_ += mass;
+    escapedParcels_++;
+}
+
+
 template<class CloudType>
 void Foam::PatchInteractionModel<CloudType>::info(Ostream& os)
 {
-    // do nothing
+    const label escapedParcels0 =
+        this->template getBaseProperty<label>("escapedParcels");
+    const label escapedParcelsTotal =
+        escapedParcels0 + returnReduce(escapedParcels_, sumOp<label>());
+
+    const scalar escapedMass0 =
+        this->template getBaseProperty<scalar>("escapedMass");
+    const scalar escapedMassTotal =
+        escapedMass0 + returnReduce(escapedMass_, sumOp<scalar>());
+
+    os  << "    Parcel fate: system (number, mass)" << nl
+        << "      - escape                      = " << escapedParcelsTotal
+        << ", " << escapedMassTotal << endl;
+
+    if (this->outputTime())
+    {
+        this->setBaseProperty("escapedParcels", escapedParcelsTotal);
+        escapedParcels_ = 0;
+
+        this->setBaseProperty("escapedMass", escapedMassTotal);
+        escapedMass_ = 0.0;
+    }
 }
 
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.H b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.H
index bfc4ae277b..66fc753f07 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/PatchInteractionModel/PatchInteractionModel.H
@@ -74,7 +74,7 @@ public:
         static wordList interactionTypeNames_;
 
 
-private:
+protected:
 
     // Private data
 
@@ -82,6 +82,16 @@ private:
         const word UName_;
 
 
+        // Counters
+
+            //- Number of parcels escaped
+            label escapedParcels_;
+
+            //- Mass of parcels escaped
+            scalar escapedMass_;
+
+
+
 public:
 
     //- Runtime type information
@@ -158,6 +168,9 @@ public:
             const tetIndices& tetIs
         ) = 0;
 
+        //- Add to escaped parcels
+        void addToEscapedParcels(const scalar mass);
+
 
         // I-O
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/SurfaceFilmModel/SurfaceFilmModel/SurfaceFilmModel.C b/src/lagrangian/intermediate/submodels/Kinematic/SurfaceFilmModel/SurfaceFilmModel/SurfaceFilmModel.C
index d2dbb92d6c..3f3d997403 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/SurfaceFilmModel/SurfaceFilmModel/SurfaceFilmModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/SurfaceFilmModel/SurfaceFilmModel/SurfaceFilmModel.C
@@ -266,8 +266,9 @@ void Foam::SurfaceFilmModel<CloudType>::info(Ostream& os)
     label nInjectTotal =
         nInject0 + returnReduce(nParcelsInjected_, sumOp<label>());
 
-    os  << "    Parcels absorbed into film      = " << nTransTotal << nl
-        << "    New film detached parcels       = " << nInjectTotal << endl;
+    os  << "    Surface film:" << nl
+        << "      - parcels absorbed            = " << nTransTotal << nl
+        << "      - parcels ejected             = " << nInjectTotal << endl;
 
     if (this->outputTime())
     {
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
index 2b4e987a73..d06d08a3f3 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
@@ -250,7 +250,7 @@ const Foam::scalarField& Foam::CompositionModel<CloudType>::Y0
 
 
 template<class CloudType>
-Foam::scalarField Foam::CompositionModel<CloudType>::X
+Foam::tmp<Foam::scalarField> Foam::CompositionModel<CloudType>::X
 (
     const label phasei,
     const scalarField& Y
@@ -284,7 +284,7 @@ Foam::scalarField Foam::CompositionModel<CloudType>::X
         {
             FatalErrorIn
             (
-                "scalarField CompositionModel<CloudType>::X"
+                "tmp<scalarField> CompositionModel<CloudType>::X"
                 "("
                     "const label, "
                     "const scalarField&"
@@ -294,9 +294,8 @@ Foam::scalarField Foam::CompositionModel<CloudType>::X
         }
     }
 
-    X /= WInv;
-
-    return X;
+    tmp<scalarField> tfld = X/(WInv + ROOTVSMALL);
+    return tfld;
 }
 
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
index b57da79fde..97aa1c5c7c 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
@@ -188,7 +188,11 @@ public:
 
                 //- Return the list of phase phaseI volume fractions fractions
                 //  based on supplied mass fractions Y
-                scalarField X(const label phaseI, const scalarField& Y) const;
+                tmp<scalarField> X
+                (
+                    const label phaseI,
+                    const scalarField& Y
+                ) const;
 
 
             // Mixture properties
diff --git a/src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.C b/src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.C
index 2d47253de9..7bec31ed52 100644
--- a/src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.C
+++ b/src/lagrangian/intermediate/submodels/Thermodynamic/SurfaceFilmModel/ThermoSurfaceFilm/ThermoSurfaceFilm.C
@@ -686,7 +686,7 @@ void Foam::ThermoSurfaceFilm<CloudType>::info(Ostream& os)
     label nSplashTotal =
         nSplash0 + returnReduce(nParcelsSplashed_, sumOp<label>());
 
-    os  << "    New film splash parcels         = " << nSplashTotal << endl;
+    os  << "      - new splash parcels          = " << nSplashTotal << endl;
 
     if (this->outputTime())
     {
diff --git a/src/lagrangian/spray/Make/options b/src/lagrangian/spray/Make/options
index 3b5a4e72e9..b45129b838 100644
--- a/src/lagrangian/spray/Make/options
+++ b/src/lagrangian/spray/Make/options
@@ -19,6 +19,7 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
@@ -45,6 +46,7 @@ LIB_LIBS = \
     -lincompressibleTransportModels \
     -lregionModels \
     -lsurfaceFilmModels \
+    -ldynamicMesh \
     -ldynamicFvMesh \
     -lsampling \
     -lfiniteVolume \
diff --git a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
index ae1651e78d..f4f562543d 100644
--- a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
+++ b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
@@ -24,6 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "SprayParcel.H"
+#include "BreakupModel.H"
 #include "CompositionModel.H"
 #include "AtomizationModel.H"
 
@@ -153,13 +154,19 @@ void Foam::SprayParcel<ParcelType>::calcAtomization
     const label cellI
 )
 {
+    typedef typename TrackData::cloudType::sprayCloudType sprayCloudType;
+    const AtomizationModel<sprayCloudType>& atomization =
+        td.cloud().atomization();
+
+    if (!atomization.active())
+    {
+        return;
+    }
+
     typedef typename TrackData::cloudType::reactingCloudType reactingCloudType;
     const CompositionModel<reactingCloudType>& composition =
         td.cloud().composition();
 
-    typedef typename TrackData::cloudType::sprayCloudType sprayCloudType;
-    const AtomizationModel<sprayCloudType>& atomization =
-        td.cloud().atomization();
 
     // Average molecular weight of carrier mix - assumes perfect gas
     scalar Wc = this->rhoc_*RR*this->Tc()/this->pc();
@@ -223,11 +230,19 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
     typedef typename TrackData::cloudType cloudType;
     typedef typename cloudType::parcelType parcelType;
     typedef typename cloudType::forceType forceType;
+    typedef typename TrackData::cloudType::sprayCloudType sprayCloudType;
+
+    BreakupModel<sprayCloudType>& breakup = td.cloud().breakup();
+
+    if (!breakup.active())
+    {
+        return;
+    }
 
     const parcelType& p = static_cast<const parcelType&>(*this);
     const forceType& forces = td.cloud().forces();
 
-    if (td.cloud().breakup().solveOscillationEq())
+    if (breakup.solveOscillationEq())
     {
         solveTABEq(td, dt);
     }
@@ -247,7 +262,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
     const scalar mass = p.mass();
     const forceSuSp Fcp = forces.calcCoupled(p, dt, mass, Re, muAv);
     const forceSuSp Fncp = forces.calcNonCoupled(p, dt, mass, Re, muAv);
-    this->tMom() = mass/(Fcp.Sp() + Fncp.Sp());
+    this->tMom() = mass/(Fcp.Sp() + Fncp.Sp() + ROOTVSMALL);
 
     const vector g = td.cloud().g().value();
 
@@ -255,7 +270,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
     scalar dChild = 0.0;
     if
     (
-        td.cloud().breakup().update
+        breakup.update
         (
             dt,
             g,
diff --git a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.H b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.H
index 4e1fab22d0..4d7f289260 100644
--- a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.H
+++ b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.H
@@ -181,6 +181,25 @@ public:
         //- Runtime type information
         TypeName("SprayParcel");
 
+        //- String representation of properties
+        AddToPropertyList
+        (
+            ParcelType,
+          + " d0"
+          + " position0"
+          + " sigma"
+          + " mu"
+          + " liquidCore"
+          + " KHindex"
+          + " y"
+          + " yDot"
+          + " tc"
+          + " ms"
+          + " injector"
+          + " tMom"
+          + " user"
+        );
+
 
     // Constructors
 
diff --git a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcelIO.C b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcelIO.C
index 7d85441b1e..6188ccf735 100644
--- a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcelIO.C
+++ b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcelIO.C
@@ -28,6 +28,11 @@ License
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
+template<class ParcelType>
+Foam::string Foam::SprayParcel<ParcelType>::propertyList_ =
+    Foam::SprayParcel<ParcelType>::propertyList();
+
+
 template<class ParcelType>
 const std::size_t Foam::SprayParcel<ParcelType>::sizeofFields_
 (
@@ -89,7 +94,7 @@ Foam::SprayParcel<ParcelType>::SprayParcel
     (
         "SprayParcel<ParcelType>::SprayParcel"
         "("
-            "const polyMesh, "
+            "const polyMesh&, "
             "Istream&, "
             "bool"
         ")"