diff --git a/src/lagrangian/intermediate/parcels/include/makeParcelInjectionModels.H b/src/lagrangian/intermediate/parcels/include/makeParcelInjectionModels.H
index cd5a78f81c..07b00f8d93 100644
--- a/src/lagrangian/intermediate/parcels/include/makeParcelInjectionModels.H
+++ b/src/lagrangian/intermediate/parcels/include/makeParcelInjectionModels.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-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -37,7 +37,7 @@ License
 #include "ManualInjection.H"
 #include "NoInjection.H"
 #include "PatchInjection.H"
-
+#include "PatchFlowRateInjection.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -53,7 +53,8 @@ License
     makeInjectionModelType(KinematicLookupTableInjection, CloudType);         \
     makeInjectionModelType(ManualInjection, CloudType);                       \
     makeInjectionModelType(NoInjection, CloudType);                           \
-    makeInjectionModelType(PatchInjection, CloudType);
+    makeInjectionModelType(PatchInjection, CloudType);                        \
+    makeInjectionModelType(PatchFlowRateInjection, CloudType);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/lagrangian/intermediate/parcels/include/makeReactingMultiphaseParcelInjectionModels.H b/src/lagrangian/intermediate/parcels/include/makeReactingMultiphaseParcelInjectionModels.H
index 126f8a8512..711f88459f 100644
--- a/src/lagrangian/intermediate/parcels/include/makeReactingMultiphaseParcelInjectionModels.H
+++ b/src/lagrangian/intermediate/parcels/include/makeReactingMultiphaseParcelInjectionModels.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-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,6 +35,7 @@ License
 #include "ManualInjection.H"
 #include "NoInjection.H"
 #include "PatchInjection.H"
+#include "PatchFlowRateInjection.H"
 #include "ReactingMultiphaseLookupTableInjection.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -49,6 +50,7 @@ License
     makeInjectionModelType(ManualInjection, CloudType);                       \
     makeInjectionModelType(NoInjection, CloudType);                           \
     makeInjectionModelType(PatchInjection, CloudType);                        \
+    makeInjectionModelType(PatchFlowRateInjection, CloudType);                \
     makeInjectionModelType(ReactingMultiphaseLookupTableInjection, CloudType);
 
 
diff --git a/src/lagrangian/intermediate/parcels/include/makeReactingParcelInjectionModels.H b/src/lagrangian/intermediate/parcels/include/makeReactingParcelInjectionModels.H
index 34fddbb50d..5fb4c6d562 100644
--- a/src/lagrangian/intermediate/parcels/include/makeReactingParcelInjectionModels.H
+++ b/src/lagrangian/intermediate/parcels/include/makeReactingParcelInjectionModels.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-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,6 +35,7 @@ License
 #include "ManualInjection.H"
 #include "NoInjection.H"
 #include "PatchInjection.H"
+#include "PatchFlowRateInjection.H"
 #include "ReactingLookupTableInjection.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -49,6 +50,7 @@ License
     makeInjectionModelType(ManualInjection, CloudType);                       \
     makeInjectionModelType(NoInjection, CloudType);                           \
     makeInjectionModelType(PatchInjection, CloudType);                        \
+    makeInjectionModelType(PatchFlowRateInjection, CloudType);                \
     makeInjectionModelType(ReactingLookupTableInjection, CloudType);
 
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.C
index 8ae15bab76..4764de1af6 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/CellZoneInjection/CellZoneInjection.C
@@ -77,7 +77,7 @@ void Foam::CellZoneInjection<CloudType>::setPositions
         for (label tetI = 1; tetI < cellTetIs.size() - 1; tetI++)
         {
             cTetVFrac[tetI] =
-                (cTetVFrac[tetI-1] + cellTetIs[tetI].tet(mesh).mag())/V[cellI];
+                cTetVFrac[tetI-1] + cellTetIs[tetI].tet(mesh).mag()/V[cellI];
         }
         cTetVFrac.last() = 1.0;
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
index 1c10870726..f1c6d2e524 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
@@ -45,7 +45,7 @@ bool Foam::InjectionModel<CloudType>::validInjection(const label parcelI)
 
 
 template<class CloudType>
