hm_read_mat102.F

Go to the documentation of this file.

Copyright>        OpenRadioss
Copyright>        Copyright (C) 1986-2025 Altair Engineering Inc.
Copyright>
Copyright>        This program is free software: you can redistribute it and/or modify
Copyright>        it under the terms of the GNU Affero General Public License as published by
Copyright>        the Free Software Foundation, either version 3 of the License, or
Copyright>        (at your option) any later version.
Copyright>
Copyright>        This program is distributed in the hope that it will be useful,
Copyright>        but WITHOUT ANY WARRANTY; without even the implied warranty of
Copyright>        MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
Copyright>        GNU Affero General Public License for more details.
Copyright>
Copyright>        You should have received a copy of the GNU Affero General Public License
Copyright>        along with this program.  If not, see <https://www.gnu.org/licenses/>.
Copyright>
Copyright>
Copyright>        Commercial Alternative: Altair Radioss Software
Copyright>
Copyright>        As an alternative to this open-source version, Altair also offers Altair Radioss
Copyright>        software under a commercial license.  Contact Altair to discuss further if the
Copyright>        commercial version may interest you: https://www.altair.com/radioss/.
!||====================================================================
!||    hm_read_mat102           ../starter/source/materials/mat/mat102/hm_read_mat102.F
!||--- called by ------------------------------------------------------
!||    hm_read_mat              ../starter/source/materials/mat/hm_read_mat.F90
!||--- calls      -----------------------------------------------------
!||    ancmsg                   ../starter/source/output/message/message.F
!||    hm_get_floatv            ../starter/source/devtools/hm_reader/hm_get_floatv.F
!||    hm_get_intv              ../starter/source/devtools/hm_reader/hm_get_intv.F
!||    hm_option_is_encrypted   ../starter/source/devtools/hm_reader/hm_option_is_encrypted.F
!||    init_mat_keyword         ../starter/source/materials/mat/init_mat_keyword.F
!||--- uses       -----------------------------------------------------
!||    elbuftag_mod             ../starter/share/modules1/elbuftag_mod.F
!||    message_mod              ../starter/share/message_module/message_mod.F
!||    submodel_mod             ../starter/share/modules1/submodel_mod.F
!||====================================================================
      SUBROUTINE hm_read_mat102(UPARAM ,MAXUPARAM,NUPARAM  ,ISRATE   ,
     .                          NUVAR  ,NFUNC    ,PARMAT   ,
     .                          UNITAB ,MAT_ID   ,TITR     ,MTAG     , LSUBMODEL,
     .                          PM     ,IPM      ,MATPARAM )
C-----------------------------------------------
C   D e s c r i p t i o n
C-----------------------------------------------
C   READ MAT LAW102 WITH HM READER ( TO BE COMPLETED )
C
C   DUMMY ARGUMENTS DESCRIPTION:
C   ===================
C
C     NAME            DESCRIPTION                         
C
C     IPM             MATERIAL ARRAY(INTEGER)
C     PM              MATERIAL ARRAY(REAL)
C     UNITAB          UNITS ARRAY
C     MAT_ID          MATERIAL ID(INTEGER)
C     TITR            MATERIAL TITLE
C     LSUBMODEL       SUBMODEL STRUCTURE   
C
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE unitab_mod
      USE elbuftag_mod            
      USE message_mod      
      USE submodel_mod , ONLY : submodel_data, nsubmod
      USE matparam_def_mod
      USE names_and_titles_mod , ONLY : nchartitle, ncharkey, ncharfield
      USE constant_mod , ONLY : pi, hundred80, four, zero, infinity, nine, one, six, three, two
      USE precision_mod , ONLY : wp
 
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
      IMPLICIT NONE
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "units_c.inc"
#include      "param_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      TYPE (UNIT_TYPE_),INTENT(IN)                       :: UNITAB
      REAL(KIND=WP), DIMENSION(NPROPM),INTENT(INOUT)     :: PM
      REAL(KIND=WP), DIMENSION(100),INTENT(INOUT)        :: PARMAT
      REAL(KIND=wp), DIMENSION(MAXUPARAM) ,INTENT(INOUT) :: uparam
      INTEGER, DIMENSION(NPROPMI),INTENT(INOUT)          :: IPM
      INTEGER, INTENT(INOUT)                             :: ISRATE,NFUNC,MAXUPARAM,NUPARAM,NUVAR
      TYPE(mlaw_tag_),INTENT(INOUT)                      :: MTAG
      INTEGER,INTENT(IN)                                 :: MAT_ID
      CHARACTER(LEN=NCHARTITLE) ,INTENT(IN)              :: TITR
      TYPE(submodel_data),INTENT(IN)                     :: LSUBMODEL(NSUBMOD)
      TYPE(matparam_struct_) ,INTENT(INOUT)              :: MATPARAM
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      REAL(KIND=wp)    :: e,nu,c,pstar,amax,g, delta,stifint,pmin,phi_deg
      REAL(KIND=wp)    :: a0,a1,a2,rho0,rhor
      INTEGER          :: IFORM
      REAL(KIND=wp)    :: phi,k,alpha
      LOGICAL          :: IS_ENCRYPTED,IS_AVAILABLE
      CHARACTER*64     :: CHAIN
