LAROMANCEStressUpdateBase.h

Go to the documentation of this file.

//* This file is part of the MOOSE framework
//* https://mooseframework.inl.gov
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#pragma once

#include "RadialReturnCreepStressUpdateBase.h"
#include "nlohmann/json.h"

template <bool is_ad>
class LAROMANCEStressUpdateBaseTempl : public RadialReturnCreepStressUpdateBaseTempl<is_ad>
{

public:
  static InputParameters validParams();

  LAROMANCEStressUpdateBaseTempl(const InputParameters & parameters);

  virtual void resetIncrementalMaterialProperties() override;
  virtual void
  storeIncrementalMaterialProperties(const unsigned int total_number_substeps) override;

protected:
  virtual void exportJSON();

  virtual bool substeppingCapabilityEnabled() override;

  enum class ROMInputTransform
  {
    LINEAR,
    LOG,
    EXP
  };

  virtual void initialSetup() override;

  // Setup unit conversion factors. The required units in the ROM are:
  // Cell dislocation density: m^-2
  // Wall dislocation density: m^-2
  // MX phase fracture: nondim.
  // stress: MPa
  // strain: nondim.
  // temperature: K
  virtual void setupUnitConversionFactors(const InputParameters & parameters);

  virtual void initQpStatefulProperties() override;
  virtual void
  computeStressInitialize(const GenericReal<is_ad> & effective_trial_stress,
                          const GenericRankFourTensor<is_ad> & elasticity_tensor) override;
  virtual GenericReal<is_ad> computeResidual(const GenericReal<is_ad> & effective_trial_stress,
                                             const GenericReal<is_ad> & scalar) override;

  virtual GenericReal<is_ad>
  computeDerivative(const GenericReal<is_ad> & /*effective_trial_stress*/,
                    const GenericReal<is_ad> & /*scalar*/) override
  {
    return _derivative;
  }

  virtual void
  computeStressFinalize(const GenericRankTwoTensor<is_ad> & plastic_strain_increment) override;
  virtual GenericReal<is_ad>
  maximumPermissibleValue(const GenericReal<is_ad> & effective_trial_stress) const override;
  virtual Real computeTimeStepLimit() override;

  void outputIterationSummary(std::stringstream * iter_output,
                              const unsigned int total_it) override;

  virtual void outputIterationStep(std::stringstream * iter_output,
                                   const GenericReal<is_ad> & effective_trial_stress,
                                   const GenericReal<is_ad> & scalar,
                                   const Real reference_residual) override;

  enum class WindowFailure
  {
    ERROR,
    EXCEPTION,
    WARN,
    IGNORE,
    DONOTHING,
    USELIMIT,
    EXTRAPOLATE
  };

  void precomputeROM(const unsigned out_index);

  GenericReal<is_ad> computeROM(const unsigned int tile,
                                const unsigned int partition,
                                const unsigned out_index,
                                const bool derivative = false);

  bool checkInTile(const unsigned int p, const unsigned int t) const;

  bool areTilesNotIdentical(const unsigned int p,
                            const unsigned int t,
                            const unsigned int tt,
                            const unsigned int in_index);

  GenericReal<is_ad> normalizeInput(const GenericReal<is_ad> & input,
                                    const ROMInputTransform transform,
                                    const Real transform_coef,
                                    const std::vector<Real> & transformed_limits,
                                    const bool derivative = false);

  void buildPolynomials(const unsigned int p,
                        const GenericReal<is_ad> & rom_input,
                        std::vector<GenericReal<is_ad>> & polynomial_inputs,
                        const GenericReal<is_ad> & drom_input = 0,
                        const bool derivative = false);

  void precomputeValues(const unsigned int p,
                        const std::vector<Real> & coefs,
                        const std::vector<std::vector<GenericReal<is_ad>>> & polynomial_inputs,
                        std::vector<GenericReal<is_ad>> & precomputed);

  GenericReal<is_ad>
  computeValues(const unsigned int p,
                const std::vector<GenericReal<is_ad>> & precomputed,
                const std::vector<std::vector<GenericReal<is_ad>>> & polynomial_inputs,
                const std::vector<GenericReal<is_ad>> & dpolynomial_inputs = {},
                const bool derivative = false);