-void Foam::InjectionModel<CloudType>::prepareForNextTimeStep
+bool Foam::InjectionModel<CloudType>::prepareForNextTimeStep
 (
     const scalar time,
     label& newParcels,
@@ -55,12 +55,13 @@ void Foam::InjectionModel<CloudType>::prepareForNextTimeStep
     // Initialise values
     newParcels = 0;
     newVolume = 0.0;
+    bool validInjection = false;
 
     // Return if not started injection event
     if (time < SOI_)
     {
         timeStep0_ = time;
-        return;
+        return validInjection;
     }
 
     // Make times relative to SOI
@@ -73,16 +74,27 @@ void Foam::InjectionModel<CloudType>::prepareForNextTimeStep
     // Volume of parcels to inject
     newVolume = this->volumeToInject(t0, t1);
 
-    // Hold previous time if no parcels, but non-zero volume fraction
-    if ((newParcels == 0) && (newVolume > 0.0))
+    if (newVolume > 0)
     {
-        // hold value of timeStep0_
+        if (newParcels > 0)
+        {
+            timeStep0_ = time;
+            validInjection = true;
+        }
+        else
+        {
+            // injection should have started, but not sufficient volume to
+            // produce (at least) 1 parcel - hold value of timeStep0_
+            validInjection = false;
+        }
     }
     else
     {
-        // advance value of timeStep0_
         timeStep0_ = time;
+        validInjection = false;
     }
+
+    return validInjection;
 }
 
 
@@ -416,7 +428,7 @@ Foam::label Foam::InjectionModel<CloudType>::parcelsToInject
         "("
             "const scalar, "
             "const scalar"
-        ") const"
+        ")"
     );
 
     return 0;
@@ -436,7 +448,7 @@ Foam::scalar Foam::InjectionModel<CloudType>::volumeToInject
         "("
             "const scalar, "
             "const scalar"
-        ") const"
+        ")"
     );
 
     return 0.0;
@@ -461,9 +473,6 @@ void Foam::InjectionModel<CloudType>::inject(TrackData& td)
     }
 
     const scalar time = this->owner().db().time().value();
-    const scalar trackTime = this->owner().solution().trackTime();
-    const polyMesh& mesh = this->owner().mesh();
-    typename TrackData::cloudType& cloud = td.cloud();
 
     // Prepare for next time step
     label parcelsAdded = 0;
@@ -471,96 +480,95 @@ void Foam::InjectionModel<CloudType>::inject(TrackData& td)
     label newParcels = 0;
     scalar newVolume = 0.0;
 
