diff --git a/MaterialLib/FractureModels/Coulomb.cpp b/MaterialLib/FractureModels/Coulomb.cpp
index 7afad4a60a5..d27f8a8a223 100644
--- a/MaterialLib/FractureModels/Coulomb.cpp
+++ b/MaterialLib/FractureModels/Coulomb.cpp
@@ -185,13 +185,9 @@ void Coulomb<DisplacementDim>::computeConstitutiveRelation(
             sigma.noalias() += sigma_prev;
         };
 
-        auto newton_solver =
-            NumLib::NewtonRaphson<decltype(linear_solver), JacobianMatrix,
-                                  decltype(update_jacobian), ResidualVectorType,
-                                  decltype(update_residual),
-                                  decltype(update_solution)>(
-                linear_solver, update_jacobian, update_residual,
-                update_solution, _nonlinear_solver_parameters);
+        auto newton_solver = NumLib::NewtonRaphson(
+            linear_solver, update_jacobian, update_residual, update_solution,
+            _nonlinear_solver_parameters);
 
         auto const success_iterations = newton_solver.solve(jacobian);
 
diff --git a/MaterialLib/SolidModels/CreepBGRa.cpp b/MaterialLib/SolidModels/CreepBGRa.cpp
index 7fa4b426f52..f1827eead35 100644
--- a/MaterialLib/SolidModels/CreepBGRa.cpp
+++ b/MaterialLib/SolidModels/CreepBGRa.cpp
@@ -113,12 +113,8 @@ CreepBGRa<DisplacementDim>::integrateStress(
     { solution += increment; };
 
     auto newton_solver =
-        NumLib::NewtonRaphson<decltype(linear_solver), JacobianMatrix,
-                              decltype(update_jacobian), ResidualVectorType,
-                              decltype(update_residual),
-                              decltype(update_solution)>(
-            linear_solver, update_jacobian, update_residual, update_solution,
-            _nonlinear_solver_parameters);
+        NumLib::NewtonRaphson(linear_solver, update_jacobian, update_residual,
+                              update_solution, _nonlinear_solver_parameters);
 
     JacobianMatrix jacobian;
     auto const success_iterations = newton_solver.solve(jacobian);
diff --git a/MaterialLib/SolidModels/Ehlers.cpp b/MaterialLib/SolidModels/Ehlers.cpp
index 160a353d7b0..e9c26e117e2 100644
--- a/MaterialLib/SolidModels/Ehlers.cpp
+++ b/MaterialLib/SolidModels/Ehlers.cpp
@@ -657,10 +657,7 @@ SolidEhlers<DisplacementDim>::integrateStress(
                     mp.G * solution.template segment<KelvinVectorSize>(0)};
             };
 