  virtual GenericReal<is_ad> convertOutput(const unsigned int p,
                                           const std::vector<Real> & old_input_values,
                                           const GenericReal<is_ad> & rom_output,
                                           const unsigned out_index,
                                           const GenericReal<is_ad> & drom_output = 0.0,
                                           const bool derivative = false);

  GenericReal<is_ad> computePolynomial(const GenericReal<is_ad> & value,
                                       const unsigned int degree,
                                       const bool derivative = false);

  GenericReal<is_ad> sigmoid(const Real lower,
                             const Real upper,
                             const GenericReal<is_ad> & val,
                             const bool derivative = false);

  void computeTileWeight(std::vector<std::vector<GenericReal<is_ad>>> & weights,
                         GenericReal<is_ad> & input,
                         const unsigned int in_index,
                         const bool derivative = false);

  virtual void computePartitionWeights(std::vector<GenericReal<is_ad>> & weights,
                                       std::vector<GenericReal<is_ad>> & dweights_dstress)
  {
    if (_num_partitions != 1)
      mooseError("Internal error: If number of partitions is not one, then computePartitionWeights "
                 "must be defined");
    weights[0] = 1.0;
    dweights_dstress[0] = 0.0;
  }

  template <typename T>
  void convertValue(T & x,
                    const ROMInputTransform transform,
                    const Real coef,
                    const bool derivative = false)
  {
    using std::exp, std::log;

    if (transform == ROMInputTransform::EXP)
    {
      mooseAssert(coef != 0, "Coefficient must not be zero.");
      if (derivative)
        x = exp(x / coef) / coef;
      else
        x = exp(x / coef);
    }
    else if (transform == ROMInputTransform::LOG)
    {
      mooseAssert(x + coef > 0, "Sum must be greater than 0.");
      if (derivative)
        x = 1.0 / (x + coef);
      else
        x = log(x + coef);
    }
    else if (transform == ROMInputTransform::LINEAR)
    {
      mooseAssert(coef == 0.0, "Coefficient cannot be supplied with linear transformation");
      if (derivative)
        x = 1.0;
    }
  }

  /*
   * Calculates and returns vector utilized in assign values
   * @param p Partition index
   * @return Vector that preallocates indexing calculations for polynomial calculation
   */
  std::vector<unsigned int> getMakeFrameHelper(const unsigned int p) const;

  /*
   * Calculates and returns the transformed limits for the ROM calculations
   * Indexes are [partition][tile][output][input].
   * Inputs ordering is
   * input[0]: cell_old
   * input[1]: wall_old
   * input[2]: trial stress,
   * input[3]: effective strain old,
   * input[4]: temperature
   * input[5]: environmental factor (optional)
   * output ordering is:
   * output[0]: cell dislocations increment
   * output[1]: wall dislocations increment
   * output[2]: strain increment
   * @param p Partition index
   * @param limits Human readable limits
   * @return Multi-dimentional vector of transformed limits
   */
  std::vector<std::vector<std::vector<std::vector<Real>>>>
  getTransformedLimits(const unsigned int p,
                       const std::vector<std::vector<std::vector<Real>>> limits);

  /*
   * Returns vector of the functions to use for the conversion of input variables.
   * Indexes are [partition][tile][output][input].
   * Inputs ordering is
   * input[0]: cell_old
   * input[1]: wall_old
   * input[2]: trial stress,
   * input[3]: effective strain old,
   * input[4]: temperature
   * input[5]: environmental factor (optional)
   * output ordering is:
   * output[0]: cell dislocations increment
   * output[1]: wall dislocations increment
   * output[2]: strain increment
   * @return vector of the functions to use for the conversion of input variables.
   */
  virtual std::vector<std::vector<std::vector<std::vector<ROMInputTransform>>>> getTransform()
  {
    checkJSONKey("transform");
    return _json["transform"]
        .template get<std::vector<std::vector<std::vector<std::vector<ROMInputTransform>>>>>();
  }

