Add ability to output the EPS strains as we still need it for some optimization capabilities

Author: rcarson3

src/option_parser.cpp

@@ -517,6 +517,7 @@ void ExaOptions::get_visualizations()
    avg_euler_strain_fname = _avg_euler_strain_fname;
    std::string _avg_pl_work_fname = toml::find_or<std::string>(table, "avg_pl_work_fname", "avg_pl_work.txt");
    avg_pl_work_fname = _avg_pl_work_fname;
+   avg_eps_fname = toml::find_or<std::string>(table, "avg_eps_fname", "avg_eps.txt");
    light_up = toml::find_or<bool>(table, "light_up", false);
    if (light_up) {
 
@@ -761,6 +762,10 @@ void ExaOptions::print_options()
       std::cout << "Additional averages being computed" << std::endl;
       std::cout << "Average deformation gradient filename: " << avg_def_grad_fname << std::endl;
       std::cout << "Average plastic work filename: " << avg_pl_work_fname << std::endl;
+      std::cout << "Average eulerian strain filename: " << avg_euler_strain_fname << std::endl;
+      std::cout << "Average equivalent plastic strain filename: " << avg_eps_fname << std::endl;
+
+
    }
    else
    {

src/option_parser.hpp

@@ -58,6 +58,7 @@ class ExaOptions {
       // average stress file name
       std::string avg_stress_fname;
       std::string avg_pl_work_fname;
+      std::string avg_eps_fname;
       std::string avg_def_grad_fname;
       std::string avg_euler_strain_fname;
       bool additional_avgs;

src/options.toml

@@ -274,6 +274,8 @@ Version = "0.8.0"
     avg_pl_work_fname = "avg_pl_work.txt"
     # Optional - the file name for our average eulerian strain file
     avg_euler_strain_fname = "avg_euler_strain.txt"
+    # Optional - the file name for our average equivalent plastic strain file
+    avg_eps_fname = "avg_eps.txt"
     # Options to drive light_up type calculations in-situ in the code
     light_up = false
     # What HKL planes we want to do the light_up measurements on

src/system_driver.cpp

@@ -116,6 +116,7 @@ SystemDriver::SystemDriver(ParFiniteElementSpace &fes,
      avg_stress_fname(options.avg_stress_fname), avg_pl_work_fname(options.avg_pl_work_fname),
      avg_def_grad_fname(options.avg_def_grad_fname),
      avg_euler_strain_fname(options.avg_euler_strain_fname),
+     avg_eps_fname(options.avg_eps_fname),
      vgrad_origin_flag(options.vgrad_origin_flag), mono_def_flag(options.mono_def_flag),
      def_grad(q_kinVars0), evec(q_evec)
 {
@@ -637,12 +638,11 @@ void SystemDriver::UpdateModel()
       Vector state_var(qstate_var->GetVDim());
       state_var = 0.0;
 
-      std::string s_pl_work = "pl_work";
-      auto qf_mapping = model->GetQFMapping();
-      auto pair = qf_mapping->find(s_pl_work)->second;
-
       exaconstit::kernel::ComputeVolAvgTensor<false>(fes, qstate_var, state_var, state_var.Size(), class_device);
 
+      mfem::Vector history(model->GetMatVars0()->GetVDim());
+      exaconstit::kernel::ComputeVolAvgTensor<true>(fes, model->GetMatVars0(), history, history.Size(), class_device);
+
       std::cout.setf(std::ios::fixed);
       std::cout.setf(std::ios::showpoint);
       std::cout.precision(8);
@@ -651,16 +651,32 @@ void SystemDriver::UpdateModel()
       MPI_Comm_rank(MPI_COMM_WORLD, &my_id);
       // Now we're going to save off the average stress tensor to a file
       if (my_id == 0) {
-         std::ofstream file;
-         file.open(avg_pl_work_fname, std::ios_base::app);
-         file << state_var[pair.first] << std::endl;
+         auto qf_mapping = model->GetQFMapping();
+         {
+            std::string s_pl_work = "pl_work";
+            auto pair = qf_mapping->find(s_pl_work)->second;
+
+            std::ofstream file;
+            file.open(avg_pl_work_fname, std::ios_base::app);
+            file << state_var[pair.first] << std::endl;
+         }
+
+         {
+            std::string s_eps = "shrEff";
+            auto pair = qf_mapping->find(s_eps)->second;
+
+            std::ofstream file;
+            file.open(avg_eps_fname, std::ios_base::app);
+            file << state_var[pair.first] << std::endl;
+         }
+
       }
-      mech_operator->CalculateDeformationGradient(def_grad);
    }
 
    if (additional_avgs)
    {
       CALI_CXX_MARK_SCOPE("extra_avgs_def_grad_computation");
+      mech_operator->CalculateDeformationGradient(def_grad);
       const QuadratureFunction *qstate_var = &def_grad;
       // Here we're getting the average stress value
       Vector dgrad(qstate_var->GetVDim());

src/system_driver.hpp

@@ -66,6 +66,7 @@ class SystemDriver
       std::string avg_pl_work_fname;
       std::string avg_def_grad_fname;
       std::string avg_euler_strain_fname;
+      std::string avg_eps_fname;
       std::string auto_dt_fname;
 
       // define a boundary attribute array and initialize to 0

test/data/voce_ea.toml

@@ -98,6 +98,8 @@ Version = "0.6.0"
     avg_pl_work_fname = "test_voce_ea_pl_work.txt"
     # Optional - the file name for our average eulerian strain file
     avg_euler_strain_fname = "test_voce_ea_euler_strain.txt"
+    # Optional - the file name for our average equivalent plastic strain file
+    avg_eps_fname = "test_voce_ea_eps.txt"
 [Solvers]
     # Option for how our assembly operation is conducted. Possible choices are
     # FULL, PA, EA

test/data/voce_ea_cs.toml

@@ -99,6 +99,7 @@ Version = "0.6.0"
     avg_def_grad_fname  = "test_voce_ea_cs_def_grad.txt"
     avg_pl_work_fname   = "test_voce_ea_cs_pl_work.txt"
     avg_euler_strain_fname = "test_voce_ea_cs_euler_strain.txt"
+    avg_eps_fname = "test_voce_ea_cs_eps.txt"
 [Solvers]
     # Option for how our assembly operation is conducted. Possible choices are
     # FULL, PA, EA

test/data/voce_ea_cs_eps.txt

@@ -0,0 +1,40 @@
+0
+3.46811e-07
+3.88961e-05
+0.000130992
+0.000235568
+0.000342301
+0.000449272
+0.000556043
+0.000662545
+0.000768796
+0.000874843
+0.000980726
+0.00108647
+0.00119208
+0.00129759
+0.001403
+0.00150833
+0.00161357
+0.00171875
+0.00182386
+0.0019289
+0.00213875
+0.00234838
+0.0025578
+0.00276702
+0.00297608
+0.00318499
+0.0036022
+0.00401895
+0.00443526
+0.00485118
+0.00505908
+0.00568185
+0.00620045
+0.00671866
+0.00749508
+0.00827084
+0.00904598
+0.0098206
+0.0108524

test/data/voce_ea_eps.txt

@@ -0,0 +1,40 @@
+0
+3.46778e-07
+3.88802e-05
+0.000130941
+0.000235474
+0.000342153
+0.00044906
+0.000555758
+0.000662176
+0.000768333
+0.000874276
+0.000980044
+0.00108566
+0.00119114
+0.0012965
+0.00140175
+0.00150691
+0.00161198
+0.00171697
+0.00182188
+0.00192672
+0.00213613
+0.00234528
+0.00255418
+0.00276285
+0.00297131
+0.00317957
+0.00359543
+0.00401064
+0.00442527
+0.00483934
+0.00504623
+0.00566586
+0.00618153
+0.00669657
+0.00746785
+0.0082379
+0.00900677
+0.00977455
+0.0107965

test/test_mechanics.py

@@ -132,12 +132,18 @@ def runExtraSystemCommands(params):
     check_stress(ans_pwd, test_pwd, test)
     cmd = 'rm ' + pwd.rstrip() + '/test_'+tresult+'_euler_strain.txt'
     subprocess.run(cmd.rstrip(), stdout=subprocess.PIPE, shell=True)
+    ans_pwd = pwd.rstrip() + '/' + ans[4]
+    tresult = test.split(".")[0]
+    test_pwd = pwd.rstrip() + '/test_'+tresult+'_eps.txt'
+    check_stress(ans_pwd, test_pwd, test)
+    cmd = 'rm ' + pwd.rstrip() + '/test_'+tresult+'_eps.txt'
+    subprocess.run(cmd.rstrip(), stdout=subprocess.PIPE, shell=True)
     return True
 
 def runExtra():
     test_cases = ["voce_ea.toml"]
 
-    test_results = [("voce_ea_stress.txt", "voce_ea_def_grad.txt", "voce_ea_pl_work.txt", "voce_ea_euler_strain.txt")]
+    test_results = [("voce_ea_stress.txt", "voce_ea_def_grad.txt", "voce_ea_pl_work.txt", "voce_ea_euler_strain.txt", "voce_ea_eps.txt")]
 
     result = subprocess.run('pwd', stdout=subprocess.PIPE)
 
@@ -149,7 +155,8 @@ def runExtra():
         cmd = 'rm ' + pwd.rstrip() + '/test_'+tresult+'_stress.txt ' + pwd.rstrip() \
             + '/test_'+tresult+'_pl_work.txt ' + pwd.rstrip() \
             + '/test_'+tresult+'_def_grad.txt' + pwd.rstrip() \
-            + '/test_'+tresult+'_euler_strain.txt' + pwd.rstrip()
+            + '/test_'+tresult+'_euler_strain.txt' + pwd.rstrip() \
+            + '/test_'+tresult+'_eps.txt' + pwd.rstrip()
         result = subprocess.run(cmd.rstrip(), stdout=subprocess.PIPE, shell=True)
 
     params =  zip(test_cases, test_results)