From 436fd6ece9bdbac59ba7dbe09f88fe2d47b515fb Mon Sep 17 00:00:00 2001 From: Mark Messner Date: Mon, 18 Dec 2023 16:01:24 -0600 Subject: [PATCH] New driver for verifying large deformation NEML2 models --- util/neml2_verification/README.md | 30 ++++- ...eml2_test_large_deformation_incremental.py | 118 ++++++++++++++++++ 2 files changed, 146 insertions(+), 2 deletions(-) create mode 100755 util/neml2_verification/run_neml2_test_large_deformation_incremental.py diff --git a/util/neml2_verification/README.md b/util/neml2_verification/README.md index 541bf54f..590232d0 100644 --- a/util/neml2_verification/README.md +++ b/util/neml2_verification/README.md @@ -1,9 +1,11 @@ -This utility runs a model for comparison to NEML2. A complete +These utilities runs a model for comparison to NEML2. A complete test case includes: 1. A NEML XML file defining the model. 2. A plaintext test file with the following format: +*For small deformations* + ``` # These are comment lines # Header: with the following space-separated entries @@ -25,9 +27,33 @@ time strain_xx strain_yy strain_zz strain_yz strain_xz strain_xy stress_xx stres ... ``` +*For large deformations with an incremental interface* +``` +# These are comment lines +# Header: with the following space-separated entries +# 1. NEML XML file name +# 2. NEML model name +# 3. NEML2 serialization file name +# 4. NEML2 model name +# 5. Either "with_temperature" or "no_temperature" +# 6. rtol for comparison +# 7. atol for comparison +viscoplastic.xml viscoplastic_linear_isotropic xx xx no_temperature 1e-5 1e-8 +# A line with a test description that could be read by a test harness +Compare Perzyna viscoplasticity with linear isotropic hardening +# Test input as a space-separated list of time strain stress and +# (optionally) temperature values +time deformation_rate_xx deformation_rate_yy deformation_rate_zz deformation_rate_yz deformation_rate_xz deformation_rate_xy vorticity_zy vorticity_xz vorticity_yx stress_xx stress_yy stress_zz stress_yz stress_xz stress_xy (temperature) + +... +``` + 3. Eventually, a serialized NEML2 model file. -This utility takes the path to a NEML XML file, the name of the model, +The small deformation utility takes the path to a NEML XML file, the name of the model, the strain rate, the maximum strain, the number of steps, and (optionally) an isothermal temperature, runs the model, and writes the test file in the same directory as the XML file. + +The large deformation utility is similar, expect it takes a spatial velocity gradient +tensor, unrolled in row major order, and a maximum time. diff --git a/util/neml2_verification/run_neml2_test_large_deformation_incremental.py b/util/neml2_verification/run_neml2_test_large_deformation_incremental.py new file mode 100755 index 00000000..0487a39d --- /dev/null +++ b/util/neml2_verification/run_neml2_test_large_deformation_incremental.py @@ -0,0 +1,118 @@ +#!/usr/bin/env python3 + +from neml import models, parse, drivers +from neml.math import nemlmath + +import numpy as np +import numpy.linalg as la + +import shutil +import argparse +import os.path + +def drive_all_L(model, L, tmax, nsteps, T): + """ + Strain controlled driver where the user gives a full + Mandel strain vector + """ + d = nemlmath.sym(0.5*(L+L.T)) + w = nemlmath.skew(0.5*(L-L.T)) + times = np.linspace(0, tmax, nsteps+1) + dt = tmax / nsteps + + h_n = model.init_store() + + stress = np.zeros((nsteps+1,6)) + d_n = np.zeros((6,)) + w_n = np.zeros((3,)) + temperature = np.zeros((nsteps+1,)) + temperature[0] = T + + u_n = 0.0 + p_n = 0.0 + t_n = 0.0 + + + for i in range(1,nsteps+1): + d_np1 = d_n + d * dt + w_np1 = w_n + w * dt + t_np1 = t_n + dt + temperature[i] = temperature[i-1] + + s_np1, h_np1, A_np1, B_np1, u_np1, p_np1 = model.update_ld_inc(d_np1, d_n, w_np1, w_n, temperature[i], temperature[i-1], t_np1, t_n, stress[i-1], h_n, u_n, p_n) + + d_n = np.copy(d_np1) + w_n = np.copy(w_np1) + stress[i] = s_np1 + h_n = np.copy(h_np1) + t_n = t_np1 + u_n = u_np1 + p_n = p_np1 + + deformation_rate = np.zeros((nsteps+1,6)) + deformation_rate[:] = d + + vorticity = np.zeros((nsteps+1,3)) + vorticity[:] = w + + return (times, + deformation_rate, + vorticity, + stress, + temperature) + + +if __name__ == "__main__": + parser = argparse.ArgumentParser( + prog = "NEML2 verification test generator", + description = "Generates large deformation verification test files to cross " + "reference against the new NEML2 implementation") + parser.add_argument("file_path") + parser.add_argument("model_name") + parser.add_argument("spatial_velocity_gradient", type = float, nargs = 9) + parser.add_argument("max_time", type = float) + parser.add_argument("nsteps", type = int) + parser.add_argument("test_file") + parser.add_argument("--temperature", type = float, + required = False) + parser.add_argument("--rtol", type = float, + default = 1e-5) + parser.add_argument("--atol", type = float, + default = 1.0e-8) + + args = parser.parse_args() + + basepath = os.path.dirname(args.file_path) + model = parse.parse_xml(args.file_path, args.model_name) + mfile = os.path.basename(args.file_path) + + if args.temperature is not None: + T = args.temperature + use_T = True + T_option = "with_temperature" + else: + T = 20.0 + use_T = False + T_option = "no_temperature" + + # Here + time, D, W, stress, temperature = drive_all_L( + model, np.array(args.spatial_velocity_gradient).reshape(3,3), args.max_time, + args.nsteps, T) + + if use_T: + data = np.hstack((time[:,None],D,W,stress,temperature[:,None])) + else: + data = np.hstack((time[:,None],D,W,stress)) + + with open(os.path.join(basepath, args.test_file), 'w') as f: + f.write("# NEML2/NEML verification test\n") + f.write("%s %s x x %s %e %e\n" % (mfile, args.model_name, T_option, args.rtol, args.atol)) + f.write("# Comment line\n") + f.write("\n") + f.write("# The following uses the Mandel convention to store the strain and stress\n") + data_header = "time deformation_rate_xx deformation_rate_yy deformation_rate_zz deformation_rate_yz deformation_rate_xz deformation_rate_xy vorticity_zy vorticity_xz vorticity_yx stress_xx stress_yy stress_zz stress_yz stress_xz stress_xy" + if use_T: + data_header += " temperature" + f.write(data_header + "\n") + np.savetxt(f, data)