48 use element_mod , only :nixs,nixc,nixtg
49
50
51
52#include "implicit_f.inc"
53
54
55
56#include "mvsiz_p.inc"
57
58
59
60#include "com04_c.inc"
61#include "param_c.inc"
62#include "units_c.inc"
63
64
65
66 INTEGER NRTS, NRTM, NINT, NTY, NOINT,NSN, NMN, IGAP,
67 . INACTI, INTTH
68 INTEGER IRECTS(4,*), IRECTM(4,*), IXS(NIXS,*), IXC(NIXC,*),
69 . NSV(*), IXTG(NIXTG,*),
70 . KNOD2ELS(*), KNOD2ELC(*), KNOD2ELTG(*), NOD2ELS(*), NOD2ELC(*),
71 . NOD2ELTG(*),IELES(*),
72 . MSR(*), ITAB(*), IKINE(*), IPARTS(*), IPARTC(*), IPARTG(*)
73 INTEGER , INTENT (IN) :: RESORT
74
75 my_real ,
INTENT(IN) :: dgapload
77 . gap,gapmin,criter, gapmax, gapscale, depth, drad, lxm, lym, lzm
79 . x(3,*), pm(npropm,*), geo(npropg,*),
80 . gap_s(*), thknod(*), stf(*), stfn(*),
81 . gap_s0(*), area_s0(*), xm0(3,*),thk_part(*),thknod0(*)
82 INTEGER ID
83 CHARACTER(LEN=NCHARTITLE) :: TITR
84
85
86
87 INTEGER NDX, I, J, II, INRT, NEL,
88 . N1,N2,N3,N4, IX, N, L, LLT, NN, IP, STAT
89 INTEGER ITMP(NUMNOD)
90
92 . dxm, gapmx, gapmn, dx,gaps1,gaps2, gapm,
93 . xxx, yyy, zzz, x0, x1, y0, y1, z0, z1
95 . x12(mvsiz),y12(mvsiz),z12(mvsiz),
96 . x13(mvsiz),y13(mvsiz),z13(mvsiz),
97 . x24(mvsiz),y24(mvsiz),z24(mvsiz),
98 . nx(mvsiz),ny(mvsiz),nz(mvsiz),aa(mvsiz)
99 my_real,
DIMENSION(:),
ALLOCATABLE :: thk_part_nods
100
101 dxm=zero
102 ndx = 0
103 gapmx=ep30
104 gapmn=ep30
105 gaps1=zero
106 gaps2=zero
107
108
109
110 DO 250 i=1,nrts
111 gapm =zero
112 inrt=i
113 CALL i4gmx3(x,irects,inrt,gapmx)
114 250 CONTINUE
115
116
117
118 IF(igap>=1)THEN
119 ALLOCATE (thk_part_nods(numnod) ,stat=stat)
120 IF (stat /= 0)
CALL ancmsg(msgid=268,anmode=aninfo,
121 . msgtype=msgerror,
122 . c1='THK_PART_NODS')
123 thk_part_nods(1:numnod) = zero
124 DO i=1,nrts
125 nel = ieles(i)
126 IF(nel<=numels) THEN
127 ip = iparts(nel)
128 DO n =1,4
129 nn = irects(n,i)
130 thk_part_nods(nn) =
max(thk_part_nods(n),thk_part(ip))
131 ENDDO
132 ELSEIF(nel<=(numels+numelc)) THEN
133 ip = ipartc(nel-numels)
134 DO n =1,4
135 nn = irects(n,i)
136 thk_part_nods(nn) =
max(thk_part_nods(n),thk_part(ip))
137 ENDDO
138 ELSE
139 ip = ipartg(nel-numels-numelc)
140 DO n =1,4
141 nn = irects(n,i)
142 thk_part_nods(nn) =
max(thk_part_nods(n),thk_part(ip))
143 ENDDO
144 ENDIF
145 ENDDO
146 ENDIF
147
148
149
150
151 IF(igap>=1)THEN
152 DO i=1,nsn
153 IF(thk_part_nods(nsv(i))/=zero) THEN
154 dx = thk_part_nods(nsv(i))*gapscale
155 ELSE
156 dx = thknod(nsv(i))*gapscale
157 ENDIF
158 gapm = half*dx
159
160 gaps2 =
max(gaps2,gapm)
161 gap_s(i)= gapm
162
163
164
165 dxm=dxm+dx
166 ndx=ndx+1
167
168 thknod0(i) = thknod(nsv(i))
169 ENDDO
170 IF (ALLOCATED(thk_part_nods)) DEALLOCATE(thk_part_nods)
171 ENDIF
172
173
174
175
176 gapmx=sqrt(gapmx)
177 IF(igap==0)THEN
178
179 IF(gap<=zero)THEN
180 DO i=1,nsn
181 dx = thknod(nsv(i))
182
183
184
185 dxm=dxm+dx
186 ndx=ndx+1
187 ENDDO
188 gap = half*dxm/ndx
189 IF (resort==0) WRITE(iout,1000)gap
190 ENDIF
191 gapmin = gap
192 gapmax = gap
193 ELSE
194
195 IF(gap>zero)gapmin=gap
196 IF (resort==0) WRITE(iout,1000)gapmin
197
198
199 IF(gapmax==zero)gapmax=ep30
200 IF (resort==0) WRITE(iout,1500)gapmax
201 gap =
min(gap,gapmax)
202 ENDIF
203
204
205
206 gap =
min(gapmax,
max(gaps2,gapmin))
207
208
209
210
211 IF (igap==0) THEN
212 criter=gap
213 ELSE
214 criter=ep30
215 DO i = 1, nsn
216 criter =
min(criter,gap_s(i))
217 ENDDO
218 criter=
max(criter,gapmin)
219 ENDIF
220
221 IF(dgapload > zero) criter=
max(criter,em01*(gap + dgapload))
222
223 IF(depth==zero)THEN
224
226
227 ELSEIF(depth<gap)THEN
228
229 depth=gap
230 END IF
231 IF (resort==0) WRITE(iout,2000)depth
232
233 criter=
max(criter,em01*depth)
234
235 IF(depth>gapmx .AND. resort==0 )THEN
237 . msgtype=msgwarning,
238 . anmode=aninfo_blind_2,
240 . c1=titr,
241 . r1=depth,
242 . r2=gapmx,
244 ENDIF
245
246 IF(intth/=0)THEN
247 IF(drad==zero)THEN
249 ELSEIF(drad<gap)THEN
250 drad=gap
251 END IF
252 IF (resort==0) WRITE(iout,2001)drad
253
254 criter=
max(criter,em01*drad)
255
256 IF(drad>gapmx .AND. resort==0)THEN
258 . msgtype=msgwarning,
259 . anmode=aninfo_blind_2,
261 . c1=titr,
262 . r1=drad ,
263 . r2=gapmx,
265 END IF
266 END IF
267
268
269
270 DO i=1,nrtm
271 stf(i)=one
272 END DO
273
274
275
276 DO i=1,nsn
277 stfn(i) = one
278 END DO
279
280 IF(igap==2)THEN
281 DO i=1,nsn
282 gap_s0(i) =
min(gap_s(i),gapmax)
283 gap_s0(i) =
max(gapmin ,gap_s0(i))
284 END DO
285
286 IF(intth == 0) THEN
287 itmp=0
288 DO i=1,nsn
289 ii=nsv(i)
290 itmp(ii)=i
291 END DO
292 DO n=1,nrts,mvsiz
293
294 llt=
min(nrts-n+1,mvsiz)
295
296 DO l=1,llt
297 i=n+l-1
298
299 n1=irects(1,i)
300 n2=irects(2,i)
301 n3=irects(3,i)
302 n4=irects(4,i)
303 IF(n4/=n3)THEN
304 x13(l)=x(1,n3)-x(1,n1)
305 y13(l)=x(2,n3)-x(2,n1)
306 z13(l)=x(3,n3)-x(3,n1)
307 x24(l)=x(1,n4)-x(1,n2)
308 y24(l)=x(2,n4)-x(2,n2)
309 z24(l)=x(3,n4)-x(3,n2)
310 nx(l)=y13(l)*z24(l)-z13(l)*y24(l)
311 ny(l)=z13(l)*x24(l)-x13(l)*z24(l)
312 nz(l)=x13(l)*y24(l)-y13(l)*x24(l)
313 aa(l)=one_over_8*sqrt(nx(l)*nx(l)+ny(l)*ny(l)+nz(l)*nz(l))
314 area_s0(itmp(n1))=area_s0(itmp(n1))+aa(l)
315 area_s0(itmp(n2))=area_s0(itmp(n2))+aa(l)
316 area_s0(itmp(n3))=area_s0(itmp(n3))+aa(l)
317 area_s0(itmp(n4))=area_s0(itmp(n4))+aa(l)
318 ELSE
319 x12(l)=x(1,n2)-x(1,n1)
320 y12(l)=x(2,n2)-x(2,n1)
321 z12(l)=x(3,n2)-x(3,n1)
322 x13(l)=x(1,n3)-x(1,n1)
323 y13(l)=x(2,n3)-x(2,n1)
324 z13(l)=x(3,n3)-x(3,n1)
325 nx(l)=y12(l)*z13(l)-z12(l)*y13(l)
326 ny(l)=z12(l)*x13(l)-x12(l)*z13(l)
327 nz(l)=x12(l)*y13(l)-y12(l)*x13(l)
328 aa(l)=one_over_6*sqrt(nx(l)*nx(l)+ny(l)*ny(l)+nz(l)*nz(l))
329 area_s0(itmp(n1))=area_s0(itmp(n1))+aa(l)
330 area_s0(itmp(n2))=area_s0(itmp(n2))+aa(l)
331 area_s0(itmp(n3))=area_s0(itmp(n3))+aa(l)
332 END IF
333 END DO
334 END DO
335 igap = 1
336 ENDIF
337 ELSE
338 IF(intth==0) THEN
339 itmp=0
340 DO i=1,nsn
341 ii=nsv(i)
342 itmp(ii)=i
343 END DO
344 DO n=1,nrts,mvsiz
345
346 llt=
min(nrts-n+1,mvsiz)
347
348 DO l=1,llt
349 i=n+l-1
350
351 n1=irects(1,i)
352 n2=irects(2,i)
353 n3=irects(3,i)
354 n4=irects(4,i)
355 IF(n4/=n3)THEN
356 x13(l)=x(1,n3)-x(1,n1)
357 y13(l)=x(2,n3)-x(2,n1)
358 z13(l)=x(3,n3)-x(3,n1)
359 x24(l)=x(1,n4)-x(1,n2)
360 y24(l)=x(2,n4)-x(2,n2)
361 z24(l)=x(3,n4)-x(3,n2)
362 nx(l)=y13(l)*z24(l)-z13(l)*y24(l)
363 ny(l)=z13(l)*x24(l)-x13(l)*z24(l)
364 nz(l)=x13(l)*y24(l)-y13(l)*x24(l)
365 aa(l)=one_over_8*sqrt(nx(l)*nx(l)+ny(l)*ny(l)+nz(l)*nz(l))
366 area_s0(itmp(n1))=area_s0(itmp(n1))+aa(l)
367 area_s0(itmp(n2))=area_s0(itmp(n2))+aa(l)
368 area_s0(itmp(n3))=area_s0(itmp(n3))+aa(l)
369 area_s0(itmp(n4))=area_s0(itmp(n4))+aa(l)
370 ELSE
371 x12(l)=x(1,n2)-x(1,n1)
372 y12(l)=x(2,n2)-x(2,n1)
373 z12(l)=x(3,n2)-x(3,n1)
374 x13(l)=x(1,n3)-x(1,n1)
375 y13(l)=x(2,n3)-x(2,n1)
376 z13(l)=x(3,n3)-x(3,n1)
377 nx(l)=y12(l)*z13(l)-z12(l)*y13(l)
378 ny(l)=z12(l)*x13(l)-x12(l)*z13(l)
379 nz(l)=x12(l)*y13(l)-y12(l)*x13(l)
380 aa(l)=one_over_6*sqrt(nx(l)*nx(l)+ny(l)*ny(l)+nz(l)*nz(l))
381 area_s0(itmp(n1))=area_s0(itmp(n1))+aa(l)
382 area_s0(itmp(n2))=area_s0(itmp(n2))+aa(l)
383 area_s0(itmp(n3))=area_s0(itmp(n3))+aa(l)
384 END IF
385 END DO
386 END DO
387 ENDIF
388 ENDIF
389
390 lxm=zero
391 lym=zero
392 lzm=zero
393 DO i=1,nrtm
394 x0=ep30
395 x1=-ep30
396 y0=ep30
397 y1=-ep30
398 z0=ep30
399 z1=-ep30
400 DO j=1,4
401 ix=msr(irectm(j,i))
402 xxx=x(1,ix)
403 yyy=x(2,ix)
404 zzz=x(3,ix)
411 END DO
415 ENDDO
416
417 RETURN
418 1000 FORMAT(2x,'GAP MIN = ',1pg20.13)
419 1500 FORMAT(2x,'GAP MAX = ',1pg20.13)
420 2000 FORMAT(2x,'DEPTH BEFORE RELEASE = ',1pg20.13)
421 2001 FORMAT(2x,'Maximum distance for radiation computation = ',
422 . 1pg20.13)
subroutine i4gmx3(x, irect, i, gapmax)
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)