  /*
   * Returns factors for the functions for the conversion functions given in getTransform
   * Indexes are [partition][tile][output][input].
   * Inputs ordering is
   * input[0]: cell_old
   * input[1]: wall_old
   * input[2]: trial stress,
   * input[3]: effective strain old,
   * input[4]: temperature
   * input[5]: environmental factor (optional)
   * output ordering is:
   * output[0]: cell dislocations increment
   * output[1]: wall dislocations increment
   * output[2]: strain increment
   * @return factors for the functions for the conversion functions given in getTransform
   */
  virtual std::vector<std::vector<std::vector<std::vector<Real>>>> getTransformCoefs()
  {
    checkJSONKey("transform_coefs");
    return _json["transform_coefs"]
        .template get<std::vector<std::vector<std::vector<std::vector<Real>>>>>();
  }

  /* Optional method that returns human-readable limits used for normalization. Default is to just
   * use the input limits.
   * Indexes are [partition][tile][input][upper/lower].
   * Inputs ordering is
   * input[0]: cell_old
   * input[1]: wall_old
   * input[2]: trial stress,
   * input[3]: effective strain old,
   * input[4]: temperature
   * input[5]: environmental factor (optional)
   * @return human-readable limits for the normalization limits
   */
  virtual std::vector<std::vector<std::vector<std::vector<Real>>>> getNormalizationLimits()
  {
    if (_json.contains("normalization_limits"))
      return _json["normalization_limits"]
          .template get<std::vector<std::vector<std::vector<std::vector<Real>>>>>();

    return getInputLimits();
  }

  /* Returns human-readable limits for the inputs.
   * Indexes are [partition][tile][input][upper/lower].
   * Inputs ordering is
   * input[0]: cell_old
   * input[1]: wall_old
   * input[2]: trial stress,
   * input[3]: effective strain old,
   * input[4]: temperature
   * input[5]: environmental factor (optional)
   * @return human-readable limits for the input limits
   */
  virtual std::vector<std::vector<std::vector<std::vector<Real>>>> getInputLimits()
  {
    checkJSONKey("input_limits");
    return _json["input_limits"]
        .template get<std::vector<std::vector<std::vector<std::vector<Real>>>>>();
  }

  /*
   * Material specific coefficients multiplied by the Legendre polynomials for each of the input
   * variables
   * @return Legendre polynomial coefficients
   */
  virtual std::vector<std::vector<std::vector<std::vector<Real>>>> getCoefs()
  {
    checkJSONKey("coefs");
    return _json["coefs"].template get<std::vector<std::vector<std::vector<std::vector<Real>>>>>();
  }

  /*
   * Material specific orientations of tiling
   * variables. Indexing is partition, then input
   * @return Vector of a vector declaring tiling orientation
   */
  virtual std::vector<std::vector<unsigned int>> getTilings()
  {
    if (_json.contains("tiling"))
      return _json["tiling"].template get<std::vector<std::vector<unsigned int>>>();

    if (_environmental)
      return {{1, 1, 1, 1, 1, 1}};
    return {{1, 1, 1, 1, 1}};
  };

  /*
   * Minimum strain value allowed by the ROM. This is material specific, and needs to be overwritten
   * by individual roms and each partition
   * @return Vector of material specific ROM low strain value for each partition
   */
  virtual std::vector<Real> getStrainCutoff()
  {
    checkJSONKey("cutoff");
    return _json["cutoff"].template get<std::vector<Real>>();
  }

  const GenericVariableValue<is_ad> & _temperature;

  const GenericMaterialProperty<Real, is_ad> * _environmental;

  /*
   * Vector of vectors WindowFailure enum that informs how to handle input that is outside of the
   * limits. Shape is number of inputs by 2 (lower and upper window enum)
   */
  std::vector<std::pair<WindowFailure, WindowFailure>> _window_failure;

  const bool _verbose;

  GenericMaterialProperty<Real, is_ad> & _cell_dislocations;
  const MaterialProperty<Real> & _cell_dislocations_old;

