Support crystal plasticity simulations with JAX-CPFEM

matflow/data/scripts/jax_cpfem/parse_stdout.py

@@ -0,0 +1,49 @@
+import os
+import re
+from pathlib import Path
+
+import numpy as np
+
+_PAT_FP = (
+    r"Fp\s*=\s*\[\[\s*([-+]?\d*\.\d*\s+[-+]?\d*\.\d*\s+[-+]?\d*\.\d*)\s*\]\s*"
+    r"\[\s*([-+]?\d*\.\d*\s+[-+]?\d*\.\d*\s+[-+]?\d*\.\d*)\s*\]\s*"
+    r"\[\s*([-+]?\d*\.\d*\s+[-+]?\d*\.\d*\s+[-+]?\d*\.\d*)\s*\]\s*\]"
+)
+_PAT_key_data_list = R"{key}:\s*\[\s*((?:\d+\.\d+\s*)+)\]"
+_PAT_STRESS_XX = _PAT_key_data_list.format(key="stress_xx")
+_PAT_STRESS_YY = _PAT_key_data_list.format(key="stress_yy")
+_PAT_STRESS_ZZ = _PAT_key_data_list.format(key="stress_zz")
+_PAT_STRESS_VM = _PAT_key_data_list.format(key="von_mises_stress")
+
+
+def parse_stdout(jax_cpfem_stdout: Path):
+
+    encoding = "utf-16" if os.name == "nt" else "utf-8"  # docker on windows nonsense
+    with jax_cpfem_stdout.open("rt", encoding=encoding) as fh:
+        contents = fh.read()
+
+    fp = []
+    for fp_res in re.findall(_PAT_FP, contents):
+        fp.append([[float(i) for i in fp_row.split()] for fp_row in fp_res])
+    fp_arr = np.asarray(fp)
+
+    stress_xx = np.array(
+        [float(i) for i in re.search(_PAT_STRESS_XX, contents).groups()[0].split()]
+    )
+    stress_yy = np.array(
+        [float(i) for i in re.search(_PAT_STRESS_YY, contents).groups()[0].split()]
+    )
+    stress_zz = np.array(
+        [float(i) for i in re.search(_PAT_STRESS_ZZ, contents).groups()[0].split()]
+    )
+    stress_von_mises = np.array(
+        [float(i) for i in re.search(_PAT_STRESS_VM, contents).groups()[0].split()]
+    )
+
+    return {
+        "Fp": fp_arr,
+        "stress_xx": stress_xx,
+        "stress_yy": stress_yy,
+        "stress_zz": stress_zz,
+        "stress_von_mises": stress_von_mises,
+    }

matflow/data/scripts/jax_cpfem/write_mesh.py

@@ -0,0 +1,14 @@
+import numpy as np
+
+
+def write_mesh(path, volume_element):
+    np.savez(
+        path,
+        num_cells=volume_element["grid_size"][:],  # i.e. Nx, Ny, Nz
+        domain=np.asarray(volume_element["size"]),  # i.e. domain_x, domain_y, domain_z
+        cell_grain_indices=volume_element["element_material_idx"][:].reshape(
+            -1
+        ),  # TODO: check order!
+        quaternions=volume_element["orientations"]["quaternions"][:],
+        grain_orientation_indices=volume_element["constituent_orientation_idx"][:],
+    )