diff --git a/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.C b/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.C
index 83bd7be2dd..7d3eebc3e4 100644
--- a/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.C
+++ b/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.C
@@ -90,7 +90,7 @@ Foam::volScalarField& Foam::functionObjects::phaseScalarTransport::Phi()
             (
                 IOobject
                 (
-                    "Phi" + s_.name(),
+                    typedName(IOobject::groupName("Phi", phaseName_)),
                     time_.name(),
                     mesh_,
                     IOobject::READ_IF_PRESENT,
@@ -112,16 +112,11 @@ Foam::volScalarField& Foam::functionObjects::phaseScalarTransport::Phi()
 Foam::tmp<Foam::surfaceScalarField>
 Foam::functionObjects::phaseScalarTransport::alphaPhi()
 {
-    // If alphaPhi exists then return it
-    if (mesh_.foundObject<surfaceScalarField>(alphaPhiName_))
+    if (!solveAlphaPhi_)
     {
         return mesh_.lookupObject<surfaceScalarField>(alphaPhiName_);
     }
 
-    // Otherwise generate it ...
-    Info<< type() << ": " << surfaceScalarField::typeName << " "
-        << alphaPhiName_ << " was not found, so generating it" << endl;
-
     const volScalarField& alpha =
         mesh_.lookupObject<volScalarField>(alphaName_);
     const surfaceScalarField& phi =
@@ -288,8 +283,6 @@ Foam::functionObjects::phaseScalarTransport::phaseScalarTransport
     fvMeshFunctionObject(name, runTime, dict),
     fieldName_(dict.lookup("field")),
     phaseName_(IOobject::group(fieldName_)),
-    nCorr_(0),
-    residualAlpha_(rootSmall),
     s_
     (
         IOobject
@@ -330,18 +323,20 @@ bool Foam::functionObjects::phaseScalarTransport::read(const dictionary& dict)
 {
     fvMeshFunctionObject::read(dict);
 
+    solveAlphaPhi_ = dict.lookupOrDefault<bool>("solveAlphaPhi", false);
+
     alphaName_ =
         dict.lookupOrDefault<word>
         (
             "alpha",
             IOobject::groupName("alpha", phaseName_)
         );
+    const word defaultAlphaPhiName =
+        IOobject::groupName("alphaPhi", phaseName_);
     alphaPhiName_ =
-        dict.lookupOrDefault<word>
-        (
-            "alphaPhi",
-            IOobject::groupName("alphaPhi", phaseName_)
-        );
+        solveAlphaPhi_
+      ? typedName(defaultAlphaPhiName)
+      : dict.lookupOrDefault<word>("alphaPhi", defaultAlphaPhiName);
     phiName_ = dict.lookupOrDefault<word>("phi", "phi");
     rhoName_ =
         dict.lookupOrDefault<word>
@@ -356,8 +351,8 @@ bool Foam::functionObjects::phaseScalarTransport::read(const dictionary& dict)
     alphaD_ = dict.lookupOrDefault<scalar>("alphaD", 1);
     alphaDt_ = dict.lookupOrDefault<scalar>("alphaDt", 1);
 
-    dict.readIfPresent("nCorr", nCorr_);
-    dict.readIfPresent("residualAlpha", residualAlpha_);
+    nCorr_ = dict.lookupOrDefault<label>("nCorr", 0);
+    residualAlpha_ = dict.lookupOrDefault<scalar>("residualAlpha", rootSmall);
     writeAlphaField_ = dict.lookupOrDefault<bool>("writeAlphaField", true);
 
     return true;
diff --git a/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.H b/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.H
index 2a4ca1450d..e851e2bd3f 100644
--- a/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.H
+++ b/src/functionObjects/solvers/phaseScalarTransport/phaseScalarTransport.H
@@ -37,27 +37,25 @@ Description
     detailed below. Note that the phase-name will be determined by stripping
     the extension from the supplied field name.
 
-    If \c alphaPhi is specified and found in the database then it will be used
-    to transport the field. This is the likely mode of operation if this
-    function is used with Euler-Euler solvers, in which \c alphaPhi -like
-    fluxes are available. If \c alphaPhi is not found, then a pressure-like
-    equation will be solved in order to construct it so that it exactly matches
-    the time-derivative of the phase volume or mass. This is likely to be
-    necessary in volume-of-fluid solvers where \c alphaPhi is not part of the
-    solution procedure. The pressure field name will be required in this case.
+    If the solver does not provide an \c alphaPhi flux, or that flux is for
+    some reason unreliable, then the \c solveAlphaPhi switch can be used to
+    make this function solve a pressure-like equation from which \c alphaPhi is
+    recovered.
 
 Usage
     \table
-        Property        | Description                        | Req'd? | Default
-        alpha           | Name of the volume-fraction field  | no\\
+        Property        | Description                       | Req'd? | Default
+        alpha           | Name of the volume-fraction field | no\\
                         | alpha.<phase-name>
-        alphaPhi        | Name of the phase-flux field       | no\\
+        alphaPhi        | Name of the phase-flux field      | no\\
                         | alphaPhi.<phase-name>
-        p               | Name of the pressure field         | no     | p
+        solveAlphaPhi   | Solve for the alphaPhi flux, rather than looking it\\
+                          up                                | no     | false
+        p               | Name of the pressure field        | no     | p
         residualAlpha   | Small volume fraction used to stabilise the solution\\
-                        | no | rootSmall
+                                                            | no     | rootSmall
         writeAlphaField | Also write out alpha multiplied by the field\\
-                        | no | true
+                                                            |        no | true
     \endtable
 
     Example specification for incompressibleVoF:
@@ -68,7 +66,6 @@ Usage
         libs            ("libsolverFunctionObjects.so");
 
         field           s.water;
-        p               p_rgh;
     }
     \endverbatim
 
@@ -123,6 +120,9 @@ class phaseScalarTransport
         //- Name of the phase in which to solve
         const word phaseName_;
 
+        //- Solve for the alphaPhi flux, rather than looking it up (optional)
+        bool solveAlphaPhi_;
+
         //- Name of phase volume-fraction field (optional)
         word alphaName_;
 
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.air.orig b/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.air.orig
new file mode 100644
index 0000000000..6b05c88fc0
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.air.orig
@@ -0,0 +1,39 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      tracer.air;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    #includeEtc "caseDicts/setConstraintTypes"
+
+    wall
+    {
+        type            zeroGradient;
+    }
+
+    atmosphere
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.water.orig b/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.water.orig
new file mode 100644
index 0000000000..2244f7678c
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/0/tracer.water.orig
@@ -0,0 +1,39 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      tracer.water;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    #includeEtc "caseDicts/setConstraintTypes"
+
+    wall
+    {
+        type            zeroGradient;
+    }
+
+    atmosphere
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/Allclean b/tutorials/incompressibleVoF/damBreakLaminar/Allclean
index 7b06dc70e3..9068aa2e07 100755
--- a/tutorials/incompressibleVoF/damBreakLaminar/Allclean
+++ b/tutorials/incompressibleVoF/damBreakLaminar/Allclean
@@ -6,10 +6,12 @@ cd ${0%/*} || exit 1    # Run from this directory
 
 cleanVoFCase
 
-rm -rf 0 system
+find 0 -type f -not -name tracer.* -delete
 
 find constant -type f -not \( \
     -name fvModels -or -name momentumTransport \
     \) -delete
 
+rm -f system/blockMeshDict
+
 #------------------------------------------------------------------------------
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/Allrun b/tutorials/incompressibleVoF/damBreakLaminar/Allrun
index 6bb259ecec..8064eb76ca 100755
--- a/tutorials/incompressibleVoF/damBreakLaminar/Allrun
+++ b/tutorials/incompressibleVoF/damBreakLaminar/Allrun
@@ -9,10 +9,6 @@ isTest "$@" && path=.. || path=$FOAM_TUTORIALS/incompressibleVoF
 cp -rn $path/damBreak/0 $path/damBreak/constant $path/damBreak/system .
 rm -f 0/alpha.water
 
-# Increase the number of limiter iterations
-runApplication foamDictionary system/fvSolution \
-    -entry solvers/"alpha.water.*"/nLimiterIter -set 5
-
 # Run
 runApplication blockMesh
 runApplication setFields
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/constant/momentumTransport b/tutorials/incompressibleVoF/damBreakLaminar/constant/momentumTransport
index c7745a3880..666fcc4955 100644
--- a/tutorials/incompressibleVoF/damBreakLaminar/constant/momentumTransport
+++ b/tutorials/incompressibleVoF/damBreakLaminar/constant/momentumTransport
@@ -16,12 +16,5 @@ FoamFile
 
 simulationType  laminar;
 
-RAS
-{
-    model           kEpsilon;
-    turbulence      on;
-    printCoeffs     on;
-}
-
 
 // ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/system/controlDict b/tutorials/incompressibleVoF/damBreakLaminar/system/controlDict
new file mode 100644
index 0000000000..8e235f7f13
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/system/controlDict
@@ -0,0 +1,62 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      controlDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+application     foamRun;
+
+solver          incompressibleVoF;
+
+startFrom       startTime;
+
+startTime       0;
+
+stopAt          endTime;
+
+endTime         1;
+
+deltaT          0.001;
+
+writeControl    adjustableRunTime;
+
+writeInterval   0.05;
+
+purgeWrite      0;
+
+writeFormat     ascii;
+
+writePrecision  6;
+
+writeCompression off;
+
+timeFormat      general;
+
+timePrecision   6;
+
+runTimeModifiable yes;
+
+adjustTimeStep  yes;
+
+maxCo           1;
+maxAlphaCo      1;
+
+maxDeltaT       1;
+
+functions
+{
+    #includeFunc phaseScalarTransport(field=tracer.water)
+    #includeFunc phaseScalarTransport(field=tracer.air, residualAlpha=1e-6)
+}
+
+// ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/system/fvSchemes b/tutorials/incompressibleVoF/damBreakLaminar/system/fvSchemes
new file mode 100644
index 0000000000..d96420f79b
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/system/fvSchemes
@@ -0,0 +1,58 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSchemes;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+ddtSchemes
+{
+    default         Euler;
+}
+
+gradSchemes
+{
+    default         Gauss linear;
+}
+
+divSchemes
+{
+    div(phi,alpha)  Gauss interfaceCompression vanLeer 1;
+
+    div(rhoPhi,U)   Gauss linearUpwind grad(U);
+
+    div(phi,k)      Gauss upwind;
+    div(phi,epsilon) Gauss upwind;
+    div(phi,omega)  Gauss upwind;
+
+    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
+
+    div(alphaPhi,tracer) Gauss vanLeer;
+}
+
+laplacianSchemes
+{
+    default         Gauss linear uncorrected;
+}
+
+interpolationSchemes
+{
+    default         linear;
+}
+
+snGradSchemes
+{
+    default         uncorrected;
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/system/fvSolution b/tutorials/incompressibleVoF/damBreakLaminar/system/fvSolution
new file mode 100644
index 0000000000..ff81240933
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/system/fvSolution
@@ -0,0 +1,90 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSolution;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+solvers
+{
+    "alpha.water.*"
+    {
+        nAlphaCorr      2;
+        nAlphaSubCycles 1;
+
+        MULESCorr       yes;
+        nLimiterIter    5;
+
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
+        tolerance       1e-8;
+        relTol          0;
+    }
+
+    "pcorr.*"
+    {
+        solver          PCG;
+        preconditioner  DIC;
+        tolerance       1e-5;
+        relTol          0;
+    }
+
+    p_rgh
+    {
+        solver          PCG;
+        preconditioner  DIC;
+        tolerance       1e-07;
+        relTol          0.05;
+    }
+
+    p_rghFinal
+    {
+        $p_rgh;
+        relTol          0;
+    }
+
+    "(U|k|epsilon|omega).*"
+    {
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
+        tolerance       1e-06;
+        relTol          0;
+        minIter         1;
+    }
+
+    "tracer.*"
+    {
+        solver          PBiCGStab;
+        preconditioner  DILU;
+        tolerance       1e-06;
+        relTol          0;
+    }
+}
+
+PIMPLE
+{
+    momentumPredictor   no;
+    nOuterCorrectors    1;
+    nCorrectors         3;
+    nNonOrthogonalCorrectors 0;
+}
+
+relaxationFactors
+{
+    equations
+    {
+        ".*"                1;
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/incompressibleVoF/damBreakLaminar/system/setFieldsDict b/tutorials/incompressibleVoF/damBreakLaminar/system/setFieldsDict
new file mode 100644
index 0000000000..936e153d7b
--- /dev/null
+++ b/tutorials/incompressibleVoF/damBreakLaminar/system/setFieldsDict
@@ -0,0 +1,53 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  dev
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      setFieldsDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+defaultFieldValues
+(
+    volScalarFieldValue alpha.water 0
+    volScalarFieldValue tracer.water 0
+    volScalarFieldValue tracer.air 0
+);
+
+regions
+(
+    boxToCell
+    {
+        box (0 0 -1) (0.1461 0.292 1);
+        fieldValues
+        (
+            volScalarFieldValue alpha.water 1
+        );
+    }
+    boxToCell
+    {
+        box (0 0 -1) (0.1461 0.146 1);
+        fieldValues
+        (
+            volScalarFieldValue tracer.water 1
+        );
+    }
+    boxToCell
+    {
+        box (0 0.292 -1) (0.1461 0.438 1);
+        fieldValues
+        (
+            volScalarFieldValue tracer.air 1
+        );
+    }
+);
+
+
+// ************************************************************************* //