#include "../MEF90/mef90.inc" module m_MEF90_DefMechAT1exp #include "petsc/finclude/petsc.h" use petscsys use m_MEF90_Parameters use m_MEF90_LinAlg use m_MEF90_DefMechAT_class use iso_c_binding implicit none(type) private public :: MEF90DefMechAT1exp_Type !!! AT1exp, a variant of AT1 model with an exponential stiffness interpolation !!! function: !!! !!! $a_b(s) = 1 + (e^{-bs} - 1) / (1 - e^-b)$ if $b /= 0$ !!! $a_0(s) = 1-s$ !!! !!! $a$ is convex if $b > 1$ and !!! $a_b'(0) = -b / (1-e^{-b}) < -2$ if $b < 3/2$ !!! type, extends(MEF90DefMechAT_Type) :: MEF90DefMechAT1exp_Type PetscReal :: b = 1.0_kr contains procedure, pass(self) :: a => aAT1exp procedure, pass(self) :: Da => DaAT1exp procedure, pass(self) :: D2a => D2aAT1exp procedure, pass(self) :: w => wAT1exp procedure, pass(self) :: Dw => DwAT1exp procedure, pass(self) :: D2w => D2wAT1exp procedure, pass(self) :: setFromOptions => MEF90DefMechAT1exp_setFromOptions procedure, pass(self) :: view => MEF90DefMechAT1exp_view end type MEF90DefMechAT1exp_Type contains #undef __FUNCT__ #define __FUNCT__ "MEF90DefMechAT1exp_view" !!! author: Blaise Bourdin (2025, bourdin@mcmaster.ca) !!! !!! MEF90DefMechAT1exp_view: the default viewer for a MEF90_DefMechAT_Type !!! subroutine MEF90DefMechAT1exp_View(self,viewer,ierr) class(MEF90DefMechAT1exp_Type), intent(in) :: self type(tPetscViewer), intent(in) :: viewer PetscErrorCode, intent(inout) :: ierr character(len=MEF90MXSTRLEN, kind=c_char) :: IOBuffer character(len=MEF90MXSTRLEN, kind=c_char) :: viewerType PetscCall(MEF90DefMechAT_View(self, viewer, ierr)) PetscCall(PetscViewerGetType(viewer, viewerType, ierr)) if (viewerType == 'ascii') then write(IOBuffer, "(A,': Options for MEF90DefMechAT1exp_type\n')") trim(self%prefix) // "damage_AT1exp" PetscCall(PetscViewerASCIIPrintf(viewer, IOBuffer, ierr)) write(IOBuffer, "(' b : ',ES12.5,' [N.m^(-1)]\n')") self%b PetscCall(PetscViewerASCIIPrintf(viewer, IOBuffer, ierr)) write(IOBuffer, "(' cw : ',ES12.5,' [N.m^(-1)]\n')") self%cw PetscCall(PetscViewerASCIIPrintf(viewer, IOBuffer, ierr)) end if end subroutine MEF90DefMechAT1exp_View #undef __FUNCT__ #define __FUNCT__ "MEF90DefMechAT1exp_setFromOptions" !!! author: Blaise Bourdin (2025, bourdin@mcmaster.ca) !!! !!! MEF90DefMechAT1exp_setFromOptions: initializes a MEF90_DefMechAT1exp_Type from options !!! subroutine MEF90DefMechAT1exp_setFromOptions(self,ierr) class(MEF90DefMechAT1exp_Type), intent(inout) :: self PetscErrorCode,intent(inout) :: ierr self%b = 1.0_kr self%cw = 0.597674895613893_Kr PetscCall(PetscOptionsBegin(self%comm, trim(self%prefix) // "damage_AT1exp_", "Options for MEF90DefMechATexp_type", "mef90DefMech", ierr)) PetscCall(PetscOptionsReal('-b', 'b parameter', '[]', self%b, self%b, PETSC_NULL_BOOL, ierr)) PetscCall(PetscOptionsReal('-cw', 'cw parameter', '[]', self%cw, self%cw, PETSC_NULL_BOOL, ierr)) PetscCall(PetscOptionsEnd(ierr)) self%aorder = 3 self%worder = 1 self%type = 'MEF90DefMechAT1exp' PetscCall(MEF90DefMechAT_setFromOptions(self, ierr)) end subroutine MEF90DefMechAT1exp_setFromOptions #undef __FUNCT__ #define __FUNCT__ "aAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! aAT1exp: the "a" function of the standard AT1exp model !!! PetscReal function aAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha if (self%b == 0.0_kr) then aAT1exp = 1.0_kr - alpha else aAT1exp = 1.0_kr + (exp(-self%b * alpha) - 1.0_kr) / (1.0_kr - exp(-self%b)) end if end function aAT1exp #undef __FUNCT__ #define __FUNCT__ "DaAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! DaAT1exp: the derivative of the "a" function of the standard AT1exp model !!! PetscReal function DaAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha if (self%b == 0.0_kr) then DaAT1exp = -1.0_kr else DaAT1exp = -self%b * exp(-self%b * alpha) / (1.0_kr - exp(-self%b)) end if end function DaAT1exp #undef __FUNCT__ #define __FUNCT__ "D2aAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! D2aAT1exp: the second derivative of the "a" function of the standard AT1exp model !!! PetscReal function D2aAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha if (self%b == 0.0_kr) then D2aAT1exp = 0.0_kr else D2aAT1exp = self%b**2 * exp(-self%b * alpha) / (1.0_kr - exp(-self%b)) end if end function D2aAT1exp #undef __FUNCT__ #define __FUNCT__ "wAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! wAT1exp: the "w" function of the standard AT1exp model, i.e. !!! PetscReal function wAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha wAT1exp = alpha end function wAT1exp #undef __FUNCT__ #define __FUNCT__ "DwAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! DwAT1exp: the derivative of the "w" function of the standard AT1exp model !!! PetscReal function DwAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha DwAT1exp = 1.0_kr end function DwAT1exp #undef __FUNCT__ #define __FUNCT__ "D2wAT1exp" !!! author: Blaise Bourdin (2020, bourdin@lsu.edu) !!! !!! D2wAT1exp: the second derivative of the "w" function of the standard AT1exp model !!! PetscReal function D2wAT1exp(self, alpha) class(MEF90DefMechAT1exp_Type), intent(IN) :: self PetscReal :: alpha D2wAT1exp = 0.0_kr end function D2wAT1exp end module m_MEF90_DefMechAT1exp