-            auto newton_solver = NumLib::NewtonRaphson<
-                decltype(linear_solver), JacobianMatrix,
-                decltype(update_jacobian), ResidualVectorType,
-                decltype(update_residual), decltype(update_solution)>(
+            auto newton_solver = NumLib::NewtonRaphson(
                 linear_solver, update_jacobian, update_residual,
                 update_solution, _nonlinear_solver_parameters);
 
diff --git a/MaterialLib/SolidModels/Lubby2.cpp b/MaterialLib/SolidModels/Lubby2.cpp
index a7617dc36ed..1d6fe281989 100644
--- a/MaterialLib/SolidModels/Lubby2.cpp
+++ b/MaterialLib/SolidModels/Lubby2.cpp
@@ -168,10 +168,7 @@ Lubby2<DisplacementDim>::integrateStress(
             local_lubby2_properties.update(sig_eff);
         };
 
-        auto newton_solver = NumLib::NewtonRaphson<
-            decltype(linear_solver), LocalJacobianMatrix,
-            decltype(update_jacobian), LocalResidualVector,
-            decltype(update_residual), decltype(update_solution)>(
+        auto newton_solver = NumLib::NewtonRaphson(
             linear_solver, update_jacobian, update_residual, update_solution,
             _nonlinear_solver_parameters);
 
diff --git a/NumLib/NewtonRaphson.h b/NumLib/NewtonRaphson.h
index 9a9affa0e32..b19dde0e860 100644
--- a/NumLib/NewtonRaphson.h
+++ b/NumLib/NewtonRaphson.h
@@ -28,11 +28,36 @@ struct NewtonRaphsonSolverParameters
 /// library.
 /// The current implementation does not update the solution itself, but calls a
 /// function for the solution's update with the current increment.
-template <typename LinearSolver, typename JacobianMatrix,
-          typename JacobianMatrixUpdate, typename ResidualVector,
+template <typename LinearSolver, typename JacobianMatrixUpdate,
           typename ResidualUpdate, typename SolutionUpdate>
 class NewtonRaphson final
 {
+private:
+    template <typename T>
+    struct FirstArgument
+    {
+        static_assert(
+            !std::is_same_v<T, T>,
+            "For the type T in the expression !std::is_same_v<T, T>: Could not "
+            "find overload for FirstArgument<T>. Expected T to be a callable "
+            "with one argument. Could be a missing overload of FirstArgument "
+            "for the type T.");
+    };
+
+    // Specialization for lambdas and functors with const operator()
+    template <typename C, typename R, typename Arg, typename... Args>
+    struct FirstArgument<R (C::*)(Arg, Args...) const>
+    {
+        using type = std::remove_reference_t<Arg>;
+    };
+
+    template <typename F>
+    using FirstArgumentType = typename FirstArgument<
+        decltype(&std::remove_reference_t<F>::operator())>::type;
+
+    using JacobianMatrix = FirstArgumentType<JacobianMatrixUpdate>;
+    using ResidualVector = FirstArgumentType<ResidualUpdate>;
+
 public:
     NewtonRaphson(LinearSolver& linear_solver,
                   JacobianMatrixUpdate jacobian_update,
diff --git a/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h b/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
index 23fc42edd85..2254f8c7d0f 100644
--- a/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
+++ b/ProcessLib/RichardsMechanics/ComputeMicroPorosity.h
@@ -198,13 +198,9 @@ MicroPorosityStateSpace<DisplacementDim> computeMicroPorosity(
     auto const update_solution = [&](ResidualVectorType const& increment)
     { solution += increment; };
 
-    auto newton_solver =
-        NumLib::NewtonRaphson<decltype(linear_solver), JacobianMatrix,
-                              decltype(update_jacobian), ResidualVectorType,
-                              decltype(update_residual),
-                              decltype(update_solution)>(
-            linear_solver, update_jacobian, update_residual, update_solution,
-            micro_porosity_parameters.nonlinear_solver_parameters);
+    auto newton_solver = NumLib::NewtonRaphson(
+        linear_solver, update_jacobian, update_residual, update_solution,
+        micro_porosity_parameters.nonlinear_solver_parameters);
 
     auto const success_iterations = newton_solver.solve(jacobian);
 
diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoMaterialProperties.cpp b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoMaterialProperties.cpp
index 5f26bac1373..8aff9cc716e 100644
--- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoMaterialProperties.cpp
+++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoMaterialProperties.cpp
@@ -136,10 +136,7 @@ bool TwoPhaseFlowWithPrhoMaterialProperties::computeConstitutiveRelation(
         const double residuum_tolerance(1.e-14);
         const double increment_tolerance(0);
 
-        auto newton_solver = NumLib::NewtonRaphson<
-            decltype(linear_solver), LocalJacobianMatrix,
-            decltype(update_jacobian), LocalResidualVector,
-            decltype(update_residual), decltype(update_solution)>(
+        auto newton_solver = NumLib::NewtonRaphson(
             linear_solver, update_jacobian, update_residual, update_solution,
             {maximum_iterations, residuum_tolerance, increment_tolerance});
 
diff --git a/Tests/NumLib/NewtonRaphson.cpp b/Tests/NumLib/NewtonRaphson.cpp
index 3ce9fc8d47c..97603cb2e2b 100644
--- a/Tests/NumLib/NewtonRaphson.cpp
+++ b/Tests/NumLib/NewtonRaphson.cpp
@@ -40,13 +40,9 @@ TEST(NumLibNewtonRaphson, Sqrt3)
     auto const update_solution = [&state](LocalResidualVector const& increment)
     { state += increment[0]; };
 
-    auto const newton_solver =
-        NumLib::NewtonRaphson<decltype(linear_solver), LocalJacobianMatrix,
-                              decltype(update_jacobian), LocalResidualVector,
-                              decltype(update_residual),
-                              decltype(update_solution)>(
-            linear_solver, update_jacobian, update_residual, update_solution,
-            {maximum_iterations, tolerance, 0.0});
+    auto const newton_solver = NumLib::NewtonRaphson(
+        linear_solver, update_jacobian, update_residual, update_solution,
+        {maximum_iterations, tolerance, 0.0});
     auto const success_iterations = newton_solver.solve(jacobian);
 
     EXPECT_TRUE(static_cast<bool>(success_iterations));