Only compute tangent on demand

src/lib.rs

@@ -1,7 +1,9 @@
 use crate::consts::*;
 use crate::interfaces::*;
 use crate::linear_elasticity::*;
+use crate::mandel::*;
 
 pub mod consts;
 pub mod interfaces;
-pub mod linear_elasticity;
+pub mod linear_elasticity;
+pub mod mandel;

src/linear_elasticity.rs

@@ -1,11 +1,10 @@
-
-use nalgebra::{SVector, SMatrix};
-use crate::interfaces::*;
 use crate::consts::*;
+use crate::interfaces::*;
+use crate::mandel::*;
 use core::ffi::c_double;
+use nalgebra::{SMatrix, SVector};
 use std::collections::HashMap;
 
-
 #[repr(C)]
 struct LinearElasticity3D();
 
@@ -24,20 +23,20 @@ impl ConstitutiveModel<6, 36, 0, 2> for LinearElasticity3D {
         let lambda = parameters[1];
         let del_strain_vec = SVector::<f64, 6>::from_column_slice(del_strain);
         let stress_vec = SVector::<f64, 6>::from_column_slice(stress);
-        let tangent_mat = SYM_ID_6_OUTER_SYM_ID_6 * lambda + (2.0 * mu) * ID_6;
-        let new_stress =  stress_vec + tangent_mat * del_strain_vec;
+        let new_stress =
+            stress_vec + (trace(&del_strain_vec) * lambda) * SYM_ID_6 + (2.0 * mu) * del_strain_vec;
         stress.copy_from_slice(new_stress.as_slice());
         if let Some(tangent) = tangent {
+            let tangent_mat = SYM_ID_6_OUTER_SYM_ID_6 * lambda + (2.0 * mu) * ID_6;
             tangent.copy_from_slice(tangent_mat.as_slice());
         }
     }
-    fn parameters() -> std::collections::HashMap<String, (usize,usize)> {
+    fn parameters() -> std::collections::HashMap<String, (usize, usize)> {
         let mut map = HashMap::new();
-        map.insert("mu".to_string(), (0,1));
-        map.insert("lambda".to_string(), (1,1));
+        map.insert("mu".to_string(), (0, 1));
+        map.insert("lambda".to_string(), (1, 1));
         map
     }
-
 }
 
 pub unsafe fn linear_elasticity3d_fn(
@@ -49,10 +48,12 @@ pub unsafe fn linear_elasticity3d_fn(
     history: *mut c_double,
     parameters: *const c_double,
 ) {
-    let del_strain=unsafe { &*(del_strain as *const [f64;6]) };
-    let stress=unsafe { &mut *(stress as *mut [f64;6]) };
-    let tangent=Some(unsafe { &mut *(tangent as *mut [f64;36]) });
-    let history=unsafe { &mut *(history as *mut [f64;0]) };
-    let parameters=unsafe { &*(parameters as *const [f64;2]) };
-    LinearElasticity3D::evaluate(time, del_time, del_strain, stress, tangent, history, parameters);
-}
+    let del_strain = unsafe { &*(del_strain as *const [f64; 6]) };
+    let stress = unsafe { &mut *(stress as *mut [f64; 6]) };
+    let tangent = Some(unsafe { &mut *(tangent as *mut [f64; 36]) });
+    let history = unsafe { &mut *(history as *mut [f64; 0]) };
+    let parameters = unsafe { &*(parameters as *const [f64; 2]) };
+    LinearElasticity3D::evaluate(
+        time, del_time, del_strain, stress, tangent, history, parameters,
+    );
+}

src/mandel.rs

@@ -0,0 +1,5 @@
+use nalgebra::{SMatrix, SVector, SVectorView};
+
+pub fn trace(vector: &SVector<f64, 6>) -> f64 {
+    vector.x + vector.y + vector.z
+}