diff --git a/libraries/porousModels/relativePermeabilityModels/BrooksAndCorey/krBrooksAndCorey.C b/libraries/porousModels/relativePermeabilityModels/BrooksAndCorey/krBrooksAndCorey.C
index e4b5b9f8..3855588f 100644
--- a/libraries/porousModels/relativePermeabilityModels/BrooksAndCorey/krBrooksAndCorey.C
+++ b/libraries/porousModels/relativePermeabilityModels/BrooksAndCorey/krBrooksAndCorey.C
@@ -157,7 +157,7 @@ void Foam::relativePermeabilityModels::krBrooksAndCorey::correctkrb(const volSca
     krb_ = krbmax_ * pow(Se_,n_);
     if (derivative)
     {
-        dkradS_ = -kramax_*n_*pow((scalar(1)-Se_),n_-1)/(Smax_- Smin_);
+        dkrbdS_ = krbmax_*n_*pow(Se_,n_-1)/(Smax_- Smin_);
     }
 }
 void Foam::relativePermeabilityModels::krBrooksAndCorey::correctkrb(const volScalarField& Sb, const label& celli)
diff --git a/solvers/impesFoam/SEqn.H b/solvers/impesFoam/SEqn.H
index da853ae6..b6588c0b 100644
--- a/solvers/impesFoam/SEqn.H
+++ b/solvers/impesFoam/SEqn.H
@@ -16,13 +16,15 @@
   
     fvScalarMatrix SbEqn
         (
-            eps*fvm::ddt(Sb) + fvc::div(phib) 
+            // eps*fvm::ddt(Sb) + fvc::div(phib)
+            eps*fvm::ddt(Sb) + fvm::div(phib, Sb) 
             ==
             // event source terms
             - sourceTerm
         );
 
     SbEqn.solve();
+    Sb.correctBoundaryConditions();
 
     Info << "Saturation b: " << " Min(Sb) = " << gMin(Sb) << " Max(Sb) = " << gMax(Sb) << endl;