38
39
40
44
45
46
47#include "implicit_f.inc"
48
49
50
51#include "param_c.inc"
52
53
54
55 CHARACTER(LEN=NCHARTITLE) :: TITR
56 INTEGER MAT_ID,IOUT
57 INTEGER ,INTENT(IN) :: NFUNC
58 INTEGER ,INTENT(IN) :: NFUNCT
59 INTEGER NPC(*), FUNC_ID(NFUNCT) ,IPM(NPROPMI)
61 . uparam(*),pld(*),pm(npropm)
62 my_real ,
INTENT(IN) :: gama_inf
63 INTEGER, DIMENSION(NFUNC):: IFUNC
64
65
66
67
68
69 INTEGER N_NETWORK,N,K,ITEST,ICHECK,II,JJ,NSTART,IC1,IC2,NOGD,NDATA,NMULA,
70 . TAB,NMUAL,IFLAG,IS,I,INDX,LAWID, ITAG,IS_ENCRYPTED
72 . mu(5),al(5),lam_min,lam_max,strain_min,strain_max,
73 . d11,d22,d12,lam1,lam2,lam3,lam12,invd1,invd2
75 . e,nu,gs,rbulk, d,young,
76 . errtol,ave_slope,g,mu_max,mu_min,dx,lam,beta,
77 . amula(2),mual(10),dy
78 my_real ,
DIMENSION(:),
ALLOCATABLE :: stress,stretch
79 INTEGER , DIMENSION(:), ALLOCATABLE :: ITAB_ON_A,INDEX
80
81
82
83! IDENTIFICATION
84!====================================================================
85
86 IF(ifunc(2) == 0 ) RETURN
87
88 nmual = 2*uparam(18)
89 iflag = nint(uparam(13))
90 icheck = nint(uparam(19))
91 nstart = nint(uparam(20))
92 is_encrypted = nint(uparam(21))
93 errtol = fiveem3
94
95 ic1 = npc(ifunc(2))
96 ic2 = npc(ifunc(2)+1)
97
98 nogd=(ic2-ic1)/2
99 ndata=nogd
100
101 ALLOCATE (stretch(nogd))
102 ALLOCATE (stress(nogd))
103
104 ave_slope = zero
105 jj=0
106 stretch=zero
107 stress=zero
108 g=zero
109 rbulk=zero
110 gs=zero
111 lam_max= zero
112 lam_min= zero
113 itag = 0
114 is = ifunc(2)
115 DO ii = ic1,ic2-2,2
116 jj=jj+1
117 stretch(jj) = pld(ii) + one
118 stress(jj) = pld(ii+1)
119 IF(pld(ii) < zero .AND. stress(jj) < zero ) itag = 1
120 IF(pld(ii) <= -one) THEN
122 . msgtype=msgerror,
123 . anmode=aninfo,
124 . i1=mat_id,
125 . c1=titr,
126 . i2=func_id(is))
128 ENDIF
129 ENDDO
130
131 DO jj =1,nogd - 1
132 dx = stretch(jj+1) - stretch(jj)
133 dy = stress(jj + 1) - stress(jj)
134 IF(dx <= zero .OR. dy <= zero) THEN
136 . msgtype=msgerror,
137 . anmode=aninfo,
138 . i1=mat_id,
139 . c1=titr,
140 . i2=func_id(is))
141 ENDIF
142 ENDDO
143 ave_slope = zero
144
145 IF (iflag == 2 .AND. nmual /= 4) THEN
147 . msgtype=msgerror,
148 . anmode=aninfo_blind_1,
149 . i1=mat_id,
150 . c1=titr,
151 . i2=func_id(is))
152 ENDIF
153 DO ii=1,5
154 mual(2*ii-1)=five
155 mual(2*ii)=two
156 ENDDO
157
158 WRITE(iout,277)
159 WRITE(iout,270)trim(titr),mat_id
160 IF(is_encrypted > 0) WRITE(iout,'(5X,A,//)')'CONFIDENTIAL DATA'
161
162 IF (iflag == 1 .OR. iflag == 2) THEN
163 lawid = iflag
164 CALL law69_nlsqf(lawid,stretch,stress,10,nmual,nogd,mual,
165 $ icheck, nstart, errtol,mat_id,titr,is_encrypted
166 ELSEIF(iflag == -1) THEN
167 lawid = 2
168 IF(itag == 0) icheck = abs(icheck)
170 $ icheck, nstart, errtol,mat_id,titr,is_encrypted)
171 uparam(18) = nmual / 2
172 ENDIF
173 DEALLOCATE (stretch)
174 DEALLOCATE (stress)
175 IF(gama_inf < one) THEN
176 DO ii=1,5
177 uparam(ii)=mual(2*ii-1)/gama_inf
178 uparam(ii+5)=mual(2*ii)
179 gs=gs + mual(2*ii-1)*mual(2*ii)
180 ENDDO
181 IF(is_encrypted == 0) WRITE(iout,100) mat_id
182 ELSE
183 DO ii=1,5
184 uparam(ii)=mual(2*ii-1)
185 uparam(ii+5)=mual(2*ii)
186 gs=gs+mual(2*ii-1)*mual(2*ii)
187 ENDDO
188 ENDIF
189 g=gs/two
190 nu = uparam(14)
191 rbulk=gs*(one+uparam(14))
192 & /
max(em30,three*(one-two*uparam(14)))
193 uparam(11)=rbulk
194 uparam(17)=gs
195
196 young = gs*(one + nu)
197 pm(20) = young
198 pm(21) = nu
199 pm(22) = gs
200 pm(24) = young/(one - nu**2)
201 pm(32) = gs
202 pm(100) = rbulk
203
204 ipm(252)= 2
205 pm(105) = gs/(rbulk + two_third*g)
206
207
208 IF(is_encrypted == 0) THEN
209 IF (lawid /= 2) THEN
210 WRITE(iout,278)
211 & uparam(1),uparam(2),uparam(3),uparam(4),
212 & uparam(5),uparam(6),uparam(7),uparam(8),
213 & uparam(9),uparam(10),g,rbulk
214 ELSE
215 WRITE(iout,279)
216 & uparam(1),uparam(2),uparam(6),uparam(7),
217 & half*uparam(1),-half*uparam(2),g,rbulk
218 ENDIF
219 ENDIF
220
221
222
223
224
225
226 ndata = 10000
227 lam_min= em3
228 lam_max= ten
229
230 DO i=1,5
231 mu(i) = uparam(i)
232 al(i) = uparam(5+i)
233 ENDDO
234
235 ALLOCATE (stretch(ndata))
236 ALLOCATE (stress(ndata))
237 ALLOCATE (itab_on_a(ndata))
238 ALLOCATE (index(ndata))
239 stretch(1)=lam_min
240 itab_on_a = 0
241 DO k= 2,ndata
242 stretch(k)=stretch(k-1) + em3
243 ENDDO
244 stress=zero
245
246 WRITE(iout,1000)mat_id
247
248 indx =0
249 DO i = 1, ndata
250 d11= zero
251 d22= zero
252 d12= zero
253 lam1 =stretch(i)
254 lam2 = one/sqrt(lam1)
255 lam3 = lam2
256
257 DO k=1,5
258 lam12 = (lam1*lam2)**(-al(k))
259 d11=d11+ al(k)*mu(k) * (lam1**al(k) + lam12 )
260 d22=d22+ al(k)*mu(k) * (lam2**al(k) + lam12 )
261 d12=d12+ al(k)*mu(k) * lam12
262 ENDDO
263
264 invd1 = d11 + d22
265 invd2 = d11*d22 - d12**2
266 IF (invd1 > 0 .AND. invd2 > 0) THEN
267 indx = indx +1
268 itab_on_a(indx) = 1
269 index(indx) = i
270 ENDIF
271 ENDDO
272 IF(indx > 0 .AND. indx < ndata) THEN
273 WRITE(iout,1010)
274 i = index(1) - 1
275 IF(i > 1) THEN
276 strain_min = stretch(i) - one
277 WRITE(iout,1100)strain_min
278 ENDIF
279 i = index(indx) + 1
280 IF(i <= ndata)THEN
281 strain_max = stretch(i) - one
282 WRITE(iout,1200)strain_max
283 ENDIF
284 ENDIF
285
286 indx =0
287 DO i = 1, ndata
288 d11= zero
289 d22= zero
290 d12= zero
291 lam1 =stretch(i)
292 lam2 = lam1
293 lam3 = one/lam1/lam1
294
295 DO k=1,5
296 lam12 = (lam1*lam2)**(-al(k))
297 d11=d11+ al(k)*mu(k) * (lam1**al(k) + lam12 )
298 d22=d22+ al(k)*mu(k) * (lam2**al(k) + lam12 )
299 d12=d12+ al(k)*mu(k) * lam12
300 ENDDO
301
302 invd1 = d11 + d22
303 invd2 = d11*d22 - d12**2
304 IF (invd1 > 0 .AND. invd2 > 0) THEN
305 indx = indx + 1
306 itab_on_a(indx) = 1
307 index(indx) = i
308 ENDIF
309 ENDDO
310 IF(indx > 0 .AND. indx < ndata) THEN
311 WRITE(iout,1020)
312 i = index(1) - 1
313 IF(i > 1) THEN
314 strain_min = stretch(i) - one
315 WRITE(iout,1100)strain_min
316 ENDIF
317 i = index(indx) + 1
318 IF(i <= ndata)THEN
319 strain_max = stretch(i) - one
320 WRITE(iout,1200)strain_max
321 ENDIF
322 ENDIF
323
324 indx =0
325 DO i = 1, ndata
326 d11= zero
327 d22= zero
328 d12= zero
329 lam1 =stretch(i)
330 lam2 = one
331 lam3 = one/lam1
332
333 DO k=1,5
334 lam12 = (lam1*lam2)**(-al(k))
335 d11=d11+ al(k)*mu(k) * (lam1**al(k) + lam12 )
336 d22=d22+ al(k)*mu(k) * (lam2**al(k) + lam12 )
337 d12=d12+ al(k)*mu(k) * lam12
338 ENDDO
339
340 invd1 = d11 + d22
341 invd2 = d11*d22 - d12**2
342 IF (invd1 > 0 .AND. invd2 > 0) THEN
343 indx = indx +1
344 itab_on_a(indx) = 1
345 index(indx) = i
346 ENDIF
347 ENDDO
348 IF(indx > 0 .AND. indx < ndata) THEN
349 WRITE(iout,1030)
350 i = index(1) - 1
351 IF(i > 1) THEN
352 strain_min = stretch(i) - one
353 WRITE(iout,1100)strain_min
354 ENDIF
355 i = index(indx) + 1
356 IF(i <= ndata)THEN
357 strain_max = stretch(i) - one
358 WRITE(iout,1200)strain_max
359 ENDIF
360 ENDIF
361 DEALLOCATE (stretch)
362 DEALLOCATE (stress)
363 DEALLOCATE (itab_on_a)
364 DEALLOCATE (index)
365
366 RETURN
367 100 FORMAT
368 & (//5x, 'MODIFIED FITTED MATERIAL RIGIDITY ' ,/,
369 & 5x, ' ---------------------------', /,
370 & 5x, 'MATERIAL LAW = LAW69 ',/,
371 & 5x, 'MATERIAL NUMBER =',i10,//)
372 277 FORMAT
373 & (//5x, 'FITTED PARAMETERS FOR HYPERELASTIC_MATERIAL LAW69' ,/,
374 & 5x, ' ----------------------------------------')
375 270 FORMAT(
376 & 5x,a,/,
377 & 5x, 'MATERIAL NUMBER =',i10,//)
378 278 FORMAT(
379 & 5x, 'OGDEN LAW PARAMETERS:',/
380 & ,5x, 'MU1 . . . . . . . . . . . . . . . . . .=',1pg20.13/
381 & ,5x, 'MU2 . . . . . . . . . . . . . . . . . .=',1pg20.13/
382 & ,5x, 'MU3 . . . . . . . . . . . . . . . . . .=',1pg20.13/
383 & ,5x, 'MU4 . . . . . . . . . . . . . . . . . .=',1pg20.13/
384 & ,5x, 'MU5 . . . . . . . . . . . . . . . . . .=',1pg20.13//
385 & ,5x, 'AL1 . . . . . . . . . . . . . . . . . .=',1pg20.13/
386 & ,5x, 'AL2 . . . . . . . . . . . . . . . . . .=',1pg20.13/
387 & ,5x, 'AL3 . . . . . . . . . . . . . . . . . .=',1pg20.13/
388 & ,5x, 'AL4 . . . . . . . . . . . . . . . . . .=',1pg20.13/
389 & ,5x, 'AL5 . . . . . . . . . . . . . . . . . .=',1pg20.13//
390 & ,5x, 'INITIAL SHEAR MODULUS . . . . . . . . .=',1pg20.13/
391 & ,5x, 'BULK MODULUS. . . . . . . . . . . . . .=',1pg20.13//)
392 279 FORMAT(
393 & 5x, 'MOONEY-RIVLIN LAW PARAMETERS:',/
394 & ,5x, 'MU1 . . . . . . . . . . . . . . . . . .=',1pg20.13/
395 & ,5x, 'MU2 . . . . . . . . . . . . . . . . . .=',1pg20.13/
396 & ,5x, 'AL1 . . . . . . . . . . . . . . . . . .=',1pg20.13/
397 & ,5x, 'AL2 . . . . . . . . . . . . . . . . . .=',1pg20.13//
398 & ,5x, 'C10 . . . . . . . . . . . . . . . . . .=',1pg20.13/
399 & ,5x, 'C01 . . . . . . . . . . . . . . . . . .=',1pg20.13/
400 & ,5x, 'INITIAL SHEAR MODULUS . . . . . . . . .=',1pg20.13/
401 & ,5x, 'BULK MODULUS. . . . . . . . . . . . . .=',1pg20.13//)
402 1000 FORMAT
403 & (//5x, 'CHECK THE DRUCKER PRAGER STABILITY CONDITIONS ' ,/,
404 & 5x, ' -----------------------------------------------', /,
405 & 5x, 'MATERIAL LAW = OGDEN (LAW69) ',/,
406 & 5x, 'MATERIAL NUMBER =',i10,//)
407 1010 FORMAT(
408 & 7x,'TEST TYPE = UNIXIAL ')
409 1020 FORMAT(//,
410 & 7x,'TEST TYPE = BIAXIAL ')
411 1030 FORMAT(//,
412 & 7x,'TEST TYPE = PLANAR (SHEAR)')
413 1100 FORMAT(
414 & 8x,'COMPRESSION: UNSTABLE AT A NOMINAL STRAIN LESS THAN ',1pg20.13)
415 1200 FORMAT(
416 & 8x,'TENSION: UNSTABLE AT A NOMINAL STRAIN LARGER THAN ',1pg20.13//)
417
subroutine law69_nlsqf_auto(lawid, stretch, y, ma, nmual, npt, mual, icheck, nstart, errtol, id, titr, is_encrypted)
integer, parameter nchartitle
subroutine law69_nlsqf(lawid, stretch, y, ma, nmual, npt, mual, icheck, nstart, errtol, id, titr, is_encrypted)
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)