-    prepareForNextTimeStep(time, newParcels, newVolume);
-
-    // Duration of injection period during this timestep
-    const scalar deltaT =
-        max(0.0, min(trackTime, min(time - SOI_, timeEnd() - time0_)));
-
-    // Pad injection time if injection starts during this timestep
-    const scalar padTime = max(0.0, SOI_ - time0_);
-
-    // Introduce new parcels linearly across carrier phase timestep
-    for (label parcelI = 0; parcelI < newParcels; parcelI++)
+    if (prepareForNextTimeStep(time, newParcels, newVolume))
     {
-        if (validInjection(parcelI))
+        const scalar trackTime = this->owner().solution().trackTime();
+        const polyMesh& mesh = this->owner().mesh();
+        typename TrackData::cloudType& cloud = td.cloud();
+
+        // Duration of injection period during this timestep
+        const scalar deltaT =
+            max(0.0, min(trackTime, min(time - SOI_, timeEnd() - time0_)));
+
+        // Pad injection time if injection starts during this timestep
+        const scalar padTime = max(0.0, SOI_ - time0_);
+
+        // Introduce new parcels linearly across carrier phase timestep
+        for (label parcelI = 0; parcelI < newParcels; parcelI++)
         {
-            // Calculate the pseudo time of injection for parcel 'parcelI'
-            scalar timeInj = time0_ + padTime + deltaT*parcelI/newParcels;
-
-            // Determine the injection position and owner cell,
-            // tetFace and tetPt
-            label cellI = -1;
-            label tetFaceI = -1;
-            label tetPtI = -1;
-
-            vector pos = vector::zero;
-
-            setPositionAndCell
-            (
-                parcelI,
-                newParcels,
-                timeInj,
-                pos,
-                cellI,
-                tetFaceI,
-                tetPtI
-            );
-
-            if (cellI > -1)
+            if (validInjection(parcelI))
             {
-                // Lagrangian timestep
-                scalar dt = time - timeInj;
+                // Calculate the pseudo time of injection for parcel 'parcelI'
+                scalar timeInj = time0_ + padTime + deltaT*parcelI/newParcels;
 
-                // Apply corrections to position for 2-D cases
-                meshTools::constrainToMeshCentre(mesh, pos);
+                // Determine the injection position and owner cell,
+                // tetFace and tetPt
+                label cellI = -1;
+                label tetFaceI = -1;
+                label tetPtI = -1;
 
-                // Create a new parcel
-                parcelType* pPtr = new parcelType
+                vector pos = vector::zero;
+
+                setPositionAndCell
                 (
-                    td.cloud().pMesh(),
+                    parcelI,
+                    newParcels,
+                    timeInj,
                     pos,
                     cellI,
                     tetFaceI,
                     tetPtI
                 );
 
-                // Check/set new parcel thermo properties
-                cloud.setParcelThermoProperties(*pPtr, dt);
+                if (cellI > -1)
+                {
+                    // Lagrangian timestep
+                    scalar dt = time - timeInj;
 
-                // Assign new parcel properties in injection model
-                setProperties(parcelI, newParcels, timeInj, *pPtr);
+                    // Apply corrections to position for 2-D cases
+                    meshTools::constrainToMeshCentre(mesh, pos);
 
-                // Check/set new parcel injection properties
-                cloud.checkParcelProperties(*pPtr, dt, fullyDescribed());
+                    // Create a new parcel
+                    parcelType* pPtr =
+                        new parcelType(mesh, pos, cellI, tetFaceI, tetPtI);
 
-                // Apply correction to velocity for 2-D cases
-                meshTools::constrainDirection
-                (
-                    mesh,
-                    mesh.solutionD(),
-                    pPtr->U()
-                );
+                    // Check/set new parcel thermo properties
+                    cloud.setParcelThermoProperties(*pPtr, dt);
 
-                // Number of particles per parcel
-                pPtr->nParticle() =
-                    setNumberOfParticles
+                    // Assign new parcel properties in injection model
+                    setProperties(parcelI, newParcels, timeInj, *pPtr);
+
+                    // Check/set new parcel injection properties
+                    cloud.checkParcelProperties(*pPtr, dt, fullyDescribed());
+
+                    // Apply correction to velocity for 2-D cases
+                    meshTools::constrainDirection
                     (
-                        newParcels,
-                        newVolume,
-                        pPtr->d(),
-                        pPtr->rho()
+                        mesh,
+                        mesh.solutionD(),
+                        pPtr->U()
                     );
 
-                if (pPtr->move(td, dt))
-                {
-                    td.cloud().addParticle(pPtr);
-                    massAdded += pPtr->nParticle()*pPtr->mass();
-                    parcelsAdded++;
-                }
-                else
-                {
-                    delete pPtr;
+                    // Number of particles per parcel
+                    pPtr->nParticle() =
+                        setNumberOfParticles
+                        (
+                            newParcels,
+                            newVolume,
+                            pPtr->d(),
+                            pPtr->rho()
+                        );
+
+                    if (pPtr->move(td, dt))
+                    {
+                        td.cloud().addParticle(pPtr);
+                        massAdded += pPtr->nParticle()*pPtr->mass();
+                        parcelsAdded++;
+                    }
+                    else
+                    {
+                        delete pPtr;
+                    }
                 }
             }
         }
@@ -625,14 +633,8 @@ void Foam::InjectionModel<CloudType>::injectSteadyState
             meshTools::constrainToMeshCentre(mesh, pos);
 
             // Create a new parcel
-            parcelType* pPtr = new parcelType
-            (
-                td.cloud().pMesh(),
-                pos,
-                cellI,
-                tetFaceI,
-                tetPtI
-            );
+            parcelType* pPtr =
+                new parcelType(mesh, pos, cellI, tetFaceI, tetPtI);
 
             // Check/set new parcel thermo properties
             cloud.setParcelThermoProperties(*pPtr, 0.0);
@@ -644,12 +646,7 @@ void Foam::InjectionModel<CloudType>::injectSteadyState
             cloud.checkParcelProperties(*pPtr, 0.0, fullyDescribed());
 
             // Apply correction to velocity for 2-D cases
-            meshTools::constrainDirection
-            (
-                mesh,
-                mesh.solutionD(),
-                pPtr->U()
-            );
+            meshTools::constrainDirection(mesh, mesh.solutionD(), pPtr->U());
 
             // Number of particles per parcel
             pPtr->nParticle() =
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
index 30a2c97b16..a3d67947b2 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.H
@@ -138,7 +138,7 @@ protected:
         virtual bool validInjection(const label parcelI);
 
         //- Determine properties for next time step/injection interval
-        virtual void prepareForNextTimeStep
+        virtual bool prepareForNextTimeStep
         (
             const scalar time,
             label& newParcels,
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.C
new file mode 100644
index 0000000000..bf4bd87b1c
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.C
@@ -0,0 +1,316 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2012 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 "PatchFlowRateInjection.H"
+#include "TimeDataEntry.H"
+#include "distributionModel.H"
+#include "mathematicalConstants.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::PatchFlowRateInjection<CloudType>::PatchFlowRateInjection
+(
+    const dictionary& dict,
+    CloudType& owner,
+    const word& modelName
+)
+:
+    InjectionModel<CloudType>(dict, owner, modelName,typeName),
+    patchName_(this->coeffDict().lookup("patchName")),
+    patchId_(owner.mesh().boundaryMesh().findPatchID(patchName_)),
+    patchArea_(0.0),
+    patchNormal_(vector::zero),
+    phiName_(this->coeffDict().template lookupOrDefault<word>("phi", "phi")),
+    rhoName_(this->coeffDict().template lookupOrDefault<word>("rho", "rho")),
+    duration_(readScalar(this->coeffDict().lookup("duration"))),
+    concentration_(readScalar(this->coeffDict().lookup("concentration"))),
+    parcelsPerSecond_
+    (
+        readScalar(this->coeffDict().lookup("parcelsPerSecond"))
+    ),
+    U0_(vector::zero),
+    sizeDistribution_
+    (
+        distributionModels::distributionModel::New
+        (
+            this->coeffDict().subDict("sizeDistribution"),
+            owner.rndGen()
+        )
+    ),
+    cellOwners_(),
+    fraction_(1.0),
+    pMeanVolume_(0.0)
+{
+    if (patchId_ < 0)
+    {
+        FatalErrorIn
+        (
+            "PatchFlowRateInjection<CloudType>::PatchFlowRateInjection"
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
+        )   << "Requested patch " << patchName_ << " not found" << nl
+            << "Available patches are: " << owner.mesh().boundaryMesh().names()
+            << nl << exit(FatalError);
+    }
+
+    const polyPatch& patch = owner.mesh().boundaryMesh()[patchId_];
+
+    duration_ = owner.db().time().userTimeToTime(duration_);
+
+    cellOwners_ = patch.faceCells();
+
+    // TODO: retrieve mean diameter from distrution model
+    scalar pMeanDiameter =
+        readScalar(this->coeffDict().lookup("meanParticleDiameter"));
+    pMeanVolume_ = constant::mathematical::pi*pow3(pMeanDiameter)/6.0;
+
+    // patch geometry
+    label patchSize = cellOwners_.size();
+    label totalPatchSize = patchSize;
+    reduce(totalPatchSize, sumOp<label>());
+    fraction_ = scalar(patchSize)/totalPatchSize;
+
+    patchArea_ = gSum(mag(patch.faceAreas()));
+    patchNormal_ = gSum(patch.faceNormals())/totalPatchSize;
+    patchNormal_ /= mag(patchNormal_);
+
+    // Re-initialise total mass/volume to inject to zero
+    // - will be reset during each injection
+    this->volumeTotal_ = 0.0;
+    this->massTotal_ = 0.0;
+}
+
+
+template<class CloudType>
+Foam::PatchFlowRateInjection<CloudType>::PatchFlowRateInjection
+(
+    const PatchFlowRateInjection<CloudType>& im
+)
+:
+    InjectionModel<CloudType>(im),
+    patchName_(im.patchName_),
+    patchId_(im.patchId_),
+    patchArea_(im.patchArea_),
+    patchNormal_(im.patchNormal_),
+    phiName_(im.phiName_),
+    rhoName_(im.rhoName_),
+    duration_(im.duration_),
+    concentration_(im.concentration_),
+    parcelsPerSecond_(im.parcelsPerSecond_),
+    U0_(im.U0_),
+    sizeDistribution_(im.sizeDistribution_().clone().ptr()),
+    cellOwners_(im.cellOwners_),
+    fraction_(im.fraction_),
+    pMeanVolume_(im.pMeanVolume_)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::PatchFlowRateInjection<CloudType>::~PatchFlowRateInjection()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::scalar Foam::PatchFlowRateInjection<CloudType>::timeEnd() const
+{
+    return this->SOI_ + duration_;
+}
+
+
+template<class CloudType>
+Foam::label Foam::PatchFlowRateInjection<CloudType>::parcelsToInject
+(
+    const scalar time0,
+    const scalar time1
+)
+{
+    if ((time0 >= 0.0) && (time0 < duration_))
+    {
+        scalar nParcels = fraction_*(time1 - time0)*parcelsPerSecond_;
+
+        cachedRandom& rnd = this->owner().rndGen();
+
+        label nParcelsToInject = floor(nParcels);
+
+        // Inject an additional parcel with a probability based on the
+        // remainder after the floor function
+        if
+        (
+            nParcelsToInject > 0
+         && (
+               nParcels - scalar(nParcelsToInject)
+             > rnd.position(scalar(0), scalar(1))
+            )
+        )
+        {
+            ++nParcelsToInject;
+        }
+
+        return nParcelsToInject;
+    }
+    else
+    {
+        return 0;
+    }
+}
+
+
+template<class CloudType>
+Foam::scalar Foam::PatchFlowRateInjection<CloudType>::volumeToInject
+(
+    const scalar time0,
+    const scalar time1
+)
+{
+    scalar volume = 0.0;
+
+    if ((time0 >= 0.0) && (time0 < duration_))
+    {
+        const polyMesh& mesh = this->owner().mesh();
+
+        const surfaceScalarField& phi =
+            mesh.lookupObject<surfaceScalarField>(phiName_);
+
+        const scalarField& phip = phi.boundaryField()[patchId_];
+
+        scalar carrierVolume = 0.0;
+        if (phi.dimensions() == dimVelocity*dimArea)
+        {
+            const scalar flowRateIn = max(0.0, -sum(phip));
+            U0_ = -patchNormal_*flowRateIn/patchArea_;
+            carrierVolume = (time1 - time0)*flowRateIn;
+        }
+        else
+        {
+            const volScalarField& rho =
+                mesh.lookupObject<volScalarField>(rhoName_);
+            const scalarField& rhop = rho.boundaryField()[patchId_];
+
+            const scalar flowRateIn = max(0.0, -sum(phip/rhop));
+            U0_ = -patchNormal_*flowRateIn/patchArea_;
+            carrierVolume = (time1 - time0)*flowRateIn;
+        }
+
+        const scalar newParticles = concentration_*carrierVolume;
+
+        volume = pMeanVolume_*newParticles;
+    }
+
+    reduce(volume, sumOp<scalar>());
+
+    this->volumeTotal_ = volume;
+    this->massTotal_ = volume*this->owner().constProps().rho0();
+
+    return fraction_*volume;
+}
+
+
+template<class CloudType>
+void Foam::PatchFlowRateInjection<CloudType>::setPositionAndCell
+(
+    const label,
+    const label,
+    const scalar,
+    vector& position,
+    label& cellOwner,
+    label& tetFaceI,
+    label& tetPtI
+)
+{
+    if (cellOwners_.size() > 0)
+    {
+        cachedRandom& rnd = this->owner().rndGen();
+
+        label cellI = rnd.position<label>(0, cellOwners_.size() - 1);
+
+        cellOwner = cellOwners_[cellI];
+
+        // The position is between the face and cell centre, which could be
+        // in any tet of the decomposed cell, so arbitrarily choose the first
+        // face of the cell as the tetFace and the first point after the base
+        // point on the face as the tetPt.  The tracking will pick the cell
+        // consistent with the motion in the firsttracking step.
+        tetFaceI = this->owner().mesh().cells()[cellOwner][0];
+        tetPtI = 1;
+
+        // position perturbed between cell and patch face centres
+        const vector& pc = this->owner().mesh().C()[cellOwner];
+        const vector& pf =
+            this->owner().mesh().Cf().boundaryField()[patchId_][cellI];
+
+        const scalar a = rnd.sample01<scalar>();
+        const vector d = pf - pc;
+        position = pc + 0.5*a*d;
+    }
+    else
+    {
+        cellOwner = -1;
+        tetFaceI = -1;
+        tetPtI = -1;
+        // dummy position
+        position = pTraits<vector>::max;
+    }
+}
+
+
+template<class CloudType>
+void Foam::PatchFlowRateInjection<CloudType>::setProperties
+(
+    const label,
+    const label,
+    const scalar,
+    typename CloudType::parcelType& parcel
+)
+{
+    // set particle velocity
+    parcel.U() = U0_;
+
+    // set particle diameter
+    parcel.d() = sizeDistribution_->sample();
+}
+
+
+template<class CloudType>
+bool Foam::PatchFlowRateInjection<CloudType>::fullyDescribed() const
+{
+    return false;
+}
+
+
+template<class CloudType>
+bool Foam::PatchFlowRateInjection<CloudType>::validInjection(const label)
+{
+    return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.H
new file mode 100644
index 0000000000..7421ad2afd
--- /dev/null
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchFlowRateInjection/PatchFlowRateInjection.H
@@ -0,0 +1,204 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2012 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::PatchFlowRateInjection
+
+Description
+    Patch injection
+    - uses patch flow rate to determine concentration and velociy
+
+    - User specifies
+      - Total mass to inject
+      - Name of patch
+      - Injection duration
+      - Initial parcel velocity
+      - Injection target concentration/carrier volume flow rate
+    - Parcel diameters obtained by distribution model
+    - Parcels injected at cell centres adjacent to patch
+
+SourceFiles
+    PatchFlowRateInjection.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef PatchFlowRateInjection_H
+#define PatchFlowRateInjection_H
+
+#include "InjectionModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+template<class Type>
+class TimeDataEntry;
+
+class distributionModel;
+
+/*---------------------------------------------------------------------------*\
+                   Class PatchFlowRateInjection Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CloudType>
+class PatchFlowRateInjection
+:
+    public InjectionModel<CloudType>
+{
+    // Private data
+
+        //- Name of patch
+        const word patchName_;
+
+        //- Id of patch
+        const label patchId_;
+
+        //- Patch area
+        scalar patchArea_;
+
+        //- Patch normal direction
+        vector patchNormal_;
+
+        //- Name of carrier (mass or volume) flux field
+        const word phiName_;
+
+        //- Name of carrier density field
+        const word rhoName_;
+
+        //- Injection duration [s]
+        scalar duration_;
+
+        //- Concentration of particles to carrier [] (particles/m3)
+        const scalar concentration_;
+
+        //- Number of parcels to introduce per second []
+        const label parcelsPerSecond_;
+
+        //- Initial parcel velocity [m/s]
+        vector U0_;
+
+        //- Parcel size distribution model
+        const autoPtr<distributionModels::distributionModel> sizeDistribution_;
+
+        //- List of cell labels corresponding to injector positions
+        labelList cellOwners_;
+
+        //- Fraction of injection controlled by this processor
+        scalar fraction_;
+
+        //- Mean particle volume TODO: temporary measure - return from PDF
+        scalar pMeanVolume_;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("patchFlowRateInjection");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        PatchFlowRateInjection
+        (
+            const dictionary& dict,
+            CloudType& owner,
+            const word& modelName
+        );
+
+        //- Construct copy
+        PatchFlowRateInjection(const PatchFlowRateInjection<CloudType>& im);
+
+        //- Construct and return a clone
+        virtual autoPtr<InjectionModel<CloudType> > clone() const
+        {
+            return autoPtr<InjectionModel<CloudType> >
+            (
+                new PatchFlowRateInjection<CloudType>(*this)
+            );
+        }
+
+
+    //- Destructor
+    virtual ~PatchFlowRateInjection();
+
+
+    // Member Functions
+
+        //- Return the end-of-injection time
+        scalar timeEnd() const;
+
+        //- Number of parcels to introduce relative to SOI
+        virtual label parcelsToInject(const scalar time0, const scalar time1);
+
+        //- Volume of parcels to introduce relative to SOI
+        virtual scalar volumeToInject(const scalar time0, const scalar time1);
+
+
+        // Injection geometry
+
+            //- Set the injection position and owner cell, tetFace and tetPt
+            virtual void setPositionAndCell
+            (
+                const label parcelI,
+                const label nParcels,
+                const scalar time,
+                vector& position,
+                label& cellOwner,
+                label& tetFaceI,
+                label& tetPtI
+            );
+
+            virtual void setProperties
+            (
+                const label parcelI,
+                const label nParcels,
+                const scalar time,
+                typename CloudType::parcelType& parcel
+            );
+
+            //- Flag to identify whether model fully describes the parcel
+            virtual bool fullyDescribed() const;
+
+            //- Return flag to identify whether or not injection of parcelI is
+            //  permitted
+            virtual bool validInjection(const label parcelI);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "PatchFlowRateInjection.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
index aae52b8e22..bc33cdc06f 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C
@@ -141,7 +141,7 @@ Foam::label Foam::PatchInjection<CloudType>::parcelsToInject
 {
     if ((time0 >= 0.0) && (time0 < duration_))
     {
-        scalar nParcels =fraction_*(time1 - time0)*parcelsPerSecond_;
+        scalar nParcels = fraction_*(time1 - time0)*parcelsPerSecond_;
 
         cachedRandom& rnd = this->owner().rndGen();
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
index c1ce95e215..2350ee62f8 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.C
@@ -39,8 +39,26 @@ Foam::LocalInteraction<CloudType>::LocalInteraction
     nEscape_(patchData_.size(), 0),
     massEscape_(patchData_.size(), 0.0),
     nStick_(patchData_.size(), 0),
-    massStick_(patchData_.size(), 0.0)
+    massStick_(patchData_.size(), 0.0),
+    writeFields_(this->coeffDict().lookupOrDefault("writeFields", false)),
+    massEscapePtr_(NULL),
+    massStickPtr_(NULL)
 {
+    if (writeFields_)
+    {
+        word massEscapeName(this->owner().name() + "::massEscape");
+        word massStickName(this->owner().name() + "::massStick");
+        Info<< "    Interaction fields will be written to " << massEscapeName
+            << " and " << massStickName << endl;
+
+        (void)massEscape();
+        (void)massStick();
+    }
+    else
+    {
+        Info<< "    Interaction fields will not be written" << endl;
+    }
+
     // check that interactions are valid/specified
     forAll(patchData_, patchI)
     {
@@ -74,7 +92,10 @@ Foam::LocalInteraction<CloudType>::LocalInteraction
     nEscape_(pim.nEscape_),
     massEscape_(pim.massEscape_),
     nStick_(pim.nStick_),
-    massStick_(pim.massStick_)
+    massStick_(pim.massStick_),
+    writeFields_(pim.writeFields_),
+    massEscapePtr_(NULL),
+    massStickPtr_(NULL)
 {}
 
 
@@ -87,6 +108,64 @@ Foam::LocalInteraction<CloudType>::~LocalInteraction()
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
+template<class CloudType>
+Foam::volScalarField& Foam::LocalInteraction<CloudType>::massEscape()
+{
+    if (!massEscapePtr_.valid())
+    {
+        const fvMesh& mesh = this->owner().mesh();
+
+        massEscapePtr_.reset
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    this->owner().name() + "::massEscape",
+                    mesh.time().timeName(),
+                    mesh,
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::AUTO_WRITE
+                ),
+                mesh,
+                dimensionedScalar("zero", dimMass, 0.0)
+            )
+        );
+    }
+
+    return massEscapePtr_();
+}
+
+
+template<class CloudType>
+Foam::volScalarField& Foam::LocalInteraction<CloudType>::massStick()
+{
+    if (!massStickPtr_.valid())
+    {
+        const fvMesh& mesh = this->owner().mesh();
+
+        massStickPtr_.reset
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    this->owner().name() + "::massStick",
+                    mesh.time().timeName(),
+                    mesh,
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::AUTO_WRITE
+                ),
+                mesh,
+                dimensionedScalar("zero", dimMass, 0.0)
+            )
+        );
+    }
+
+    return massStickPtr_();
+}
+
+
 template<class CloudType>
 bool Foam::LocalInteraction<CloudType>::correct
 (
@@ -97,14 +176,13 @@ bool Foam::LocalInteraction<CloudType>::correct
     const tetIndices& tetIs
 )
 {
-    vector& U = p.U();
-
-    bool& active = p.active();
-
     label patchI = patchData_.applyToPatch(pp.index());
 
     if (patchI >= 0)
     {
+        vector& U = p.U();
+        bool& active = p.active();
+
         typename PatchInteractionModel<CloudType>::interactionType it =
             this->wordToInteractionType
             (
@@ -115,20 +193,36 @@ bool Foam::LocalInteraction<CloudType>::correct
         {
             case PatchInteractionModel<CloudType>::itEscape:
             {
+                scalar dm = p.mass()*p.nParticle();
+
                 keepParticle = false;
                 active = false;
                 U = vector::zero;
                 nEscape_[patchI]++;
-                massEscape_[patchI] += p.mass()*p.nParticle();
+                massEscape_[patchI] += dm;
+                if (writeFields_)
+                {
+                    label pI = pp.index();
+                    label fI = pp.whichFace(p.face());
+                    massEscape().boundaryField()[pI][fI] += dm;
+                }
                 break;
             }
             case PatchInteractionModel<CloudType>::itStick:
             {
+                scalar dm = p.mass()*p.nParticle();
+
                 keepParticle = true;
                 active = false;
                 U = vector::zero;
                 nStick_[patchI]++;
-                massStick_[patchI] += p.mass()*p.nParticle();
+                massStick_[patchI] += dm;
+                if (writeFields_)
+                {
+                    label pI = pp.index();
+                    label fI = pp.whichFace(p.face());
+                    massStick().boundaryField()[pI][fI] += dm;
+                }
                 break;
             }
             case PatchInteractionModel<CloudType>::itRebound:
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.H b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.H
index a3a28121f0..3a8c776776 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/PatchInteractionModel/LocalInteraction/LocalInteraction.H
@@ -34,6 +34,7 @@ Description
 
 #include "PatchInteractionModel.H"
 #include "patchInteractionDataList.H"
+#include "Switch.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -69,6 +70,16 @@ class LocalInteraction
             List<scalar> massStick_;
 
 
+        //- Flag to output data as fields
+        Switch writeFields_;
+
+        //- Mass escape field
+        autoPtr<volScalarField> massEscapePtr_;
+
+        //- Mass stick field
+        autoPtr<volScalarField> massStickPtr_;
+
+
 public:
 
     //- Runtime type information
@@ -99,6 +110,12 @@ public:
 
     // Member Functions
 
+        //- Return access to the massEscape field
+        volScalarField& massEscape();
+
+        //- Return access to the massStick field
+        volScalarField& massStick();
+
         //- Apply velocity correction
         //  Returns true if particle remains in same cell
         virtual bool correct
diff --git a/src/lagrangian/spray/parcels/include/makeSprayParcelInjectionModels.H b/src/lagrangian/spray/parcels/include/makeSprayParcelInjectionModels.H
index dc5bd6225e..dba5fa2e2d 100644
--- a/src/lagrangian/spray/parcels/include/makeSprayParcelInjectionModels.H
+++ b/src/lagrangian/spray/parcels/include/makeSprayParcelInjectionModels.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-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,6 +36,7 @@ License
 #include "ManualInjection.H"
 #include "NoInjection.H"
 #include "PatchInjection.H"
+#include "PatchFlowRateInjection.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -50,6 +51,7 @@ License
     makeInjectionModelType(InflationInjection, CloudType);                    \
     makeInjectionModelType(ManualInjection, CloudType);                       \
     makeInjectionModelType(NoInjection, CloudType);                           \
+    makeInjectionModelType(PatchFlowRateInjection, CloudType);                \
     makeInjectionModelType(PatchInjection, CloudType);