From 3c3aea838f0bb84d283eeee9918ebb7002ab6e95 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <github@naumov.de>
Date: Thu, 15 Sep 2016 00:42:34 +0000
Subject: [PATCH] [MatL]computeConstitutiveRelation returns boolean.
On error the function returs false to indicate that the local iterations
did not succeed.
---
.../SolidModels/LinearElasticIsotropic.h | 3 +-
MaterialLib/SolidModels/Lubby2-impl.h | 6 ++--
MaterialLib/SolidModels/Lubby2.h | 2 +-
MaterialLib/SolidModels/MechanicsBase.h | 28 +++++++++++--------
.../SmallDeformation/SmallDeformationFEM.h | 7 +++--
5 files changed, 28 insertions(+), 18 deletions(-)
diff --git a/MaterialLib/SolidModels/LinearElasticIsotropic.h b/MaterialLib/SolidModels/LinearElasticIsotropic.h
index 063f1f07433..bf2d729c2fc 100644
--- a/MaterialLib/SolidModels/LinearElasticIsotropic.h
+++ b/MaterialLib/SolidModels/LinearElasticIsotropic.h
@@ -79,7 +79,7 @@ public:
{
}
- void computeConstitutiveRelation(
+ bool computeConstitutiveRelation(
double const t,
ProcessLib::SpatialPosition const& x,
double const /*dt*/,
@@ -97,6 +97,7 @@ public:
C.noalias() += 2 * _mp.mu(t, x) * KelvinMatrix::Identity();
sigma.noalias() = sigma_prev + C * (eps - eps_prev);
+ return true;
}
private:
diff --git a/MaterialLib/SolidModels/Lubby2-impl.h b/MaterialLib/SolidModels/Lubby2-impl.h
index 392fd020ede..909029fd66c 100644
--- a/MaterialLib/SolidModels/Lubby2-impl.h
+++ b/MaterialLib/SolidModels/Lubby2-impl.h
@@ -17,7 +17,7 @@ namespace MaterialLib
namespace Solids
{
template <int DisplacementDim>
-void Lubby2<DisplacementDim>::computeConstitutiveRelation(
+bool Lubby2<DisplacementDim>::computeConstitutiveRelation(
double const t,
ProcessLib::SpatialPosition const& x,
double const dt,
@@ -115,7 +115,7 @@ void Lubby2<DisplacementDim>::computeConstitutiveRelation(
auto const success_iterations = newton_solver.solve(K_loc);
if (!success_iterations)
- OGS_FATAL("The local Newton method did not succeed.");
+ return false;
// If the Newton loop didn't run, the linear solver will not be
// initialized.
@@ -142,6 +142,8 @@ void Lubby2<DisplacementDim>::computeConstitutiveRelation(
auto const& P_sph = Invariants::spherical_projection;
C.noalias() = state.GM * dzdE * P_dev + 3. * state.KM * P_sph;
+
+ return true;
}
template <int DisplacementDim>
diff --git a/MaterialLib/SolidModels/Lubby2.h b/MaterialLib/SolidModels/Lubby2.h
index e2d8d50672e..2d41b6fc954 100644
--- a/MaterialLib/SolidModels/Lubby2.h
+++ b/MaterialLib/SolidModels/Lubby2.h
@@ -127,7 +127,7 @@ public:
{
}
- void computeConstitutiveRelation(
+ bool computeConstitutiveRelation(
double const t,
ProcessLib::SpatialPosition const& x_position,
double const dt,
diff --git a/MaterialLib/SolidModels/MechanicsBase.h b/MaterialLib/SolidModels/MechanicsBase.h
index 8875291d6fc..7040c53aed8 100644
--- a/MaterialLib/SolidModels/MechanicsBase.h
+++ b/MaterialLib/SolidModels/MechanicsBase.h
@@ -53,7 +53,9 @@ struct MechanicsBase
/// Dynamic size Kelvin vector and matrix wrapper for the polymorphic
/// constitutive relation compute function.
- void computeConstitutiveRelation(
+ /// Returns false in case of errors in the computation if Newton iterations
+ /// did not converge, for example.
+ bool computeConstitutiveRelation(
double const t,
ProcessLib::SpatialPosition const& x,
double const dt,
@@ -76,25 +78,29 @@ struct MechanicsBase
KelvinVector sigma_{sigma};
KelvinMatrix C_{C};
- computeConstitutiveRelation(t,
- x,
- dt,
- eps_prev_,
- eps_,
- sigma_prev_,
- sigma_,
- C_,
- material_state_variables);
+ bool const result =
+ computeConstitutiveRelation(t,
+ x,
+ dt,
+ eps_prev_,
+ eps_,
+ sigma_prev_,
+ sigma_,
+ C_,
+ material_state_variables);
sigma = sigma_;
C = C_;
+ return result;
}
/// Computation of the constitutive relation for specific material model.
/// This should be implemented in the derived model. Fixed Kelvin vector and
/// matrix size version; for dynamic size arguments there is an overloaded
/// wrapper function.
- virtual void computeConstitutiveRelation(
+ /// Returns false in case of errors in the computation if Newton iterations
+ /// did not converge, for example.
+ virtual bool computeConstitutiveRelation(
double const t,
ProcessLib::SpatialPosition const& x,
double const dt,
diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index eed4503422d..01cb3ab818c 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
@@ -234,9 +234,10 @@ public:
Eigen::Map<typename BMatricesType::NodalForceVectorType const>(
local_x.data(), ShapeFunction::NPOINTS * DisplacementDim);
- _ip_data[ip]._solid_material.computeConstitutiveRelation(
- t, x_position, _process_data.dt, eps_prev, eps, sigma_prev,
- sigma, C, material_state_variables);
+ if (!_ip_data[ip]._solid_material.computeConstitutiveRelation(
+ t, x_position, _process_data.dt, eps_prev, eps, sigma_prev,
+ sigma, C, material_state_variables))
+ OGS_FATAL("Computation of local constitutive relation failed.");
local_b.noalias() -= B.transpose() * sigma * detJ * wp.getWeight();
local_Jac.noalias() +=
--
GitLab