umat.C

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

#include <stdexcept>

#include "omi_for_c.h"
#include "build_stiffness_matrix.h"

extern "C" void FOR_NAME(umat, UMAT)(double * stress,
                                     double * /*statev*/,
                                     double * ddsdde,
                                     double * /*sse*/,
                                     double * /*spd*/,
                                     double * /*scd*/,
                                     double * /*rpl*/,
                                     double * /*ddsddt*/,
                                     double * /*drplde*/,
                                     double * /*drpldt*/,
                                     double * stran,
                                     double * dstran,
                                     double * /*time*/,
                                     double * /*dtime*/,
                                     double * /*temp*/,
                                     double * /*dtemp*/,
                                     double * /*predef*/,
                                     double * /*dpred*/,
                                     char * /*cmname*/,
                                     int * /*ndi*/,
                                     int * /*nshr*/,
                                     int * /*ntens*/,
                                     int * /*nstatv*/,
                                     double * props,
                                     int * nprops,
                                     double * /*coords*/,
                                     double * /*drot*/,
                                     double * /*pnewdt*/,
                                     double * /*celent*/,
                                     double * /*dfgrd0*/,
                                     double * /*dfgrd1*/,
                                     int * /*noel*/,
                                     int * /*npt*/,
                                     int * /*layer*/,
                                     int * /*kspt*/,
                                     int * /*kstep*/,
                                     int * /*kinc*/,
                                     short /*cmname_len*/)
{
  if (*nprops != 2)
    throw std::invalid_argument("This UMAT requires exactly two properties.");

  double E = props[0];
  double nu = props[1];
  double G = E / 2.0 / (1.0 + nu);
  double lambda = 2.0 * G * nu / (1.0 - 2.0 * nu);
  double eps[6];

  // Build stress as in
  // https://github.com/michael-schw/Abaqus-UMAT-Cpp-Subroutine/blob/main/umat.cpp
  for (int i = 0; i < 6; i++)
    eps[i] = stran[i] + dstran[i];

  auto eps_trace = eps[0] + eps[1] + eps[2];

  for (int i = 0; i < 3; i++)
  {
    stress[i] = lambda * eps_trace + 2.0 * G * eps[i];
    stress[i + 3] = G * eps[i + 3];
  }

  Eigen::Matrix<double, 6, 6> C;
  buildStiffnessMatrix(C, G, lambda);

  for (int i = 0; i < 6; i++)
    for (int j = 0; j < 6; j++)
      ddsdde[6 * i + j] = 0.0;

  for (int i = 0; i < 3; i++)
  {
    ddsdde[6 * i + 0] = C(0, i);
    ddsdde[6 * i + 1] = C(1, i);
    ddsdde[6 * i + 2] = C(2, i);
    ddsdde[6 * (i + 3) + (i + 3)] = C(i + 3, i + 3);
  }
}