diff --git a/opm/autodiff/BlackoilModelEbos.hpp b/opm/autodiff/BlackoilModelEbos.hpp
index 9c4558fa1..edff9a97e 100644
--- a/opm/autodiff/BlackoilModelEbos.hpp
+++ b/opm/autodiff/BlackoilModelEbos.hpp
@@ -463,19 +463,6 @@ namespace Opm {
 
             auto& ebosJac = ebosSimulator_.model().linearizer().jacobian();
             auto& ebosResid = ebosSimulator_.model().linearizer().residual();
-            // J = [A, B; C, D], where A is the reservoir equations, B and C the interaction of well
-            // with the reservoir and D is the wells itself.
-            // The full system is reduced to a number of cells X number of cells system via Schur complement
-            // A -= B^T D^-1 C
-            // If matrix_add_well_contribution is false, the Ax operator is modified. i.e Ax -= B^T D^-1 C x in the WellModelMatrixAdapter
-            // instead of A.
-            // The residual is modified similarly.
-            // r = [r, r_well], where r is the residual and r_well the well residual.
-            // r -= B^T * D^-1 r_well
-            wellModel().apply(ebosResid);
-            if (param_.matrix_add_well_contributions_) {
-                wellModel().addWellContributions(ebosJac.istlMatrix());
-            }
 
             // set initial guess
             x = 0.0;
diff --git a/opm/autodiff/BlackoilWellModel_impl.hpp b/opm/autodiff/BlackoilWellModel_impl.hpp
index fff1a6218..aca7821c2 100644
--- a/opm/autodiff/BlackoilWellModel_impl.hpp
+++ b/opm/autodiff/BlackoilWellModel_impl.hpp
@@ -101,18 +101,23 @@ namespace Opm {
         if (!localWellsActive())
             return;
 
-        // we don't what to add the schur complement
-        // here since it affects the getConvergence method
-        /*
+        if (!param_.matrix_add_well_contributions_) {
+            // if the well contributions are not supposed to be included explicitly in
+            // the matrix, we only apply the vector part of the Schur complement here.
+            for (const auto& well: well_container_) {
+                // r = r - duneC_^T * invDuneD_ * resWell_
+                well->apply(res);
+            }
+            return;
+        }
+
         for (const auto& well: well_container_) {
-            if (param_.matrix_add_well_contributions_)
-                well->addWellContributions(mat);
+            well->addWellContributions(mat.istlMatrix());
 
             // applying the well residual to reservoir residuals
             // r = r - duneC_^T * invDuneD_ * resWell_
             well->apply(res);
         }
-        */
     }