From 06ebd074a4ccbffb6a28390193b2d15963fc4c1a Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Tue, 25 Jun 2024 15:34:59 +0200 Subject: [PATCH] [NL] Remove Matrix/Vector types from NR interface The types for the jacobian matrix and the residual vector are not needed in the NewtonRaphson class interface, only in the interior and can be (mostly) deduced from the JacobianMatrixUpdate and ResidualUpdate first function argument. --- MaterialLib/FractureModels/Coulomb.cpp | 10 ++----- MaterialLib/SolidModels/CreepBGRa.cpp | 8 ++--- MaterialLib/SolidModels/Ehlers.cpp | 5 +--- MaterialLib/SolidModels/Lubby2.cpp | 5 +--- NumLib/NewtonRaphson.h | 29 +++++++++++++++++-- .../RichardsMechanics/ComputeMicroPorosity.h | 10 ++----- ...TwoPhaseFlowWithPrhoMaterialProperties.cpp | 5 +--- Tests/NumLib/NewtonRaphson.cpp | 10 ++----- 8 files changed, 41 insertions(+), 41 deletions(-) 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::computeConstitutiveRelation( sigma.noalias() += sigma_prev; }; - auto newton_solver = - NumLib::NewtonRaphson( - 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::integrateStress( { solution += increment; }; auto newton_solver = - NumLib::NewtonRaphson( - 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::integrateStress( mp.G * solution.template segment(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::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 class NewtonRaphson final { +private: + template + struct FirstArgument + { + static_assert( + !std::is_same_v, + "For the type T in the expression !std::is_same_v: Could not " + "find overload for FirstArgument. 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 + struct FirstArgument + { + using type = std::remove_reference_t; + }; + + template + using FirstArgumentType = typename FirstArgument< + decltype(&std::remove_reference_t::operator())>::type; + + using JacobianMatrix = FirstArgumentType; + using ResidualVector = FirstArgumentType; + 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 computeMicroPorosity( auto const update_solution = [&](ResidualVectorType const& increment) { solution += increment; }; - auto newton_solver = - NumLib::NewtonRaphson( - 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( - 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(success_iterations));