  const Real _max_cell_increment;

  const Function * const _cell_function;

  GenericReal<is_ad> _cell_dislocation_increment;

  GenericMaterialProperty<Real, is_ad> & _wall_dislocations;
  const MaterialProperty<Real> & _wall_dislocations_old;

  const Real _max_wall_increment;

  const Function * const _wall_function;

  const Function * const _stress_function;

  GenericReal<is_ad> _wall_dislocation_increment;

  const unsigned int _cell_input_index;

  const unsigned int _wall_input_index;

  const unsigned int _stress_input_index;

  const unsigned int _old_strain_input_index;

  const unsigned int _temperature_input_index;

  const unsigned int _environmental_input_index;

  const unsigned int _cell_output_index;

  const unsigned int _wall_output_index;

  const unsigned int _strain_output_index;

  const Function * const _creep_strain_old_forcing_function;

  unsigned int _num_partitions;

  std::vector<unsigned int> _num_tiles;

  unsigned int _num_inputs;

  unsigned int _num_outputs;

  std::vector<unsigned int> _degree;

  std::vector<unsigned int> _num_coefs;

  std::vector<std::vector<std::vector<std::vector<ROMInputTransform>>>> _transform;

  std::vector<std::vector<std::vector<std::vector<Real>>>> _transform_coefs;

  std::vector<std::vector<std::vector<std::vector<Real>>>> _input_limits;

  std::vector<std::vector<std::vector<std::vector<Real>>>> _normalization_limits;

  std::vector<std::vector<std::vector<std::vector<Real>>>> _coefs;

  std::vector<std::vector<std::vector<std::vector<std::vector<Real>>>>>
      _transformed_normalization_limits;

  std::vector<std::vector<unsigned int>> _makeframe_helper;

  GenericMaterialProperty<Real, is_ad> & _creep_rate;

  GenericMaterialProperty<Real, is_ad> & _cell_rate;

  GenericMaterialProperty<Real, is_ad> & _wall_rate;

  MaterialProperty<Real> & _extrapolation;

  GenericMaterialProperty<Real, is_ad> & _second_partition_weight;

  GenericReal<is_ad> _derivative;

  std::vector<GenericReal<is_ad>> _input_values;

  std::vector<Real> _old_input_values;

  std::vector<std::vector<std::vector<GenericReal<is_ad>>>> _rom_inputs;

  std::vector<std::vector<std::vector<std::vector<GenericReal<is_ad>>>>> _polynomial_inputs;

  std::vector<std::vector<std::vector<GenericReal<is_ad>>>> _precomputed_vals;

  std::vector<std::pair<Real, Real>> _global_limits;

  std::vector<std::vector<GenericReal<is_ad>>> _non_stress_weights;

  std::vector<std::vector<GenericReal<is_ad>>> _weights;

  std::vector<GenericReal<is_ad>> _partition_weights;

  std::vector<GenericReal<is_ad>> _dpartition_weight_dstress;

  std::vector<std::vector<unsigned int>> _tiling;

  std::vector<Real> _cutoff;

  Real _stress_ucf;

  GenericMaterialProperty<Real, is_ad> & _wall_dislocations_step;
  GenericMaterialProperty<Real, is_ad> & _cell_dislocations_step;

  RankTwoTensor _plastic_strain_increment;

  MaterialProperty<Real> & _number_of_substeps;

  std::vector<std::string> _index_name;

  void checkJSONKey(const std::string & key);

  nlohmann::json _json;

  usingTransientInterfaceMembers;
  using Material::_name;
  using Material::_q_point;
  using Material::_qp;
  using Material::coupledGenericValue;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::computeResidual;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::computeDerivative;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::_apply_strain;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::initQpStatefulProperties;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::outputIterationStep;
  using RadialReturnCreepStressUpdateBaseTempl<is_ad>::outputIterationSummary;
};

typedef LAROMANCEStressUpdateBaseTempl<false> LAROMANCEStressUpdateBase;
typedef LAROMANCEStressUpdateBaseTempl<true> ADLAROMANCEStressUpdateBase;