OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
hm_read_mat36.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| hm_read_mat36 ../starter/source/materials/mat/mat036/hm_read_mat36.F
25!||--- called by ------------------------------------------------------
26!|| hm_read_mat ../starter/source/materials/mat/hm_read_mat.F90
27!||--- calls -----------------------------------------------------
28!|| ancmsg ../starter/source/output/message/message.F
29!|| hm_get_float_array_index ../starter/source/devtools/hm_reader/hm_get_float_array_index.F
30!|| hm_get_float_array_index_dim ../starter/source/devtools/hm_reader/hm_get_float_array_index_dim.F
31!|| hm_get_floatv ../starter/source/devtools/hm_reader/hm_get_floatv.F
32!|| hm_get_floatv_dim ../starter/source/devtools/hm_reader/hm_get_floatv_dim.F
33!|| hm_get_int_array_index ../starter/source/devtools/hm_reader/hm_get_int_array_index.F
34!|| hm_get_intv ../starter/source/devtools/hm_reader/hm_get_intv.F
35!|| hm_option_is_encrypted ../starter/source/devtools/hm_reader/hm_option_is_encrypted.f
36!|| init_mat_keyword ../starter/source/materials/mat/init_mat_keyword.F
37!||--- uses -----------------------------------------------------
38!|| elbuftag_mod ../starter/share/modules1/elbuftag_mod.F
39!|| hm_option_read_mod ../starter/share/modules1/hm_option_read_mod.f
40!|| message_mod ../starter/share/message_module/message_mod.F
41!|| submodel_mod ../starter/share/modules1/submodel_mod.F
42!||====================================================================
43 SUBROUTINE hm_read_mat36(UPARAM ,MAXUPARAM,NUPARAM ,NUVAR ,NVARTMP ,
44 . IFUNC ,MAXFUNC ,MFUNC ,PARMAT ,UNITAB ,
45 . ID ,MTAG ,TITR ,LSUBMODEL,PM ,
46 . ISRATE ,MATPARAM )
47C-----------------------------------------------
48C D e s c r i p t i o n
49C-----------------------------------------------
50C
51C DUMMY ARGUMENTS DESCRIPTION:
52C ===================
53C
54C NAME DESCRIPTION
55C
56C IPM MATERIAL ARRAY(INTEGER)
57C PM MATERIAL ARRAY(REAL)
58C UNITAB UNITS ARRAY
59C ID MATERIAL ID(INTEGER)
60C TITR MATERIAL TITLE
61C LSUBMODEL SUBMODEL STRUCTURE
62C
63C-----------------------------------------------
64C M o d u l e s
65C-----------------------------------------------
66 USE unitab_mod
67 USE elbuftag_mod
68 USE message_mod
69 USE submodel_mod
70 USE matparam_def_mod
71 USE hm_option_read_mod
72 USE names_and_titles_mod , ONLY : nchartitle
73C-----------------------------------------------
74C I m p l i c i t T y p e s
75C-----------------------------------------------
76#include "implicit_f.inc"
77C-----------------------------------------------
78C C o m m o n B l o c k s
79C-----------------------------------------------
80#include "units_c.inc"
81#include "param_c.inc"
82C-----------------------------------------------
83C D u m m y A r g u m e n t s
84C-----------------------------------------------
85 TYPE (UNIT_TYPE_),INTENT(IN) ::UNITAB
86 my_real, DIMENSION(NPROPM) ,INTENT(INOUT) :: PM
87 my_real, DIMENSION(100) ,INTENT(INOUT) :: PARMAT
88 my_real, DIMENSION(MAXUPARAM) ,INTENT(INOUT) :: uparam
89 INTEGER, DIMENSION(MAXFUNC) ,INTENT(INOUT) :: IFUNC
90 INTEGER, INTENT(INOUT) :: MFUNC,NUPARAM,NUVAR,NVARTMP,ISRATE
91 TYPE(mlaw_tag_),INTENT(INOUT) :: MTAG
92 INTEGER,INTENT(IN) :: ID,MAXFUNC,MAXUPARAM
93 CHARACTER(LEN=NCHARTITLE) ,INTENT(IN) :: TITR
94 TYPE(submodel_data),INTENT(IN) :: LSUBMODEL(*)
95 TYPE(matparam_struct_) ,INTENT(INOUT) :: MATPARAM
96C-----------------------------------------------
97C L o c a l V a r i a b l e s
98C-----------------------------------------------
99 INTEGER :: NBMAT, MAT_ID ! Number of declared materials
100 INTEGER :: I,J,VP,YLDCHECK
101 INTEGER :: RHOFLAG,ICOMP,NRATE1,NRATE,IPFUN,IFUNCE,ISRAT,ISMOOTH,
102 . nbline,nbread,ifail,opte,ilaw,nfunc
103 my_real :: rho0, rhor,e,nu,g,c1,soundsp, epsmax,epsr1,epsr2,epsf,fisokin,fcut,
104 . pscal_unit,pscale,einf,ce ,
105 . yfac(maxfunc),rate(max(1,maxfunc)),strainrate_unit(maxfunc),yfac_unit(maxfunc)
106 LOGICAL :: IS_AVAILABLE,IS_ENCRYPTED
107C-----------------------------------------------
108C S o u r c e L i n e s
109C-----------------------------------------------
110 is_encrypted = .false.
111 is_available = .false.
112 rate(1) = zero
113C--------------------------------------------------
114C EXTRACT DATA (IS OPTION CRYPTED)
115C--------------------------------------------------
116 CALL hm_option_is_encrypted(is_encrypted)
117C-----------------------------------------------
118 ilaw = 36
119 CALL hm_get_floatv('MAT_RHO' ,rho0 ,is_available, lsubmodel, unitab)
120 CALL hm_get_floatv('Refer_Rho',rhor ,is_available, lsubmodel, unitab)
121C-----------------------------------------------
122Card1
123 CALL hm_get_floatv('MAT_E' ,e ,is_available, lsubmodel, unitab)
124 CALL hm_get_floatv('MAT_NU' ,nu ,is_available, lsubmodel, unitab)
125 CALL hm_get_floatv('MAT_EPS' ,epsmax ,is_available, lsubmodel, unitab)
126 CALL hm_get_floatv('MAT_EPST1',epsr1 ,is_available, lsubmodel, unitab)
127 CALL hm_get_floatv('MAT_EPST2',epsr2 ,is_available, lsubmodel, unitab)
128C-----------------------------------------------
129Card2
130 CALL hm_get_intv ('NFUNC' ,nrate ,is_available,lsubmodel)
131 CALL hm_get_intv ('Fsmooth' ,ismooth ,is_available,lsubmodel)
132 CALL hm_get_floatv('MAT_HARD' ,fisokin ,is_available, lsubmodel, unitab)
133 CALL hm_get_floatv('Fcut' ,fcut ,is_available, lsubmodel, unitab)
134 CALL hm_get_floatv('MAT_Epsilon_F',epsf ,is_available, lsubmodel, unitab)
135 CALL hm_get_intv ('Vflag' ,vp ,is_available,lsubmodel)
136Card3
137 CALL hm_get_intv ('Xr_fun' ,ipfun ,is_available,lsubmodel)
138 CALL hm_get_floatv('MAT_FScale' ,pscale ,is_available, lsubmodel, unitab)
139 CALL hm_get_intv ('Yr_fun' ,ifunce ,is_available,lsubmodel)
140 CALL hm_get_floatv('MAT_EFIB' ,einf ,is_available, lsubmodel, unitab)
141 CALL hm_get_floatv('MAT_C' ,ce ,is_available, lsubmodel, unitab)
142C-----------------------------------------------
143 ! Poisson's ratio check
144 IF (nu < zero .OR. nu >= half) THEN
145 CALL ancmsg(msgid=49,msgtype=msgerror,anmode=aninfo_blind_2,r1=nu,i1=id,c1=titr)
146 ENDIF
147c
148 IF(nrate > 100)THEN
149 CALL ancmsg(msgid=215, msgtype=msgerror, anmode=aninfo,i1=36,i2=id,c1=titr)
150 ELSEIF (nrate <= 0) THEN
151 CALL ancmsg(msgid=740, msgtype=msgerror, anmode=aninfo,i1=id,c1=titr)
152 ENDIF
153c
154 IF (ipfun == 0) THEN
155 pscale = zero
156 ELSEIF (pscale == zero) THEN
157 !units
158 CALL hm_get_floatv_dim('MAT_FScale' ,pscal_unit ,is_available, lsubmodel, unitab)
159 pscale = one * pscal_unit
160 ELSE
161 pscale = one /pscale
162 ENDIF
163c------------------------
164 IF (nrate > 0) THEN
165 DO j=1,nrate
166 CALL hm_get_int_array_index ('FUN_LOAD',ifunc(j),j,is_available,lsubmodel)
167 ENDDO
168 DO j=1,nrate
169 CALL hm_get_float_array_index ('SCALE_LOAD',yfac(j) ,j,is_available,lsubmodel,unitab)
170 IF(yfac(j) == zero) THEN
171 CALL hm_get_float_array_index_dim('SCALE_LOAD',yfac_unit(j),j,is_available,lsubmodel,unitab)
172 yfac(j)=one * yfac_unit(j)
173 ENDIF
174 ENDDO
175
176 rate(1:maxfunc) = zero
177 DO j=1,nrate
178 CALL hm_get_float_array_index ('STRAINRATE_LOAD',rate(j),j,is_available,lsubmodel,unitab)
179 ENDDO
180
181 DO i=1,nrate-1
182 IF (rate(i) >= rate(i+1)) THEN
183 CALL ancmsg(msgid=478, msgtype=msgerror, anmode=aninfo_blind_1,i1=id,c1=titr)
184 EXIT
185 ENDIF
186 ENDDO
187 DO i=1,nrate
188 IF (ifunc(i) == 0) THEN
189 CALL ancmsg(msgid=126, msgtype=msgerror, anmode=aninfo_blind_1,i1=id,c1=titr,i2=ifunc(i))
190 ENDIF
191 ENDDO
192 ENDIF
193c-------------------------
194 IF (nrate == 1) THEN
195 nfunc = 1
196 ismooth= 0
197 israt = 0 ! strain rate is always calculated for output
198 fcut = zero
199 vp = 0 !!! no plastic strain rate dependency with single static curve
200 ELSE ! NRATE > 1
201 israt = 1
202 IF (rate(1) == zero) THEN
203 nfunc = nrate
204 ELSE
205 nfunc = nrate+1
206 DO j=nrate,1,-1
207 ifunc(j+1) =ifunc(j)
208 rate(j+1) =rate(j)
209 yfac(j+1) =yfac(j)
210 ENDDO
211 rate(1) = zero
212 ENDIF
213c
214 IF (fcut == zero .or. vp == 1) THEN
215 fcut = 10000.0d0*unitab%FAC_T_WORK
216 END IF
217 ENDIF
218 israte = max(israte,israt)
219c
220 IF (nu == half) nu = zep499
221 mfunc = nfunc + 1
222 ifunc(mfunc) = ipfun
223c-----------------------------------------------
224 IF (fisokin > one .OR. fisokin < zero) THEN
225 CALL ancmsg(msgid=912, msgtype=msgerror, anmode=aninfo_blind_1,i1=id,c1='36',c2=titr)
226 END IF
227
228 IF (epsr1 == zero .AND. epsr2 == zero .AND. epsf == zero) THEN
229 IF (epsmax == zero) THEN
230 ifail = 0
231 ELSE
232 ifail = 1
233 END IF
234 ELSE
235 ifail = 2
236 ENDIF
237c IFAIL = 0 => no failure at all inside material
238c IFAIL = 1 => only failure vs max plastic strain
239c IFAIL = 2 => failure + damage vs principal tensile strain
240c
241 IF (ifail > 0) THEN
242 mtag%G_DMG = 1
243 mtag%L_DMG = 1
244 ENDIF
245 IF (epsmax== zero) epsmax= infinity
246 IF (epsr1 == zero) epsr1 = infinity
247 IF (epsr2 == zero) epsr2 = two*infinity
248 IF (epsf == zero) epsf = three*infinity
249c Limit max failure values
250 epsmax = min(epsmax ,infinity)
251 epsr1 = min(epsr1 ,infinity)
252 epsr2 = min(epsr2 ,two*infinity)
253 epsf = min(epsf ,three*infinity)
254c
255 IF (epsr1 /= zero .AND. epsr2 /= zero) THEN
256 IF (epsr1 >= epsr2) THEN
257 CALL ancmsg(msgid=480, msgtype=msgerror, anmode=aninfo_blind_1,i1=id,c1=titr)
258 ENDIF
259 ENDIF
260
261 IF(e <= zero)THEN
262 CALL ancmsg(msgid=276,msgtype=msgerror,anmode=aninfo,i1=36,i2=id,c1=titr)
263 e=zero
264 ENDIF
265
266C------------------------------
267 g = half*e/(one+nu)
268 c1= e/three/(one - two*nu)
269 soundsp = sqrt((c1 + four_over_3*g)/rho0)
270 yldcheck = 0 ! check if yld function decreases to zero (set in law36_upd.F)
271 opte = 0
272C------------------------------
273C------UPARAM STORAGE----------
274C------------------------------
275 uparam(1)= nfunc
276 uparam(2)= e
277 uparam(3)= e/(one - nu*nu)
278 uparam(4)= nu*uparam(3)
279 uparam(5)= g
280 uparam(6)= nu
281 DO j=1,nfunc
282 uparam(6 + j)= rate(j)
283 ENDDO
284 DO j=1,nfunc
285 uparam(nfunc + 6+j)= yfac(j)
286 ENDDO
287 uparam(2*nfunc + 7) = epsmax
288 uparam(2*nfunc + 8) = epsr1
289 uparam(2*nfunc + 9) = epsr2
290 uparam(2*nfunc + 10)= two*g
291 uparam(2*nfunc + 11)= three*g
292 uparam(2*nfunc + 12)= c1
293 uparam(2*nfunc + 13)= soundsp ! soundspeed solids
294 uparam(2*nfunc + 14)= fisokin
295 uparam(2*nfunc + 15)= epsf
296 IF (ipfun == 0) THEN
297 uparam(2*nfunc + 16) = 0
298 ELSE
299 uparam(2*nfunc + 16) = mfunc
300 ENDIF
301 uparam(2*nfunc + 17) = pscale
302c sound speed (shell)
303 uparam(2*nfunc + 18) = sqrt(e/(one - nu*nu)/rho0) ! soundspeed shells
304 uparam(2*nfunc + 19) = nu / (one-nu)
305 uparam(2*nfunc + 20) = three / (one+nu)
306 uparam(2*nfunc + 21) = one / (one-nu)
307c -----------------------
308 IF (ifunce > 0 ) opte = 1
309 mfunc = mfunc + 1
310 ifunc(mfunc) = ifunce
311 uparam(2*nfunc + 22) = mfunc
312 uparam(2*nfunc + 23) = opte
313 uparam(2*nfunc + 24) = einf
314 uparam(2*nfunc + 25) = ce
315 uparam(2*nfunc + 26) = vp
316 uparam(2*nfunc + 27) = ifail
317 uparam(2*nfunc + 28) = yldcheck
318 uparam(2*nfunc + 29) = ismooth
319c
320 nuparam = 2*nfunc + 29
321c--------------------------------
322 IF (rhor == zero) rhor=rho0
323 pm(1) = rhor
324 pm(89)= rho0
325 pm(27)= sqrt(e/rho0) ! Sound speed for beam elements
326c--------------------------------
327 parmat(1) = c1
328 parmat(2) = e
329 parmat(3) = nu
330 parmat(4) = israte
331 parmat(5) = fcut
332 !PARMAT(6) = NOFAIL
333 parmat(7) = epsr1
334 parmat(8) = epsr2
335 parmat(9) = epsf
336C Formulation for solid elements time step computation.
337 parmat(16) = 2
338 parmat(17) = two*g/(c1+four_over_3*g) ! == (1-2*nu)/(1-nu)
339C----------------
340 nuvar = 0
341 IF (vp == 1) THEN
342 nuvar = 3
343 ENDIF
344 nvartmp = 2+ nfunc
345C-----------------------
346 mtag%G_EPSD = 1
347 mtag%L_EPSD = 1
348 mtag%G_PLA = 1
349 mtag%L_PLA = 1
350 IF (fisokin /= zero) THEN
351 mtag%L_SIGB = 6
352 ENDIF
353C-----------------------
354 CALL init_mat_keyword(matparam,"ELASTO_PLASTIC")
355 CALL init_mat_keyword(matparam,"INCREMENTAL")
356 CALL init_mat_keyword(matparam,"LARGE_STRAIN")
357 ! Properties compatibility
358 CALL init_mat_keyword(matparam,"SHELL_ISOTROPIC")
359 CALL init_mat_keyword(matparam,"solid_isotropic")
360 CALL INIT_MAT_KEYWORD(MATPARAM,"sph")
361 CALL INIT_MAT_KEYWORD(MATPARAM,"beam_integrated")
362c
363 ! Material compatibility with /EOS option
364 CALL INIT_MAT_KEYWORD(MATPARAM,"eos")
365C-----------------------
366 WRITE(IOUT,1001) TRIM(TITR),ID,36
367 WRITE(IOUT,1000)
368 IF (IS_ENCRYPTED)THEN
369 WRITE(IOUT,'(5X,A,//)')'CONFIDENTIAL DATA'
370 ELSE
371 WRITE(IOUT,1002) RHO0
372 WRITE(IOUT,1100) E,NU,EPSMAX,EPSR1,EPSR2,EPSF,FISOKIN,ISMOOTH,FCUT,VP
373 WRITE(IOUT,1200)(IFUNC(J),YFAC(J),RATE(J),J=1,NFUNC)
374 WRITE(IOUT,1300) IPFUN,PSCALE, IFUNCE,EINF,CE
375 WRITE(IOUT,*)' '
376 ENDIF
377C-----------
378 RETURN
379C-----------
380 1000 FORMAT(
381 & 5X,' TABULATED ELASTIC PLASTIC LAW 36 ',/,
382 & 5X,' -------------------------------- ' ,//)
383 1001 FORMAT(/
384 & 5X,A,/,
385 & 5X,'MATERIAL NUMBER . . . . . . . . . . . . . .=',I10/,
386 & 5X,'MATERIAL LAW. . . . . . . . . . . . . . . .=',I10/)
387 1002 FORMAT(
388 & 5X,'INITIAL DENSITY . . . . . . . . . . . . . .=',1PG20.13/)
389 1100 FORMAT(
390 & 5X,'YOUNG MODULUS . . . . . . . . . . . . . . .=',1PG20.13/
391 & 5X,'POISSON RATIO . . . . . . . . . . . . . . .=',1PG20.13/
392 & 5X,'MAXIMUM PLASTIC STRAIN . . . . . . . . . .=',1PG20.13/
393 & 5X,'TENSION FAILURE STRAIN 1 . . . . . . . . .=',1PG20.13/
394 & 5X,'TENSION FAILURE STRAIN 2 . . . . . . . . .=',1PG20.13/
395 & 5X,'MAXIMUM TENSION FAILURE STRAIN . . . . . .=',1PG20.13/
396 & 5X,'ISO-KINEMATIC HARDENING FACTOR. . . . . . .=',1PG20.13/
397 & 5X,'SMOOTH STRAIN RATE OPTION . . . . . . . . .=',I10/
398 & 5X,' 0 -> NO SMOOTHING ',/,
399 & 5X,' 1 -> SMOOTH + LINEAR INTERPOLATION ',/,
400 & 5X,' 2 -> SMOOTH + LOG_N INTERPOLATION ',/
401 & 5X,'STRAIN RATE CUTTING FREQUENCY . . . . . . .=',1PG20.13/
402 & 5X,'PLASTIC STRAIN RATE DEPENDENCY FLAG . . . .=',I10/
403 & 5X,' FLAG_PL = 0 -> TOTAL SR DEPENDENCY ',/,
404 & 5X,' FLAG_PL = 1 -> PLASTIC SR DEPENDENCY ',/,
405 & 5X,'STRAIN RATE INTERPOLATION FLAG. . . . . . .=',I10/)
406 1200 FORMAT(
407 & 5X,'YIELD STRESS FUNCTION NUMBER. . . . . . . .=',I10/
408 & 5X,'YIELD SCALE FACTOR. . . . . . . . . . . . .=',1PG20.13/
409 & 5X,'STRAIN RATE . . . . . . . . . . . . . . . .=',1PG20.13)
410 1300 FORMAT(
411 & 5X,'PRESSURE DEPENDENT YIELD FUNCTION . . . . .=',I10/
412 & 5X,'PRESSURE SCALE FACTOR . . . . . . . . . . .=',1PG20.13/
413 & 5X,'YOUNG MODULUS SCALE FACTOR FUNCTION . . . .=',I10/
414 & 5X,'YOUNG MODULUS EINF. . . . . . . . . . . . .=',1PG20.13/
415 & 5X,'PARAMETER CE. . . . . . . . . . . . . . . .=',1PG20.13)
416c-----------
417 RETURN
418 END
hm_get_float_array_index
subroutine hm_get_float_array_index(name, rval, index, is_available, lsubmodel, unitab)
Definition hm_get_float_array_index.F:105
hm_get_float_array_index_dim
subroutine hm_get_float_array_index_dim(name, dim_fac, index, is_available, lsubmodel, unitab)
Definition hm_get_float_array_index_dim.F:47
hm_get_floatv
subroutine hm_get_floatv(name, rval, is_available, lsubmodel, unitab)
Definition hm_get_floatv.F:415
hm_get_floatv_dim
subroutine hm_get_floatv_dim(name, dim_fac, is_available, lsubmodel, unitab)
Definition hm_get_floatv_dim.F:159
hm_get_int_array_index
subroutine hm_get_int_array_index(name, ival, index, is_available, lsubmodel)
Definition hm_get_int_array_index.F:164
hm_get_intv
subroutine hm_get_intv(name, ival, is_available, lsubmodel)
Definition hm_get_intv.F:523
hm_option_is_encrypted
subroutine hm_option_is_encrypted(is_encrypted)
Definition hm_option_is_encrypted.F:284
hm_read_mat36
subroutine hm_read_mat36(uparam, maxuparam, nuparam, nuvar, nvartmp, ifunc, maxfunc, mfunc, parmat, unitab, id, mtag, titr, lsubmodel, pm, israte, matparam)
Definition hm_read_mat36.F:47
init_mat_keyword
subroutine init_mat_keyword(matparam, keyword)
Definition init_mat_keyword.F:150
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
elbuftag_mod
Definition elbuftag_mod.F:227
hm_option_read_mod
Definition hm_option_read_mod.F:408
message_mod
Definition message_mod.F:1249
names_and_titles_mod
Definition names_and_titles_mod.F:997
names_and_titles_mod::nchartitle
integer, parameter nchartitle
Definition names_and_titles_mod.F:1003
submodel_mod
Definition submodel_mod.F:712
unitab_mod
Definition unitab_mod.F:513
ancmsg
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:889
starter
program starter
Definition starter.F:39
elbuftag_mod::mlaw_tag_
Definition elbuftag_mod.F:232
submodel_mod::submodel_data
Definition submodel_mod.F:723