C-----------------------------------------------
C   S o u r c e   L i n e s
C-----------------------------------------------
      is_encrypted = .false.
      is_available = .false.
      k = 0.0d0
 
      CALL hm_option_is_encrypted(is_encrypted)
 
      CALL hm_get_intv  ('IFORM',iform  ,is_available, lsubmodel)
      
      CALL hm_get_floatv('MAT_RHO'      ,rho0  ,is_available, lsubmodel, unitab)
      
      CALL hm_get_floatv('MAT_E'        ,e     ,is_available, lsubmodel, unitab)
      CALL hm_get_floatv('MAT_NU'       ,nu    ,is_available, lsubmodel, unitab)
      
      CALL hm_get_floatv('MAT102_C'     ,c     ,is_available, lsubmodel, unitab)
      CALL hm_get_floatv('MAT102_PHI'   ,phi   ,is_available, lsubmodel, unitab)
      CALL hm_get_floatv('MAT102_AMAX'  ,amax  ,is_available, lsubmodel, unitab)
 
      CALL hm_get_floatv('MAT102_PMIN'  ,pmin    ,is_available, lsubmodel, unitab) 
 
      !   C   : Mohr-Coulomb cohesion (intercept) (in Pa)
      ! PHI   : Mohr-Coulomb angle of internal friction (angle, degree)
      ! Iform : Formulation Flag (circumscribed,inscribed,middle)
 
C----------------------------------------------- 
C     UNITS   
C----------------------------------------------- 
      phi_deg = phi
      phi     = phi*pi/hundred80
C----------------------------------------------- 
C     DEFAULTS   
C----------------------------------------------- 
      IF(iform<=0 .OR. iform>=4)iform=2
C----------------------------------------------- 
C     YIELD PARAMETERS   
C-----------------------------------------------
      g=e/two/(one+nu) 
      SELECT CASE(iform)
        CASE(1)
          k     = six*c*cos(phi)/sqrt(three)/(three-sin(phi))
          alpha = two*sin(phi)/sqrt(three)/(three-sin(phi))
        CASE(2)
          k     = six*c*cos(phi)/sqrt(three)/(three+sin(phi))
          alpha = two*sin(phi)/sqrt(three)/(three+sin(phi))
        CASE(3)
          k     = three*c*cos(phi)/sqrt(nine+three*sin(phi)*sin(phi))
          alpha = sin(phi)/sqrt(nine+three*sin(phi)*sin(phi))
      END SELECT
      a0 = k*k
      a1 = six*k*alpha
      a2 = nine*alpha*alpha
 
      IF(e<=zero)THEN
        chain='YOUNG MODULUS MUST BE DEFINED                           '
        CALL ancmsg(msgid=829, msgtype=msgerror, anmode=aninfo, i1=10, i2=mat_id, c1='ERROR', c2=titr, c3=chain)
      ENDIF 
 
      IF(nu<=zero)THEN
        chain='POISSON RATIO MUST BE DEFINED                           '
        CALL ancmsg(msgid=829, msgtype=msgerror, anmode=aninfo, i1=10, i2=mat_id, c1='ERROR', c2=titr, c3=chain)
      ENDIF 
 
      pstar = -infinity
      IF(a2==zero .AND. a1/=zero)THEN 
        pstar=-a0/a1 
      ELSEIF(a2/=zero)THEN
        delta = a1*a1-four*a0*a2
        !If intersection with the axis
        IF(delta >= zero)THEN
          delta=sqrt(delta)
          pstar = (-a1+delta)/two/a2
        
        ELSE
          pstar = -a1/two/a2 ! extremum          
    !      chain='FITTED YIELD FUNCTION HAS NO ROOT. CHECK INPUT PARAMETER        '      
    !      CALL ANCMSG(  MSGID=829, MSGTYPE=MSGWARNING, ANMODE=ANINFO, I1=10, I2=MAT_ID, C1='WARNING', C2=TITR, C3=chain)
        ENDIF
      ELSE
        !do nothing let user do what he wants
        pstar = -infinity
      ENDIF
C
      IF(amax==zero) amax = infinity
      IF(pmin==zero) pmin =-infinity
 
      rhor=rho0
      pm(1) = rhor
      pm(89)= rho0
      pm(37)= pmin
      
      israte=0
 
