diff --git a/MaterialLib/SolidModels/Ehlers-impl.h b/MaterialLib/SolidModels/Ehlers-impl.h
index 441fd72957ad8483795cb543863252375686fa06..882d7a8fc750ea45556a6641a0a906b71c35b44a 100644
--- a/MaterialLib/SolidModels/Ehlers-impl.h
+++ b/MaterialLib/SolidModels/Ehlers-impl.h
@@ -36,7 +36,7 @@
#include <boost/math/special_functions/pow.hpp>
#include <logog/include/logog.hpp>
#include "MaterialLib/SolidModels/KelvinVector.h"
-#include "NewtonRaphson.h"
+#include "NumLib/NewtonRaphson.h"
namespace MaterialLib
{
@@ -82,7 +82,7 @@ struct OnePlusGamma_pTheta final
double const m_p)
: value{1 + gamma_p * theta},
pow_m_p{std::pow(value, m_p)},
- pow_m_p1{pow_m_p/value}
+ pow_m_p1{pow_m_p / value}
{
}
@@ -104,16 +104,14 @@ double plasticFlowVolumetricPart(
}
template <int DisplacementDim>
-typename SolidEhlers<DisplacementDim>::KelvinVector
-plasticFlowDeviatoricPart(
+typename SolidEhlers<DisplacementDim>::KelvinVector plasticFlowDeviatoricPart(
PhysicalStressWithInvariants<DisplacementDim> const& s,
OnePlusGamma_pTheta const& one_gt, double const sqrtPhi,
- typename SolidEhlers<DisplacementDim>::KelvinVector const&
- dtheta_dsigma,
+ typename SolidEhlers<DisplacementDim>::KelvinVector const& dtheta_dsigma,
double const gamma_p, double const m_p)
{
- return (one_gt.pow_m_p * (s.D +
- s.J_2 * m_p * gamma_p * dtheta_dsigma / one_gt.value)) /
+ return (one_gt.pow_m_p *
+ (s.D + s.J_2 * m_p * gamma_p * dtheta_dsigma / one_gt.value)) /
(2 * sqrtPhi);
}
template <int DisplacementDim>
@@ -194,8 +192,8 @@ void calculatePlasticResidual(
s, one_gt, sqrtPhi, dtheta_dsigma, gamma_p, m_p);
KelvinVector const lambda_flow_D = lambda * flow_D;
- residual.template segment<KelvinVectorSize>(KelvinVectorSize)
- .noalias() = eps_p_D_dot - lambda_flow_D;
+ residual.template segment<KelvinVectorSize>(KelvinVectorSize).noalias() =
+ eps_p_D_dot - lambda_flow_D;
// plastic volume strain
{
@@ -494,8 +492,7 @@ void SolidEhlers<DisplacementDim>::MaterialProperties::
template <int DisplacementDim>
typename SolidEhlers<DisplacementDim>::KelvinVector predict_sigma(
double const G, double const K,
- typename SolidEhlers<DisplacementDim>::KelvinVector const&
- sigma_prev,
+ typename SolidEhlers<DisplacementDim>::KelvinVector const& sigma_prev,
typename SolidEhlers<DisplacementDim>::KelvinVector const& eps,
typename SolidEhlers<DisplacementDim>::KelvinVector const& eps_prev,
double const eps_V)
@@ -512,8 +509,8 @@ typename SolidEhlers<DisplacementDim>::KelvinVector predict_sigma(
// dimensioness hydrostatic stress increment
double const pressure = pressure_prev - K / G * (eps_V - e_prev);
// dimensionless deviatoric initial stress
- typename SolidEhlers<DisplacementDim>::KelvinVector const
- sigma_D_prev = P_dev * sigma_prev / G;
+ typename SolidEhlers<DisplacementDim>::KelvinVector const sigma_D_prev =
+ P_dev * sigma_prev / G;
// dimensionless deviatoric stress
typename SolidEhlers<DisplacementDim>::KelvinVector const sigma_D =
sigma_D_prev + 2 * P_dev * (eps - eps_prev);
@@ -599,9 +596,8 @@ bool SolidEhlers<DisplacementDim>::computeConstitutiveRelation(
auto const update_solution = [&](
ResidualVectorType const& increment) {
- sigma.noalias() +=
- increment.template segment<KelvinVectorSize>(
- KelvinVectorSize * 0);
+ sigma.noalias() += increment.template segment<KelvinVectorSize>(
+ KelvinVectorSize * 0);
s = PhysicalStressWithInvariants<DisplacementDim>{G * sigma};
_state.eps_p_D.noalias() +=
increment.template segment<KelvinVectorSize>(
@@ -618,13 +614,12 @@ bool SolidEhlers<DisplacementDim>::computeConstitutiveRelation(
int const maximum_iterations(100);
double const tolerance(1e-14);
- auto newton_solver =
- NewtonRaphson<decltype(linear_solver), JacobianMatrix,
- decltype(update_jacobian), ResidualVectorType,
- decltype(update_residual),
- decltype(update_solution)>(
- linear_solver, update_jacobian, update_residual,
- update_solution, maximum_iterations, tolerance);
+ 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, maximum_iterations, tolerance);
auto const success_iterations = newton_solver.solve(jacobian);
diff --git a/MaterialLib/SolidModels/Lubby2-impl.h b/MaterialLib/SolidModels/Lubby2-impl.h
index 909029fd66c591b202cbad97f6a0dd75e00600ec..703fee922b3d6b9f18b9c67342b188d6cb5f8ba4 100644
--- a/MaterialLib/SolidModels/Lubby2-impl.h
+++ b/MaterialLib/SolidModels/Lubby2-impl.h
@@ -10,7 +10,7 @@
#ifndef MATERIALLIB_SOLIDMODELS_LUBBY2_IMPL_H_
#define MATERIALLIB_SOLIDMODELS_LUBBY2_IMPL_H_
-#include "NewtonRaphson.h"
+#include "NumLib/NewtonRaphson.h"
namespace MaterialLib
{
@@ -105,12 +105,12 @@ bool Lubby2<DisplacementDim>::computeConstitutiveRelation(
const int maximum_iterations(20);
const double tolerance(1.e-10);
- auto newton_solver =
- 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);
+ auto 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);
auto const success_iterations = newton_solver.solve(K_loc);
diff --git a/MaterialLib/SolidModels/NewtonRaphson.h b/NumLib/NewtonRaphson.h
similarity index 92%
rename from MaterialLib/SolidModels/NewtonRaphson.h
rename to NumLib/NewtonRaphson.h
index f59366340f657ed33948883b69cad9474eb21aba..531e769d7254e7aeb5248554191caa5d65284c66 100644
--- a/MaterialLib/SolidModels/NewtonRaphson.h
+++ b/NumLib/NewtonRaphson.h
@@ -7,17 +7,15 @@
*
*/
-#ifndef MATERIALLIB_SOLIDMODELS_NEWTONRAPHSON_H_
-#define MATERIALLIB_SOLIDMODELS_NEWTONRAPHSON_H_
+#ifndef NUMLIB_NEWTONRAPHSON_H_
+#define NUMLIB_NEWTONRAPHSON_H_
#include <boost/optional.hpp>
#include <logog/include/logog.hpp>
#include <Eigen/Dense>
-namespace MaterialLib
-{
-namespace Solids
+namespace NumLib
{
/// Newton-Raphson solver for system of equations using an Eigen linear solvers
/// library.
@@ -92,8 +90,6 @@ private:
const int _maximum_iterations;
const double _tolerance_squared;
};
+} // namespace NumLib
-} // namespace Solids
-} // namespace MaterialLib
-
-#endif // MATERIALLIB_SOLIDMODELS_NEWTONRAPHSON_H_
+#endif // NUMLIB_NEWTONRAPHSON_H_
diff --git a/Tests/MaterialLib/NewtonRaphson.cpp b/Tests/NumLib/NewtonRaphson.cpp
similarity index 91%
rename from Tests/MaterialLib/NewtonRaphson.cpp
rename to Tests/NumLib/NewtonRaphson.cpp
index 5a5e429abb75eff5697b414c72bd830de4f57e77..46412fc2b090abac0f7e36d2c28694314a042702 100644
--- a/Tests/MaterialLib/NewtonRaphson.cpp
+++ b/Tests/NumLib/NewtonRaphson.cpp
@@ -9,9 +9,8 @@
#include <gtest/gtest.h>
#include <limits>
-#include "MaterialLib/SolidModels/NewtonRaphson.h"
-
-TEST(MaterialLibNewtonRaphson, Sqrt3)
+#include "NumLib/NewtonRaphson.h"
+TEST(NumLibNewtonRaphson, Sqrt3)
{
static const int N = 1; // Problem's size.
@@ -40,7 +39,7 @@ TEST(MaterialLibNewtonRaphson, Sqrt3)
auto const update_solution = [&state](
LocalResidualVector const& increment) { state += increment[0]; };
- auto const newton_solver = MaterialLib::Solids::NewtonRaphson<
+ auto const newton_solver = NumLib::NewtonRaphson<
decltype(linear_solver), LocalJacobianMatrix, decltype(update_jacobian),
LocalResidualVector, decltype(update_residual),
decltype(update_solution)>(linear_solver, update_jacobian,