59
60
61
62 USE elbufdef_mod
66 USE matparam_def_mod
68 use glob_therm_mod
69 use s20temp_mod
70
71
72
73#include "implicit_f.inc"
74
75
76
77#include "mvsiz_p.inc"
78
79
80
81#include "com04_c.inc"
82#include "param_c.inc"
83#include "scr12_c.inc"
84#include "scr17_c.inc"
85#include "scry_c.inc"
86#include "vect01_c.inc"
87
88
89
90 INTEGER IXS(NIXS,*), IPARG(*),IPARTS(*),IGEO(NPROPGI,*),
91 . IXS16(8,*), IPART(LIPART1,*),IPM(NPROPMI,*), PTSOL(*),
92 . NPF(*),STRSGLOB(*),STRAGLOB(*),FAIL_INI(*),PERTURB(NPERTURB)
93 INTEGER NEL,NSIGI,IUSER,NSIGS
95 . mas(*), pm(npropm,*), x(*), geo(npropg,*),
96 . veul(lveul,*), dtelem(*),sigi(nsigs,*),skew(lskew,*),stifn(*),
97 . partsav(20,*), v(*), mss(8,*), mssx(12,*), sigsp(nsigi, *),
98 . volnod(*), bvolnod(*), vns(8,*), bns(8,*),
99 . vnsx(12,*), bnsx(12,*),bufmat(*),rnoise(nperturb,*),
100 . mcp(*), mcps(8,*),mcpsx(12,*), temp(*), tf(*)
101 TYPE(ELBUF_STRUCT_), TARGET :: ELBUF_STR
102 INTEGER,INTENT(IN) :: ILOADP(SIZLOADP,*)
103 my_real,
INTENT(IN) :: facload(lfacload,*)
104 TYPE(DETONATORS_STRUCT_)::DETONATORS
105 TYPE(t_ale_connectivity), INTENT(INOUT) ::
106 TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(INOUT) :: MAT_PARAM
107 type (glob_therm_) ,intent(in) :: glob_therm
108
109
110
111 INTEGER NF1,IBID,I,IGTYP,IP,NF2,NPTR,NPTS,NPTT,NLAY,IL,IR,IS,IT,
112 . N, NUVAR,IINT, NCC,IDEF,JHBE,IPID1,L_PLA,L_SIGB
113 INTEGER NC(MVSIZ,16),MAT(MVSIZ), PID(MVSIZ), NGL(MVSIZ),RBID(1)
114 INTEGER ,PARAMETER :: NPE=16
115 CHARACTER(LEN=NCHARTITLE)::TITR1
117 . bid, fv,wi,aa,bb
119 . mass(mvsiz),
120 . volp(mvsiz,5), sti(mvsiz),deltax(mvsiz),deltax2(mvsiz),
121 . xx(mvsiz,16), yy(mvsiz,16), zz(mvsiz,16),
122 . vx(mvsiz,16), vy(mvsiz,16), vz(mvsiz,16),
123 . px(mvsiz,16), py(mvsiz,16), pz(mvsiz,16),
124 . rx(mvsiz),ry(mvsiz),rz(mvsiz),
125 . sx(mvsiz),sy(mvsiz),sz(mvsiz),volg(mvsiz),
126 . tx(mvsiz),ty(mvsiz),tz(mvsiz),ul(mvsiz,16),
127 . ni(mvsiz,16),dnidr(mvsiz,16),dnids(mvsiz,16),dnidt(mvsiz,16),
128 . dtx(mvsiz),stin(mvsiz,16), rhocp(mvsiz),temp0(mvsiz), aire(mvsiz)
130 TYPE(L_BUFEL_) ,POINTER :: LBUF
131 TYPE(G_BUFEL_) ,POINTER :: GBUF
132 TYPE(BUF_MAT_) ,POINTER :: MBUF
133
135 . w_gauss(9,9),a_gauss(9,9),w_lobatto(9,9),a_lobatto(9,9),
136 . w_newton(9,9),a_newton(9,9)
137
138 DATA w_gauss /
139 1 2. ,0. ,0. ,
140 1 0. ,0. ,0. ,
141 1 0. ,0. ,0. ,
142 2 1. ,1. ,0. ,
143 2 0. ,0. ,0. ,
144 2 0. ,0. ,0. ,
145 3 0.555555555555556,0.888888888888889,0.555555555555556,
146 3 0. ,0. ,0. ,
147 3 0. ,0. ,0. ,
148 4 0.347854845137454,0.652145154862546,0.652145154862546,
149 4 0.347854845137454,0. ,0. ,
150 4 0. ,0. ,0. ,
151 5 0.236926885056189,0.478628670499366,0.568888888888889,
152 5 0.478628670499366,0.236926885056189,0. ,
153 5 0. ,0. ,0. ,
154 6 0.171324492379170,0.360761573048139,0.467913934572691,
155 6 0.467913934572691,0.360761573048139,0.171324492379170,
156 6 0. ,0. ,0. ,
157 7 0.129484966168870,0.279705391489277,0.381830050505119,
158 7 0.417959183673469,0.381830050505119,0.279705391489277,
159 7 0.129484966168870,0. ,0. ,
160 8 0.101228536290376,0.222381034453374,0.313706645877887,
161 8 0.362683783378362,0.362683783378362,0.313706645877887,
162 8 0.222381034453374,0.101228536290376,0. ,
163 9 0.081274388361574,0.180648160694857,0.260610696402935,
164 9 0.312347077040003,0.330239355001260,0.312347077040003,
165 9 0.260610696402935,0.180648160694857,0.081274388361574/
166 DATA a_gauss /
167 1 0. ,0. ,0. ,
168 1 0. ,0. ,0. ,
169 1 0. ,0. ,0. ,
170 2 -.577350269189626,0.577350269189626,0. ,
171 2 0. ,0. ,0. ,
172 2 0. ,0. ,0. ,
173 3 -.774596669241483,0. ,0.774596669241483,
174 3 0. ,0. ,0. ,
175 3 0. ,0. ,0. ,
176 4 -.861136311594053,-.339981043584856,0.339981043584856,
177 4 0.861136311594053,0. ,0. ,
178 4 0. ,0. ,0. ,
179 5 -.906179845938664,-.538469310105683,0. ,
180 5 0.538469310105683,0.906179845938664,0. ,
181 5 0. ,0. ,0. ,
182 6 -.932469514203152,-.661209386466265,-.238619186083197,
183 6 0.238619186083197,0.661209386466265,0.932469514203152,
184 6 0. ,0. ,0. ,
185 7 -.949107912342759,-.741531185599394,-.405845151377397,
186 7 0. ,0.405845151377397,0.741531185599394,
187 7 0.949107912342759,0. ,0. ,
188 8 -.960289856497536,-.796666477413627,-.525532409916329,
189 8 -.183434642495650,0.183434642495650,0.525532409916329,
190 8 0.796666477413627,0.960289856497536,0. ,
191 9 -.968160239507626,-.836031107326636,-.613371432700590,
192 9 -.324253423403809,0. ,0.324253423403809,
193 9 0.613371432700590,0.836031107326636,0.968160239507626/
194
195 DATA w_lobatto /
196 1 2. ,0. ,0. ,
197 1 0. ,0. ,0. ,
198 1 0. ,0. ,0. ,
199 2 1. ,1. ,0. ,
200 2 0. ,0. ,0. ,
201 2 0. ,0. ,0. ,
202 3 0.333333333333333,1.333333333333333,0.333333333333333,
203 3 0. ,0. ,0. ,
204 3 0. ,0. ,0. ,
205 4 0.166666666666667,0.833333333333333,0.833333333333333,
206 4 0.166666666666667,0. ,0. ,
207 4 0. ,0. ,0. ,
208 5 0.1 ,0.544444444444444,0.711111111111111,
209 5 0.544444444444444,0.1 ,0. ,
210 5 0. ,0. ,0. ,
211 6 0.066666666666667,0.37847496 ,0.55485838 ,
212 6 0.55485838 ,0.37847496 ,0.066666666666667,
213 6 0. ,0. ,0. ,
214 7 0.04761904 ,0.27682604 ,0.43174538 ,
215 7 0.48761904 ,0.43174538 ,0.27682604 ,
216 7 0.04761904 ,0. ,0. ,
217 8 0.03571428 ,0.21070422 ,0.34112270 ,
218 8 0.41245880 ,0.41245880 ,0.34112270 ,
219 8 0.21070422 ,0.03571428 ,0. ,
220 9 0.027777777777778,0.1654953616 ,0.2745387126 ,
221 9 0.3464285110 ,0.3715192744 ,0.3464285110 ,
222 9 0.2745387126 ,0.1654953616 ,0.027777777777778/
223 DATA a_lobatto /
224 1 0. ,0. ,0. ,
225 1 0. ,0. ,0. ,
226 1 0. ,0. ,0. ,
227 2 -1. ,1. ,0. ,
228 2 0. ,0. ,0. ,
229 2 0. ,0. ,0. ,
230 3 -1. ,0. ,1. ,
231 3 0. ,0. ,0. ,
232 3 0. ,0. ,0. ,
233 4 -1. ,-.44721360 ,0.44721360 ,
234 4 1. ,0. ,0. ,
235 4 0. ,0. ,0. ,
236 5 -1. ,-.65465367 ,0. ,
237 5 0.65465367 , 1. ,0. ,
238 5 0. ,0. ,0. ,
239 6 -1. ,-.76505532 ,-.28523152 ,
240 6 0.28523152 ,0.76505532 , 1. ,
241 6 0. ,0. ,0. ,
242 7 -1. ,-.83022390 ,-.46884879 ,
243 7 0. ,0.46884879 ,0.83022390 ,
244 7 1. ,0. ,0. ,
245 8 -1. ,-.87174015 ,-.59170018 ,
246 8 -.20929922 ,0.20929922 ,0.59170018 ,
247 8 0.87174015 , 1. ,0. ,
248 9 -1. ,-.8997579954 ,-.6771862795 ,
249 9 -.3631174638 ,0. ,0.3631174638 ,
250 9 0.6771862795 ,0.8997579954 , 1. /
251
252
253 DATA w_newton /
254 1 2. ,0. ,0. ,
255 1 0. ,0. ,0. ,
256 1 0. ,0. ,0. ,
257 2 1. ,1. ,0. ,
258 2 0. ,0. ,0. ,
259 2 0. ,0. ,0. ,
260 3 0.5 ,1. ,0.5 ,
261 3 0. ,0. ,0. ,
262 3 0. ,0. ,0. ,
263 4 0.166666666666667,0.833333333333333,0.833333333333333,
264 4 0.166666666666667,0. ,0. ,
265 4 0. ,0. ,0. ,
266 5 0.25 ,0.5 ,0.5 ,
267 5 0.5 ,0.25 ,0. ,
268 5 0. ,0. ,0. ,
269 6 0.066666666666667,0.37847496 ,0.55485838 ,
270 6 0.55485838 ,0.37847496 ,0.066666666666667,
271 6 0. ,0. ,0. ,
272 7 0.04761904 ,0.27682604 ,0.43174538 ,
273 7 0.48761904 ,0.43174538 ,0.27682604 ,
274 7 0.04761904 ,0. ,0. ,
275 8 0.03571428 ,0.21070422 ,0.34112270 ,
276 8 0.41245880 ,0.41245880 ,0.34112270 ,
277 8 0.21070422 ,0.03571428 ,0. ,
278 9 0.027777777777778,0.1654953616 ,0.2745387126 ,
279 9 0.3464285110 ,0.3715192744 ,0.3464285110 ,
280 9 0.2745387126 ,0.1654953616 ,0.027777777777778/
281 DATA a_newton /
282 1 0. ,0. ,0. ,
283 1 0. ,0. ,0. ,
284 1 0. ,0. ,0. ,
285 2 -1. ,1. ,0. ,
286 2 0. ,0. ,0. ,
287 2 0. ,0. ,0. ,
288 3 -1. ,0. ,1. ,
289 3 0. ,0. ,0. ,
290 3 0. ,0. ,0. ,
291 4 -1. ,-.44721360 ,0.44721360 ,
292 4 1. ,0. ,0. ,
293 4 0. ,0. ,0. ,
294 5 -1. ,-.5 ,0. ,
295 5 0.5 , 1. ,0. ,
296 5 0. ,0. ,0. ,
297 6 -1. ,-.76505532 ,-.28523152 ,
298 6 0.28523152 ,0.76505532 , 1. ,
299 6 0. ,0. ,0. ,
300 7 -1. ,-.83022390 ,-.46884879 ,
301 7 0. ,0.46884879 ,0.83022390 ,
302 7 1. ,0. ,0. ,
303 8 -1. ,-.87174015 ,-.59170018 ,
304 8 -.20929922 ,0.20929922 ,0.59170018 ,
305 8 0.87174015 , 1. ,0. ,
306 9 -1. ,-.8997579954 ,-.6771862795 ,
307 9 -.3631174638 ,0. ,0.3631174638 ,
308 9 0.6771862795 ,0.8997579954 , 1. /
309
310
311
312 gbuf => elbuf_str%GBUF
313 nptr = elbuf_str%NPTR
314 npts = elbuf_str%NPTS
315 nptt = elbuf_str%NPTT
316 nlay = elbuf_str%NLAY
317
318 jhbe = iparg(23)
319 iint = iparg(36)
320 igtyp = iparg(38)
321 idef = 0
322 nf1=nft+1
323 nf2=nf1-(numels8+numels10+numels20)
324 ibid = 0
325 rbid = zero
326
327 DO i=lft,llt
328 rhocp(i) = pm(69,ixs(1,nft+i))
329 temp0(i) = pm(79,ixs(1,nft+i))
330 ENDDO
331
333 1 x ,v ,ixs(1,nf1) ,ixs16(1,nf2),xx ,
334 2 yy ,zz ,vx ,vy ,vz ,
335 3 nc ,ngl ,mat ,pid ,mass ,
336 4 dtelem(nft+1),sti ,gbuf%SIG ,gbuf%EINT ,gbuf%RHO,
337 5 gbuf%QVIS ,temp0 ,temp ,nel ,glob_therm%NINTEMP)
338
339 DO n=1,16
340 DO i=lft,llt
341 ul(i,n) = zero
342 ENDDO
343 ENDDO
344 DO i=lft,llt
345 volg(i) = zero
346 ENDDO
347
348
349
350
351 IF(jthe /=0)
CALL atheri(mat,pm,gbuf%TEMP)
352
353
354
355 is = 1
356 DO it=1,nptt
357 DO ir=1,nptr
358 DO il=1,nlay
359
360 lbuf => elbuf_str%BUFLY(il)%LBUF(ir,is,it)
361 mbuf => elbuf_str%BUFLY(il)%MAT(ir,is,it)
362 l_pla = elbuf_str%BUFLY(il)%L_PLA
363 l_sigb = elbuf_str%BUFLY(il)%L_SIGB
364 ip = ir + ( (il-1) + (it-1)*nlay )*nptr
365
366 IF (iint == 1) THEN
367
368 wi = w_gauss(ir,nptr)*w_gauss(il,nlay)*w_gauss(it,nptt)
369
371 1 a_gauss(ir,nptr),a_gauss(il,nlay),a_gauss(it,nptt),ni ,
372 2 dnidr ,dnids ,dnidt )
373
375 1 a_gauss(ir,nptr),a_gauss(il,nlay),a_gauss(it,nptt),wi,
376 2 dnidr ,dnids ,dnidt ,rx ,ry ,rz
377 3 sx ,sy ,sz ,tx ,ty ,tz ,
378 4 xx ,yy ,zz ,px ,py ,pz ,
379 5 lbuf%VOL,deltax ,stin ,ni ,volg ,ul ,lbuf%VOL0DP)
380 ELSEIF (iint == 2) THEN
381
382 wi = w_gauss(ir,nptr)*w_lobatto(il,nlay)*w_gauss(it,nptt)
383
385 1 a_gauss(ir,nptr),a_lobatto(il,nlay),a_gauss(it,nptt),ni ,
386 2 dnidr ,dnids ,dnidt )
387
389 1 a_gauss(ir,nptr),a_lobatto(il,nlay),a_gauss(it,nptt),wi,
390 2 dnidr ,dnids ,dnidt ,rx ,ry ,rz ,
391 3 sx ,sy ,sz ,tx ,ty ,tz ,
392 4 xx ,yy ,zz ,px ,py ,pz ,
393 5 lbuf%VOL,deltax ,stin ,ni ,volg ,ul ,lbuf%VOL0DP )
394 ENDIF
395
396 IF (jthe == 0 .and. glob_therm%NINTEMP > 0) THEN
397 CALL s20temp(nel,numnod,mvsiz,npe, nc,ni(1,ip),temp,tempel)
398 ELSE
399 tempel(1:nel) = temp0(1:nel)
400 ENDIF
401
402 CALL matini(pm ,ixs ,nixs ,x ,
403 . geo ,ale_connectivity ,detonators ,iparg ,
404 . sigi ,nel ,skew ,igeo(1,1) ,
405 . ipart ,iparts ,
406 . mat ,ipm ,nsigs ,numsol ,ptsol ,
407 . ip ,ngl ,npf ,tf ,bufmat ,
408 . gbuf ,lbuf ,mbuf ,elbuf_str ,iloadp ,
409 . facload, deltax ,tempel )
410
411
412
413 CALL s20msi(lbuf%RHO ,mass ,lbuf%VOL ,dtelem(nft+1),sti ,
414 . lbuf%OFF ,lbuf%SIG ,lbuf%EINT ,dtx ,nel ,
415 . gbuf%OFF ,gbuf%SIG ,gbuf%EINT ,gbuf%RHO ,wi/eight)
416
417
418
419
420 IF (mtn >= 28) THEN
421 nuvar = ipm(8,ixs(1,nft+1))
422 idef =1
423 ELSE
424 nuvar = 0
425 IF(mtn == 14 .OR. mtn == 12)THEN
426 idef =1
427 ELSEIF(mtn == 24)THEN
428 idef =1
429 ELSEIF(istrain == 1)THEN
430 IF(mtn == 1)THEN
431 idef =1
432 ELSEIF(mtn == 2)THEN
433 idef =1
434 ELSEIF(mtn == 4)THEN
435 idef =1
436 ELSEIF(mtn == 3.OR.mtn == 6.OR.mtn ==10.OR.
437 . mtn == 21.OR.mtn == 22.OR.
438 . mtn == 23.OR.mtn == 49)THEN
439 idef =1
440 ENDIF
441 ENDIF
442 ENDIF
443 CALL sigin20b(lbuf%SIG,pm ,lbuf%VOL,sigsp ,
444 . sigi ,lbuf%EINT,lbuf%RHO,mbuf%VAR,lbuf%STRA,
445 . ixs ,nixs ,nsigi ,ip ,nuvar ,
446 . nel ,iuser ,idef ,nsigs ,strsglob ,
447 . straglob,jhbe
448 . mat ,lbuf%PLA ,l_pla ,ptsol ,lbuf%SIGB,
449 . l_sigb ,ipm ,bufmat ,lbuf%VOL0DP)
450 ENDDO
451 ENDDO
452 ENDDO
453
454 DO i=lft,llt
455 aa =
max(ul(i,1),ul(i,2),ul(i,3),ul(i,4),
456 . ul(i,5),ul(i,6),ul(i,7),ul(i,8))
457 bb =
max(ul(i,9) ,ul(i,10),ul(i,11),ul(i,12),ul(i,13),ul(i,14),
458 . ul(i,15),ul(i,16))
459 deltax2(i) = aa/
max(aa,bb)
460 aa = aa*thirty2
461 bb = bb*thirty2*third
462 deltax(i) = sqrt(two*volg(i)/
max(aa,bb))
463 gbuf%VOL(i) = volg(i)
464 ENDDO
465
466 aire(:) = zero
467 CALL dtmain(geo ,pm ,ipm ,pid ,mat ,fv ,
468 . gbuf%EINT ,gbuf%TEMP ,gbuf%DELTAX ,gbuf%RK ,gbuf%RE ,bufmat, deltax, aire,
469 . gbuf%VOL, dtx, igeo,igtyp)
470
472 1 mass ,mas,partsav,iparts(nf1),mss(1,nf1),volg,
473 2 xx ,yy ,zz ,vx ,vy ,vz ,
474 3 nc ,sti,stifn ,deltax2 ,gbuf%RHO ,dtx ,
475 4 dtelem(nft+1),mssx(1,nf1),rhocp, mcp, mcps(1,nf1) ,
476 5 mcpsx(1,nf1) ,gbuf%FILL )
477
478
479
480 CALL failini(elbuf_str,nptr,npts,nptt,nlay,
481 . ipm,sigsp,nsigi,fail_ini ,
482 . sigi,nsigs,ixs,nixs,ptsol,
483 . rnoise,perturb,mat_param)
484
485
486
487
488 IF (i7stifs /= 0) THEN
489 ncc=16
490 CALL sbulk3(volg ,nc ,ncc,mat,pm ,
491 2 volnod,bvolnod,vns(1,nf1),bns(1,nf1),
492 3 vnsx(1,nf1),bnsx(1,nf1) ,gbuf%FILL)
493 ENDIF
494
495 DO i=lft,llt
496 IF(ixs(10,i+nft) /= 0) THEN
497 IF (igtyp/=0 .AND. igtyp /= 14 .AND. igtyp/=20 .AND.
498 . igtyp/=21) THEN
499 ipid1=ixs(nixs-1,i+nft)
500 CALL fretitl2(titr1,igeo(npropgi-ltitr+1,ipid1),ltitr)
502 . msgtype=msgerror,
503 . anmode=aninfo_blind_1,
504 . i1=igeo(1,ipid1),
505 . c1=titr1,
506 . i2=igtyp)
507 ENDIF
508 ENDIF
509 ENDDO
510
511 RETURN
subroutine atheri(mat, pm, temp)
subroutine dtmain(geo, pm, ipm, pid, mat, fv, eint, temp, deltax, rk, re, bufmat, ddeltax, aire, vol, dtx, igeo, igtyp)
subroutine failini(elbuf_str, nptr, npts, nptt, nlay, ipm, sigsp, nsigi, fail_ini, sigi, nsigs, ix, nix, pt, rnoise, perturb, mat_param)
subroutine matini(pm, ix, nix, x, geo, ale_connectivity, detonators, iparg, sigi, nel, skew, igeo, ipart, ipartel, mat, ipm, nsig, nums, pt, ipt, ngl, npf, tf, bufmat, gbuf, lbuf, mbuf, elbuf_str, iloadp, facload, ddeltax, tempel)
integer, parameter nchartitle
subroutine s16coor3(x, v, ixs, ixs16, xx, yy, zz, vx, vy, vz, nc, ngl, mxt, ngeo, mass, dtelem, sti, sigg, eintg, rhog, qg, temp0, temp, nel, nintemp)
subroutine s16mass3(mass, ms, partsav, ipart, mss, volg, xx, yy, zz, vx, vy, vz, nc, sti, stifn, deltax2, rho, dtx, dtelem, mssx, rhocp, mcp, mcps, mcpsx, fill)
subroutine s20msi(rho, mass, volu, dtelem, sti, off, sig, eint, dtx, nel, offg, sigg, eintg, rhog, wip)
subroutine sigin20b(sig, pm, vol, sigsp, sigi, eint, rho, uvar, eps, ix, nix, nsigi, ipt, nuvar, nel, iuser, idef, nsigs, strsglob, straglob, jhbe, igtyp, x, bufgama, mat, epsp, l_pla, pt, sigb, l_sigb, ipm, bufmat, voldp)
subroutine sbulk3(volu, nc, nnc, mat, pm, volnod, bvolnod, vns, bns, vnsx, bnsx, fill)
subroutine s16rst(r, s, t, ni, dnidr, dnids, dnidt)
subroutine s16deri3(ngl, off, r, s, t, w, dnidr, dnids, dnidt, dxdr, dydr, dzdr, dxds, dyds, dzds, dxdt, dydt, dzdt, xx, yy, zz, px, py, pz, vol, deltax, kxx, ni, volg, ul, voldp)
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)