diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
index e654b56..f97328b 100644
--- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
@@ -181,19 +181,18 @@ const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity
 
     // Calculate the integration coefficients
     const vector acp = (Fcp.Sp()*td.Uc() + Fcp.Su())/massEff;
-    const vector ancp = (Fncp.Sp()*td.Uc() + Fncp.Su() + Su)/massEff;
+    const vector ancp = (Fncp.Su() + Su)/massEff;
     const scalar bcp = Fcp.Sp()/massEff;
-    const scalar bncp = Fncp.Sp()/massEff;
+    //const scalar bncp = Fncp.Sp()/massEff;
 
     // Integrate to find the new parcel velocity
-    const scalar dtEff = cloud.UIntegrator().dtEff(dt, bcp + bncp);
-    const vector deltaUcp = integrationScheme::delta(U_, dtEff, acp, bcp);
-    const vector deltaUncp = integrationScheme::delta(U_, dtEff, ancp, bncp);
+    const scalar dtEff = cloud.UIntegrator().dtEff(dt, bcp);
+    const vector deltaU = integrationScheme::delta(U_, dtEff, acp+ancp, bcp);
 
     // Calculate the new velocity and the momentum transfer terms
-    vector Unew = U_ + deltaUcp + deltaUncp;
+    vector Unew = U_ + deltaU;
 
-    dUTrans -= massEff*deltaUcp;
+    dUTrans -= massEff*(deltaU-ancp*dt);
 
     Spu = dt*Fcp.Sp();
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
index d44526e..3f2d064 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
@@ -287,14 +287,13 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer
 
     // Integrate to find the new parcel temperature
     const scalar dtEff = cloud.TIntegrator().dtEff(dt, bcp);
-    const scalar deltaTcp = integrationScheme::delta(T_, dtEff, acp, bcp);
-    const scalar deltaTncp = integrationScheme::delta(T_, dtEff, ancp, 0);
+    const scalar deltaT = integrationScheme::delta(T_, dtEff, acp+ancp, bcp);
 
     // Calculate the new temperature and the enthalpy transfer terms
-    scalar Tnew = T_ + deltaTcp + deltaTncp;
+    scalar Tnew = T_ + deltaT;
     Tnew = min(max(Tnew, cloud.constProps().TMin()), cloud.constProps().TMax());
 
-    dhsTrans -= m*Cp_*deltaTcp;
+    dhsTrans -= m*Cp_*(deltaT-ancp*dt);
 
     Sph = dt*m*Cp_*bcp;