C----------------------------------------------- 
C     USER MATERIAL PARAMETERS DEFINITION   
C----------------------------------------------- 
      nuparam   = 11
      uparam(1) = c
      uparam(2) = phi
      uparam(3) = pstar
      uparam(4) = a0
      uparam(5) = a1
      uparam(6) = a2
      uparam(7) = amax
      uparam(8) = g
      uparam(9) = iform
      uparam(10)= e
      uparam(11)= nu
      nuvar     = 0 
      nfunc     = 0 
      stifint   = e
      parmat(1) = stifint/three !
      parmat(2) = stifint
      parmat(3) = nu
c---------------------------   
      mtag%G_PLA  = 1
      mtag%L_PLA  = 1
!---------------------------------------------------------   
      !MATPARAM keywords
      matparam%IEOS = 18  ! Linear EOS is used by default
      ipm(4)        = 18
      pm(32)        = e / three/(one - two*nu) ! Bulk modulus for default EOS
!---------------------------------------------------------   
 
      CALL init_mat_keyword(matparam,"elasto_plastic")
 
      ! Material compatibility with /EOS option
      CALL INIT_MAT_KEYWORD(MATPARAM,"eos")
 
      ! EOS/Thermo keyword for pressure treatment in elements
      CALL INIT_MAT_KEYWORD(MATPARAM,"hydro_eos")
 
      ! Properties compatibility
      CALL INIT_MAT_KEYWORD(MATPARAM,"solid_isotropic")
      CALL INIT_MAT_KEYWORD(MATPARAM,"sph")
C------------------------------------------------- 
C     LISTING OUTPUT 
C------------------------------------------------- 
      WRITE(IOUT,1001) TRIM(TITR),MAT_ID,102
      WRITE(IOUT,1000)
      IF(IS_ENCRYPTED)THEN
        WRITE(IOUT,'(5X,A,//)')'CONFIDENTIAL DATA'
      ELSE
        WRITE(IOUT,1002)RHO0,RHOR    
        WRITE(IOUT,1100)E,NU,C,PHI_DEG,PMIN
        WRITE(IOUT,1200)IFORM
        SELECT CASE(IFORM)
          CASE(1)
            WRITE(IOUT,1201)
          CASE(2)
            WRITE(IOUT,1202)
          CASE(3)
            WRITE(IOUT,1203)
          CASE(4)
            WRITE(IOUT,1204)
        END SELECT
        WRITE(IOUT,1300)A0,A1,A2,AMAX,PSTAR
      ENDIF
C
 1000 FORMAT(
     & 5X,'  EXTENDED DRUCKER-PRAGER MATERIAL (DPRAG2) ',/,
     & 5X,'  ----------------------------------------- ')
 1001 FORMAT(/
     & 5X,A,/,
     & 5X,'MATERIAL NUMBER . . . . . . . . . . . .=',I10/,
     & 5X,'MATERIAL LAW. . . . . . . . . . . . . .=',I10/)
 1002 FORMAT(
     & 5X,'INITIAL DENSITY . . . . . . . . . . . .=',1PG20.13/,
     & 5X,'REFERENCE DENSITY . . . . . . . . . . .=',1PG20.13/)  
 1100 FORMAT(
     & 5X,'YOUNG MODULUS . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'POISSON RATIO . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'COHESION. . . . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'ANGLE OF INTERNAL FRICTION. . . . . . .=',1PG20.13/
     & 5X,'MINIMUM PRESSURE. . . . . . . . . . . .=',1PG20.13)   
 1200 FORMAT(
     & 5X,'DRUCKER-PRAGER MATERIAL CRITERION DEFINED FROM MOHR-COULOMB PARAMETERS',/,
     & 5X,'FORMULATION FLAG. . . . . . . . . . . .=',I10)     
 1201 FORMAT(
     & 5X,'-> CIRCUMBSCRIBED CRITERIA')   
 1202 FORMAT(
     & 5X,'-> MIDDLE CRITERIA') 
 1203 FORMAT(
     & 5X,'-> INSCRIBED CRITERIA') 
 1204 FORMAT(
     & 5X,'-> ORIGINAL MOHR-COULOMB CRITERIA')   
 1300 FORMAT(
     & 5X,'PARAMETERS USED TO DEFINE CRITERIA',/,
     & 5X,'A0. . . . . . . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'A1. . . . . . . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'A2. . . . . . . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'AMAX. . . . . . . . . . . . . . . . . .=',1PG20.13/
     & 5X,'YIELD FUNCTION PRESSURE ROOT. . . . . .=',1PG20.13//)        
 
      RETURN
      END