45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
68 USE matparam_def_mod
70
71
72
73#include "implicit_f.inc"
74
75
76
77#include "units_c.inc"
78#include "param_c.inc"
79
80
81
82 TYPE (UNIT_TYPE_),INTENT(IN) :: UNITAB
83 my_real,
INTENT(INOUT) :: pm(npropm),parmat(100),uparam(maxuparam)
84 INTEGER, INTENT(INOUT) :: IPM(NPROPMI),ISRATE,IFUNC(MAXFUNC),NFUNC,
85 . MAXFUNC,MAXUPARAM,NUPARAM,NUVAR,IMATVIS ,
86 . NVARTMP
87 TYPE(MLAW_TAG_),INTENT(INOUT) :: MTAG
88 INTEGER,INTENT(IN) :: MAT_ID
89 CHARACTER(LEN=NCHARTITLE) ,INTENT(IN) :: TITR
90 TYPE(SUBMODEL_DATA),INTENT(IN) :: LSUBMODEL(*)
91 TYPE(MATPARAM_STRUCT_) ,INTENT(INOUT) :: MATPARAM
92
93
94
95 INTEGER J,NRATE,I,ILAW,VP
97 . e11,e22,e33,nu12,nu23,nu13,g12,g13,g23,qr1,qr2,cr1,cr2,
98 . sigy,r11,r22,r33,r12,r13,r23,a1,a2,a3,hh,ff,gg,ll,mm,nn,
99 . d11,d22,d33,d12,d13,d23,a11,a22,a12,c11,c22,c33,c12,c13,
100 . c23,nu21,nu31,nu32,detc,fac,yfac(100),rate(100),dmin,dmax,
101 . yscale_unit,rho0,rhor,asrate
102 LOGICAL :: IS_AVAILABLE,IS_ENCRYPTED
103
104
105
106 is_encrypted = .false.
107 is_available = .false.
108 ilaw = 93
109
110
112
113
114
115 CALL hm_get_floatv(
'MAT_RHO',rho0 ,is_available, lsubmodel, unitab)
116
118 CALL hm_get_floatv(
'LAW93_E22' ,e22 ,is_available, lsubmodel, unitab)
119 CALL hm_get_floatv(
'LAW93_E33' ,e33 ,is_available, lsubmodel, unitab
120 CALL hm_get_floatv(
'LAW93_G12' ,g12 ,is_available, lsubmodel, unitab)
121 CALL hm_get_floatv(
'LAW93_Nu12',nu12 ,is_available, lsubmodel, unitab)
122
123 CALL hm_get_floatv(
'LAW93_G13' ,g13 ,is_available, lsubmodel, unitab)
124 CALL hm_get_floatv(
'LAW93_G23' ,g23 ,is_available, lsubmodel, unitab)
125 CALL hm_get_floatv(
'LAW93_Nu13',nu13 ,is_available, lsubmodel, unitab)
126 CALL hm_get_floatv(
'LAW93_Nu23',nu23 ,is_available, lsubmodel, unitab)
127
128 CALL hm_get_intv (
'LAW93_NL' ,nrate ,is_available, lsubmodel)
129 CALL hm_get_floatv(
'FCUT' ,asrate ,is_available, lsubmodel, unitab)
130 CALL hm_get_intv (
'VP' ,vp ,is_available, lsubmodel
131
132 IF (nrate > 0) THEN
133 DO i=1,nrate
137 IF (yfac(i) == zero) THEN
138 CALL hm_get_floatv_dim(
'LAW93_arr2' ,yscale_unit ,is_available, lsubmodel, unitab)
139 yfac(i) = one * yscale_unit
140 ENDIF
141 ENDDO
142 ENDIF
143
144 CALL hm_get_floatv(
'LAW93_Sigma_y',sigy ,is_available, lsubmodel, unitab)
145 CALL hm_get_floatv(
'LAW93_QR1' ,qr1 ,is_available, lsubmodel, unitab)
146 CALL hm_get_floatv(
'LAW93_CR1' ,cr1 ,is_available, lsubmodel, unitab)
148 CALL hm_get_floatv(
'LAW93_CR2' ,cr2 ,is_available, lsubmodel, unitab)
149
150 CALL hm_get_floatv(
'LAW93_R11' ,r11 ,is_available, lsubmodel, unitab)
151 CALL hm_get_floatv('law93_r22
' ,R22 ,IS_AVAILABLE, LSUBMODEL, UNITAB)
152 CALL HM_GET_FLOATV('law93_r12' ,R12 ,IS_AVAILABLE, LSUBMODEL, UNITAB)
153
154 CALL HM_GET_FLOATV('law93_r33' ,R33 ,IS_AVAILABLE, LSUBMODEL, UNITAB)
155 CALL HM_GET_FLOATV('law93_r13' ,R13 ,IS_AVAILABLE, LSUBMODEL, UNITAB)
156 CALL HM_GET_FLOATV('law93_r23' ,R23 ,IS_AVAILABLE, LSUBMODEL, UNITAB)
157
158 !========== DEFAULT VALUES=============!
159
160 ! Default value for functions
161 NFUNC = NRATE
162 IF (NRATE > 1) THEN
163 IF (RATE(1) == ZERO) THEN
164 NFUNC = NRATE
165 ELSE
166 NFUNC = NRATE + 1
167 DO J = NRATE,1,-1
168 IFUNC(J+1) = IFUNC(J)
169 RATE(J+1) = RATE(J)
170 YFAC(J+1) = YFAC(J)
171 ENDDO
172 RATE(1) = ZERO
173 ENDIF
174 ENDIF
175
176 ! Yield stresses
177 IF(SIGY == ZERO) SIGY = INFINITY
178 IF(R11 == ZERO) R11 = ONE
179 IF(R22 == ZERO) R22 = ONE
180 IF(R33 == ZERO) R33 = ONE
181 IF(R12 == ZERO) R12 = ONE
182 IF(R23 == ZERO) R23 = ONE
183 IF(R13 == ZERO) R13 = ONE
184
185 ! Young modulus
186 IF (E22 == ZERO) E22 = E11
187 IF (E33 == ZERO) E33 = E22
188 ! Shear modulus
189 IF (G13 == ZERO) G13 = G12
190 IF (G23 == ZERO) G23 = G12
191 ! Remaining Poisson's ratio
192 nu21 = nu12*e22/e11
193 nu31 = nu13*e33/e11
194 nu32 = nu23*e33/e22
195
196
197 if(nu12*nu21 >= one ) then
199 . msgtype=msgerror,
200 . anmode=aninfo_blind_2,
201 . i1=mat_id ,
202 . c1=titr)
203 else if(nu13*nu31 >= one ) then
205 . msgtype=msgerror,
206 . anmode=aninfo_blind_2,
207 . i1=mat_id ,
208 . c1=titr)
209 else if(nu23*nu32 >= one ) then
211 . msgtype=msgerror,
212 . anmode=aninfo_blind_2,
213 . i1=mat_id ,
214 . c1=titr)
215 endif
216
217
218 a1 = one/r11/r11
219 a2 = one/r22/r22
220 a3 = one/r33/r33
221 ff = half*(a2 + a3 - a1)
222 gg = half*(a3 + a1 - a2)
223 hh = half*(a1 + a2 - a3
224 ll = three_half/r23/r23
225 mm = three_half/r13/r13
226 nn = three_half/r12/r12
227
228
229 fac = one/(one - nu12*nu21)
230 a11 = e11*fac
231 a12 = nu21*a11
232 a22 = e22*fac
233
234 c11 = one/e11
235 c22 = one/e22
236 c33 = one/e33
237 c12 =-nu12/e11
238 c13 =-nu31/e33
239 c23 =-nu23/e22
240
241 detc= c11*c22*c33-c11*c23*c23-c12*c12*c33+c12*c13*c23
242 + +c13*c12*c23-c13*c22*c13
243 IF(detc<=zero) THEN
245 . msgtype=msgerror,
246 . anmode=aninfo,
247 . i1=mat_id,
248 . c1=titr)
249 ENDIF
250
251 d11 = (c22*c33-c23*c23)/detc
252 d12 =-(c12*c33-c13*c23)/detc
253 d13 = (c12*c23-c13*c22)/detc
254 d22 = (c11*c33-c13*c13)/detc
255 d23 =-(c11*c23-c13*c12)/detc
256 d33 = (c11*c22-c12*c12)/detc
257 dmin =
min(d11*d22 -d12**2, d11*d33 - d13**2, d22*d33 - d23**2 )
258 dmax =
max(d11,d22,d33)
259
260
261 IF (nfunc > 1) THEN
262 israte = 1
263 IF (vp == 0) vp = 2
264 IF (vp == 1) THEN
265 asrate = 1.0d4*unitab%FAC_T_WORK
266 ELSE
267 IF (asrate == zero) asrate = 1.0d4*unitab%FAC_T_WORK
268 ENDIF
269 ELSE
270 israte = 0
271 asrate = zero
272 ENDIF
273
274
275 rhor = zero
276 pm(1) = rhor
277 pm(89) = rho0
278
279
280 parmat(1) =
max(a11,a22,d11,d22,d33)
281 parmat(2) =
max(e11,e22,e33)
282 parmat(3) =
max(nu12,nu13,nu23)
283 parmat(4) = israte
284 parmat(5) = asrate
285 parmat(16) = 1
286 parmat(17) = dmin/dmax/dmax
287
288
289 mtag%G_PLA = 1
290 mtag%L_PLA = 1
291 mtag%G_SEQ = 1
292 mtag%L_SEQ = 1
293 mtag%G_EPSD = 1
294 mtag%L_EPSD = 1
295
296
302
303
307
308
309 imatvis = 0
310
311
312 IF ((nrate > 1).AND.(vp /= 2)) THEN
313 nuvar = 1
314 ELSE
315 nuvar = 0
316 ENDIF
317
318
319 nuparam = 30 + 2*nfunc
320 nvartmp = nfunc
321
322
323
324 uparam(1) = a11
325 uparam(2) = a22
326 uparam(3) = a12
327 uparam(4) = d11
328 uparam(5) = d12
329 uparam(6) = d13
330 uparam(7) = d22
331 uparam(8) = d23
332 uparam(9) = d33
333 uparam(10) = g12
334 uparam(11) = g13
335 uparam(12) = g23
336 uparam(13) = e11
337 uparam(14) = e22
338 uparam(15) = e33
339 uparam(16) = nu12
340 uparam(17) = nu13
341 uparam(18) = nu23
342
343 uparam(19) = ff
344 uparam(20) = gg
345 uparam(21) = hh
346 uparam(22) = ll
347 uparam(23) = mm
348 uparam(24) = nn
349
350 uparam(25) = sigy
351 uparam(26) = qr1
352 uparam(27) = cr1
353 uparam(28) = qr2
354 uparam(29) = cr2
355
356 uparam(30) = vp
357
358 IF (nfunc > 0) THEN
359 DO j=1,nfunc
360 uparam(30 + j) = rate(j)
361 uparam(30 + nfunc + j) = yfac(j)
362 ENDDO
363 ENDIF
364
365
366
367
368 WRITE(iout,1001) trim(titr),mat_id,ilaw
369 WRITE(iout,1000)
370 IF(is_encrypted)THEN
371 WRITE(iout,'(5X,A,//)')'CONFIDENTIAL DATA'
372 ELSE
373 WRITE(iout,1002) rho0
374 WRITE(iout,1300) e11,e22,e33,g12,g13,g23,nu12,nu13,nu23
375 IF (nrate == 0) THEN
376 WRITE(iout,1450)
377 WRITE(iout,1400) sigy,qr1,cr1,qr2,cr2
378 ELSE
379 WRITE(iout,1550)
380 DO j=1,nfunc
381 WRITE(iout,1500) ifunc(j),yfac(j),rate(j)
382 ENDDO
383 IF (nrate > 1) THEN
384 WRITE(iout,1575) asrate,vp
385 ENDIF
386 ENDIF
387 WRITE(iout,1600) r11,r22,r33,r12,r13,r23
388 ENDIF
389
390 1000 FORMAT(
391 & 5x,' ORTHOTROPIC ELASTIC + HILL CRITERION '/,
392 & 5x,' ------------------------------------ '//)
393 1001 FORMAT(
394 & 5x,a,/,
395 & 5x,'MATERIAL NUMBER . . . . . . . . . . .=',i10/,
396 & 5x,'MATERIAL LAW. . . . . . . . . . . . .=',i10/)
397 1002 FORMAT(
398 & 5x,'INITIAL DENSITY . . . . . . . . . . .=',1pg20.13/)
399 1300 FORMAT(
400 & 5x,'YOUNG MODULUS IN 11 DIRECTION . . . .=',1pg20.13/,
401 & 5x,'YOUNG MODULUS IN 22 DIRECTION . . . .=',1pg20.13/,
402 & 5x,'YOUNG MODULUS IN 33 DIRECTION . . . .=',1pg20.13/,
403 & 5x,'SHEAR MODULUS IN 12 DIRECTION . . . .=',1pg20.13/,
404 & 5x,'SHEAR MODULUS IN 13 DIRECTION . . . .=',1pg20.13/,
405 & 5x,'SHEAR MODULUS IN 23 DIRECTION . . . .=',1pg20.13/,
406 & 5x,'POISSON RATIO 12. . . . . . . . . . .=',1pg20.13/,
407 & 5x,'POISSON RATIO 13. . . . . . . . . . .=',1pg20.13/,
408 & 5x,'POISSON RATIO 23. . . . . . . . . . .=',1pg20.13//)
409 1550 FORMAT(
410 & 5x,'--------------------------------------'/,
411 & 5x,'TABULATED YIELD STRESS '/,
412 & 5x,'--------------------------------------'//)
413 1500 FORMAT(
414 & 5x,'YIELD STRESS FUNCTION NUMBER. . . . .=',i10/,
415 & 5x,'YIELD SCALE FACTOR. . . . . . . . . .=',1pg20.13/,
416 & 5x,'STRAIN RATE . . . . . . . . . . . . .=',1pg20.13/)
417 1575 FORMAT(
418 & 5x,'STRAIN RATE CUTTING FREQUENCY . . . .=',1pg20.13/
419 & 5x,'STRAIN RATE CHOICE FLAG . . . . . . .=',i10/
420 & 5x,' VP=1 EQUIVALENT PLASTIC STRAIN RATE'/
421 & 5x,' VP=2 TOTAL STRAIN RATE (DEFAULT)'/
422 & 5x,' VP=3 DEVIATORIC STRAIN RATE'/)
423 1450 FORMAT(
424 & 5x,'--------------------------------------'/,
425 & 5x,'CONTINUOUS YIELD STRESS '/,
426 & 5x,'--------------------------------------'//)
427 1400 FORMAT(
428 & 5x,'INITIAL YIELD STRESS. . . . . . . . .=',1pg20.13/,
429 & 5x,'PARAMETER QR1 OF HARDENING . . . . .=',1pg20.13/,
430 & 5x,'PARAMETER CR1 OF HARDENING . . . . .=',1pg20.13/,
431 & 5x,'PARAMETER QR2 OF HARDENING . . . . .=',1pg20.13/,
432 & 5x,'PARAMETER CR2 OF HARDENING . . . . .=',1pg20.13/,
433 & 5x,'REFERENCE STRAIN. . . . . . . . . . .=',1pg20.13//)
434 1600 FORMAT(
435 & 5x,'RATIO YIELD PARAMETER R11 . . . . . .=',1pg20.13/,
436 & 5x,'RATIO YIELD PARAMETER R22 . . . . . .=',1pg20.13/,
437 & 5x,'RATIO YIELD PARAMETER R33 . . . . . .=',1pg20.13/,
438 & 5x,'RATIO YIELD PARAMETER R12 . . . . . .=',1pg20.13/,
439 & 5x,'RATIO YIELD PARAMETER R13 . . . . . .=',1pg20.13/,
440 & 5x,'RATIO YIELD PARAMETER R23 . . . . . .=',1pg20.13/)
441 RETURN
subroutine hm_get_float_array_index(name, rval, index, is_available, lsubmodel, unitab)
subroutine hm_get_floatv(name, rval, is_available, lsubmodel, unitab)
subroutine hm_get_floatv_dim(name, dim_fac, is_available, lsubmodel, unitab)
subroutine hm_get_int_array_index(name, ival, index, is_available, lsubmodel)
subroutine hm_get_intv(name, ival, is_available, lsubmodel)
subroutine hm_option_is_encrypted(is_encrypted)
subroutine init_mat_keyword(matparam, keyword)
integer, parameter nchartitle
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)