diff --git a/opm/autodiff/BlackoilModelEbos.hpp b/opm/autodiff/BlackoilModelEbos.hpp
index d20a76773..af57eadd9 100644
--- a/opm/autodiff/BlackoilModelEbos.hpp
+++ b/opm/autodiff/BlackoilModelEbos.hpp
@@ -565,11 +565,8 @@ namespace Opm {
           {
             A_.mv( x, y );
 
-            if ( ! matrix_add_well_contributions_ )
-            {
-              // add well model modification to y
-              wellMod_.apply(x, y );
-            }
+            // add well model modification to y
+            wellMod_.apply(x, y );
 
 #if HAVE_MPI
             if( comm_ )
@@ -582,11 +579,8 @@ namespace Opm {
           {
             A_.usmv(alpha,x,y);
 
-            if ( ! matrix_add_well_contributions_ )
-            {
-              // add scaled well model modification to y
-              wellMod_.applyScaleAdd( alpha, x, y );
-            }
+            // add scaled well model modification to y
+            wellMod_.applyScaleAdd( alpha, x, y );
 
 #if HAVE_MPI
             if( comm_ )
diff --git a/opm/autodiff/BlackoilWellModel.hpp b/opm/autodiff/BlackoilWellModel.hpp
index 3af0d06d5..cba3142c3 100644
--- a/opm/autodiff/BlackoilWellModel.hpp
+++ b/opm/autodiff/BlackoilWellModel.hpp
@@ -157,12 +157,6 @@ namespace Opm {
 
             const SimulatorReport& lastReport() const;
 
-            /// \! brief Modifies matrix to include influences of the well perforations.
-            ///
-            /// \param mat The linear system with the assembled mass balance
-            ///            equations
-            void addWellContributions(Mat& mat) const;
-
 
         protected:
 
diff --git a/opm/autodiff/BlackoilWellModel_impl.hpp b/opm/autodiff/BlackoilWellModel_impl.hpp
index 6219a0769..75ae1418b 100644
--- a/opm/autodiff/BlackoilWellModel_impl.hpp
+++ b/opm/autodiff/BlackoilWellModel_impl.hpp
@@ -294,18 +294,6 @@ namespace Opm {
         }
     }
 
-
-    template<typename TypeTag>
-    void
-    BlackoilWellModel<TypeTag>::
-    addWellContributions(Mat& mat) const
-    {
-        for(const auto& well: well_container_)
-        {
-            well->addWellContributions(mat);
-        }
-    }
-
     // applying the well residual to reservoir residuals
     // r = r - duneC_^T * invDuneD_ * resWell_
     template<typename TypeTag>
@@ -333,7 +321,8 @@ namespace Opm {
     apply(const BVector& x, BVector& Ax) const
     {
         // TODO: do we still need localWellsActive()?
-        if ( ! localWellsActive() ) {
+        if ( ! localWellsActive() ||
+             well_container_[0]->jacobianContainsWellContributions() ) {
             return;
         }
 
@@ -352,7 +341,8 @@ namespace Opm {
     BlackoilWellModel<TypeTag>::
     applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const
     {
-        if ( ! localWellsActive() ) {
+        if ( ! localWellsActive() ||
+             well_container_[0]->jacobianContainsWellContributions() ) {
             return;
         }
 
diff --git a/opm/autodiff/StandardWell.hpp b/opm/autodiff/StandardWell.hpp
index ac1e18502..e741c236b 100644
--- a/opm/autodiff/StandardWell.hpp
+++ b/opm/autodiff/StandardWell.hpp
@@ -151,7 +151,14 @@ namespace Opm
 
         virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
                                                  const WellState& well_state); // should be const?
+
         virtual void  addWellContributions(Mat& mat) const;
+
+        /// \brief Wether the Jacobian will also have well contributions in it.
+        virtual bool jacobianContainsWellContributions() const
+        {
+            return param_.matrix_add_well_contributions_;
+        }
     protected:
 
         // protected functions from the Base class
diff --git a/opm/autodiff/WellInterface.hpp b/opm/autodiff/WellInterface.hpp
index 3842968ce..d8ef760df 100644
--- a/opm/autodiff/WellInterface.hpp
+++ b/opm/autodiff/WellInterface.hpp
@@ -203,7 +203,11 @@ namespace Opm
         virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
                                                  const WellState& well_state) = 0; // should be const?
 
-        virtual void addWellContributions(Mat& mat) const;
+        /// \brief Wether the Jacobian will also have well contributions in it.
+        virtual bool jacobianContainsWellContributions() const
+        {
+            return false;
+        }
 
         // updating the voidage rates in well_state when requested
         void calculateReservoirRates(WellState& well_state) const;
diff --git a/opm/autodiff/WellInterface_impl.hpp b/opm/autodiff/WellInterface_impl.hpp
index 149b8d994..e449628ff 100644
--- a/opm/autodiff/WellInterface_impl.hpp
+++ b/opm/autodiff/WellInterface_impl.hpp
@@ -844,14 +844,6 @@ namespace Opm
         return 1.0;
     }
 
-    template<typename TypeTag>
-    void
-    WellInterface<TypeTag>::addWellContributions(Mat& mat) const
-    {
-        OPM_THROW(NotImplemented, "This well class does not support adding well contributions"
-                  << "to the matrix");     
-    }
-