diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C
index ac60cf0d04..68d3edcf77 100644
--- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C
+++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C
@@ -53,6 +53,7 @@ Foam::porosityModels::solidification::solidification
 :
     porosityModel(name, modelType, mesh, dict, cellZoneName),
     TName_(coeffs_.lookupOrDefault<word>("T", "T")),
+    alphaName_(coeffs_.lookupOrDefault<word>("alpha", "none")),
     rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho")),
     D_(Function1<scalar>::New("D", coeffs_))
 {}
diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H
index 7f86849f64..119ba56805 100644
--- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H
+++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H
@@ -31,12 +31,13 @@ Description
     is represented as a porous blockage with the drag-coefficient evaluated from
 
         \f[
-            S = - \rho D(T) U
+            S = - \alpha \rho D(T) U
         \f]
 
     where
     \vartable
-        D(T) | User-defined drag-coefficient as function of temperature
+        \alpha  | Optional phase-fraction of solidifying phase
+        D(T)    | User-defined drag-coefficient as function of temperature
     \endvartable
 
     Note that the latent heat of solidification is not included and the
@@ -55,6 +56,9 @@ Description
                 (330.5     0)     // Liquid above 330.5K
             );
 
+            // Optional phase-fraction of solidifying phase
+            alpha alpha.liquid;
+
             // Solidification porosity is isotropic
             // use the global coordinate system
             coordinateSystem
@@ -103,6 +107,9 @@ class solidification
         //- Name of temperature field, default = "T"
         word TName_;
 
+        //- Name of optional phase-fraction field, default = "none"
+        word alphaName_;
+
         //- Name of density field, default = "rho"
         word rhoName_;
 
@@ -112,6 +119,27 @@ class solidification
 
     // Private Member Functions
 
+        //- Apply resistance
+        template<class AlphaFieldType, class RhoFieldType>
+        void apply
+        (
+            scalarField& Udiag,
+            const scalarField& V,
+            const AlphaFieldType& alpha,
+            const RhoFieldType& rho,
+            const volVectorField& U
+        ) const;
+
+        //- Apply resistance
+        template<class AlphaFieldType, class RhoFieldType>
+        void apply
+        (
+            tensorField& AU,
+            const AlphaFieldType& alpha,
+            const RhoFieldType& rho,
+            const volVectorField& U
+        ) const;
+
         //- Apply resistance
         template<class RhoFieldType>
         void apply
diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C
index 71b290af17..d9ea04d59e 100644
--- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C
+++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C
@@ -24,14 +24,16 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "volFields.H"
+#include "geometricOneField.H"
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-template<class RhoFieldType>
+template<class AlphaFieldType, class RhoFieldType>
 void Foam::porosityModels::solidification::apply
 (
     scalarField& Udiag,
     const scalarField& V,
+    const AlphaFieldType& alpha,
     const RhoFieldType& rho,
     const volVectorField& U
 ) const
@@ -48,12 +50,66 @@ void Foam::porosityModels::solidification::apply
         forAll(cells, i)
         {
             const label celli = cells[i];
-            Udiag[celli] += V[celli]*rho[celli]*D_->value(T[celli]);
+            Udiag[celli] +=
+                V[celli]*alpha[celli]*rho[celli]*D_->value(T[celli]);
         }
     }
 }
 
 
+template<class AlphaFieldType, class RhoFieldType>
+void Foam::porosityModels::solidification::apply
+(
+    tensorField& AU,
+    const AlphaFieldType& alpha,
+    const RhoFieldType& rho,
+    const volVectorField& U
+) const
+{
+    const volScalarField& T = mesh_.lookupObject<volScalarField>
+    (
+        IOobject::groupName(TName_, U.group())
+    );
+
+    forAll(cellZoneIDs_, zoneI)
+    {
+        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+
+        forAll(cells, i)
+        {
+            const label celli = cells[i];
+            AU[celli] +=
+                tensor::I*alpha[celli]*rho[celli]*D_->value(T[celli]);
+        }
+    }
+}
+
+
+template<class RhoFieldType>
+void Foam::porosityModels::solidification::apply
+(
+    scalarField& Udiag,
+    const scalarField& V,
+    const RhoFieldType& rho,
+    const volVectorField& U
+) const
+{
+    if (alphaName_ == "none")
+    {
+        return apply(Udiag, V, geometricOneField(), rho, U);
+    }
+    else
+    {
+        const volScalarField& alpha = mesh_.lookupObject<volScalarField>
+        (
+            IOobject::groupName(alphaName_, U.group())
+        );
+
+        return apply(Udiag, V, alpha, rho, U);
+    }
+}
+
+
 template<class RhoFieldType>
 void Foam::porosityModels::solidification::apply
 (
@@ -62,20 +118,18 @@ void Foam::porosityModels::solidification::apply
     const volVectorField& U
 ) const
 {
-    const volScalarField& T = mesh_.lookupObject<volScalarField>
-    (
-        IOobject::groupName(TName_, U.group())
-    );
-
-    forAll(cellZoneIDs_, zoneI)
+    if (alphaName_ == "none")
     {
-        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
+        return apply(AU, geometricOneField(), rho, U);
+    }
+    else
+    {
+        const volScalarField& alpha = mesh_.lookupObject<volScalarField>
+        (
+            IOobject::groupName(alphaName_, U.group())
+        );
 
-        forAll(cells, i)
-        {
-            const label celli = cells[i];
-            AU[celli] += tensor::I*rho[celli]*D_->value(T[celli]);
-        }
+        return apply(AU, alpha, rho, U);
     }
 }