material_verification_tools.hpp

Go to the documentation of this file.

// Copyright (c) Lawrence Livermore National Security, LLC and
// other Serac Project Developers. See the top-level LICENSE file for
// details.
//
// SPDX-License-Identifier: (BSD-3-Clause)
 
#pragma once
 
#include <functional>
 
#include "serac/numerics/functional/tensor.hpp"
#include "serac/numerics/functional/tuple.hpp"
 
namespace serac {
 
template <typename MaterialType, typename StateType, typename... parameter_types>
auto uniaxial_stress_test(double t_max, size_t num_steps, const MaterialType material, const StateType initial_state,
                          std::function<double(double)> epsilon_xx, const parameter_types... parameter_functions)
{
  double t = 0;
  const double dt = t_max / double(num_steps - 1);
 
  auto state = initial_state;
 
  auto sigma_yy_and_zz = [&](auto x) {
    auto epsilon_yy = x[0];
    auto epsilon_zz = x[1];
    using T = decltype(epsilon_yy);
    tensor<T, 3, 3> du_dx{};
    du_dx[0][0] = epsilon_xx(t);
    du_dx[1][1] = epsilon_yy;
    du_dx[2][2] = epsilon_zz;
    auto state_copy = state;
    auto stress = material(state_copy, du_dx, parameter_functions(t)...);
    return tensor{{stress[1][1], stress[2][2]}};
  };
 
  std::vector<tuple<double, tensor<double, 3, 3>, tensor<double, 3, 3>, StateType> > output_history;
  output_history.reserve(num_steps);
 
  tensor<double, 3, 3> dudx{};
  for (size_t i = 0; i < num_steps; i++) {
    auto initial_guess = tensor<double, 2>{dudx[1][1], dudx[2][2]};
    auto epsilon_yy_and_zz = find_root(sigma_yy_and_zz, initial_guess);
    dudx[0][0] = epsilon_xx(t);
    dudx[1][1] = epsilon_yy_and_zz[0];
    dudx[2][2] = epsilon_yy_and_zz[1];
 
    auto stress = material(state, dudx, parameter_functions(t)...);
    output_history.push_back(tuple{t, dudx, stress, state});
 
    t += dt;
  }
 
  return output_history;
}
 
template <typename MaterialType, typename StateType, typename... parameter_types>
auto uniaxial_stress_test_rate_dependent(double t_max, size_t num_steps, const MaterialType material,
                                         const StateType initial_state, std::function<double(double)> epsilon_xx,
                                         const parameter_types... parameter_functions)
{
  const double dt = t_max / double(num_steps - 1);
  auto mat_with_dt = [&material, dt](auto& state, auto du_dx, parameter_types... parameters) {
    return material(state, dt, du_dx, parameters...);
  };
  return uniaxial_stress_test(t_max, num_steps, mat_with_dt, initial_state, epsilon_xx, parameter_functions...);
}
 
template <typename MaterialType, typename StateType, typename... functions>
auto single_quadrature_point_test(double t_max, size_t num_steps, const MaterialType material,
                                  const StateType initial_state, const functions... f)
{
  double t = 0;
  const double dt = t_max / double(num_steps - 1);
  auto state = initial_state;
 
  using output_type = decltype(std::tuple{t, state, f(0.0)..., decltype(material(state, dt, f(0.0)...)){}});
  std::vector<output_type> history;
  history.reserve(num_steps);
 
  for (size_t i = 0; i < num_steps; i++) {
    auto material_output = material(state, dt, f(t)...);
    history.push_back(std::tuple{t, state, f(t)..., material_output});
    t += dt;
  }
 
  return history;
}
 
}  // namespace serac