gradient model in terms of total nonlocal strain

srosenbu · Feb 14, 2024 · a15bf25 · a15bf25
1 parent 1fdf069
commit a15bf25
Show file tree

Hide file tree

Showing 3 changed files with 20 additions and 13 deletions.
diff --git a/python/comfe/cdm.py b/python/comfe/cdm.py
@@ -361,7 +361,7 @@ class CDMNonlocal(NonlocalInterface):
     M: df.fem.Function
     # fields: dict[str, df.fem.Function]
     # q_fields: dict[str, df.fem.Function | ufl.core.expr.Expr]
-    rate_evaluator: QuadratureEvaluator
+    strain_evaluator: QuadratureEvaluator
     parameters: dict[str, float]
     form: df.fem.FormMetaClass
     t: float
@@ -379,17 +379,16 @@ def __init__(
         parameters: dict[str, float],
         quadrature_rule: QuadratureRule,
         q_fields_local: dict[str, df.fem.Function],
-        q_field_nonlocal_rate: df.fem.Function,
+        q_field_nonlocal: df.fem.Function,
         Q_local_damage: str | None = None,
     ):
         q_fields = {Q_local: q_fields_local[Q_local]}
 
         if Q_nonlocal_rate[-4:] != "rate":
             raise ValueError("Q_nonlocal_rate must be some rate for CDM, you provided " + Q_nonlocal_rate)
 
-        q_fields[Q_nonlocal_rate] = q_field_nonlocal_rate
-
         Q_nonlocal = Q_nonlocal_rate[:-5]
+        q_fields[Q_nonlocal] = q_field_nonlocal
 
         fields = {
             Q_nonlocal_rate: df.fem.Function(function_space, name=Q_nonlocal_rate),
@@ -419,7 +418,7 @@ def __init__(
 
         f_form = df.fem.form(f_ufl)
 
-        rate_evaluator = QuadratureEvaluator(fields[Q_nonlocal_rate], function_space.mesh, quadrature_rule)
+        strain_evaluator = QuadratureEvaluator(fields[Q_nonlocal], function_space.mesh, quadrature_rule)
         super().__init__(
             Q_local=Q_local,
             Q_nonlocal=Q_nonlocal,
@@ -431,7 +430,7 @@ def __init__(
             q_fields=q_fields,
             fields=fields,
             form=f_form,
-            rate_evaluator=rate_evaluator,
+            strain_evaluator=strain_evaluator,
         )
 
     def step(self, h: float) -> None:
@@ -454,8 +453,8 @@ def step(self, h: float) -> None:
         self.fields[self.Q_nonlocal].vector.array[:] += h * self.fields[self.Q_nonlocal_rate].vector.array
         self.fields[self.Q_nonlocal].x.scatter_forward()
 
-        self.rate_evaluator(self.q_fields[self.Q_nonlocal_rate])
-        self.q_fields[self.Q_nonlocal_rate].x.scatter_forward()
+        self.strain_evaluator(self.q_fields[self.Q_nonlocal])
+        self.q_fields[self.Q_nonlocal].x.scatter_forward()
 
         self.t += h
 

diff --git a/src/gradient_jh2.rs b/src/gradient_jh2.rs
@@ -51,7 +51,12 @@ impl ConstitutiveModel for GradientJH23D {
         d_eps_vol *= -1.0;
 
         let sigma_0 = input.get_vector::<{ Q::MandelStress.size() }>(Q::MandelStress, ip);
-        let del_lambda_nonlocal = del_t * input.get_scalar(Q::EqNonlocalPlasticStrainRate, ip).max(0.0);
+        //let del_lambda_nonlocal = del_t * input.get_scalar(Q::EqNonlocalPlasticStrainRate, ip).max(0.0);
+        let del_lambda_nonlocal = (input.get_scalar(Q::HistoryMaximum, ip) - input.get_scalar(Q::EqNonlocalPlasticStrain, ip)).max(0.0);
+        if del_lambda_nonlocal > 0.0 {
+            output.set_scalar(Q::HistoryMaximum, ip, input.get_scalar(Q::EqNonlocalPlasticStrain, ip));
+        }
+
         let mut del_lambda = 0.0;
 
 
@@ -74,9 +79,7 @@ impl ConstitutiveModel for GradientJH23D {
         let damage_0 = input.get_scalar(Q::Damage, ip);
         let mut damage_1 = 0.0;
         if self.parameters.E_F > 0.0 {
-            let lambda_old = - self.parameters.E_F * (1. - damage_0).ln();
-            let lambda_new = lambda_old + del_lambda_nonlocal;
-            damage_1 = 1. - (-lambda_new / self.parameters.E_F).exp();
+            damage_1 = 1. - (-input.get_scalar(Q::EqNonlocalPlasticStrain, ip) / self.parameters.E_F).exp();
         } else {
             damage_1 = (damage_0 + del_lambda_nonlocal / e_p_f).min(self.parameters.DMAX);
         }
@@ -229,7 +232,8 @@ impl ConstitutiveModel for GradientJH23D {
     fn define_input(&self) -> HashMap<Q, QDim> {
         HashMap::from([
             (Q::VelocityGradient, QDim::SquareTensor(3)),
-            (Q::EqNonlocalPlasticStrainRate, QDim::Scalar),
+            //(Q::EqNonlocalPlasticStrainRate, QDim::Scalar),
+            (Q::EqNonlocalPlasticStrain, QDim::Scalar),
             (Q::CellDiameter, QDim::Scalar),
         ])
     }
@@ -244,6 +248,7 @@ impl ConstitutiveModel for GradientJH23D {
             (Q::Density, QDim::Scalar),
             (Q::EqPlasticStrain, QDim::Scalar),
             (Q::BulkViscosity, QDim::Scalar),
+            (Q::HistoryMaximum, QDim::Scalar),
         ])
     }
 

diff --git a/src/interfaces.rs b/src/interfaces.rs
@@ -63,6 +63,8 @@ pub enum Q {
     BulkViscosity,
     #[strum(serialize = "CellDiameter", serialize = "cell_diameter")]
     CellDiameter,
+    #[strum(serialize = "HistoryMaximum", serialize = "history_maximum")]
+    HistoryMaximum,
     #[strum(serialize = "_LAST", serialize = "_last")]
     _LAST,
 }
@@ -234,6 +236,7 @@ impl Q {
             Q::InternalPlasticEnergyRate => QDim::Scalar,
             Q::BulkViscosity => QDim::Scalar,
             Q::CellDiameter => QDim::Scalar,
+            Q::HistoryMaximum => QDim::Scalar,
             Q::_LAST => QDim::Scalar,
         }
     }