60
61
62
63 USE my_alloc_mod
65 USE intbuf_fric_mod
68
69
70
71#include "implicit_f.inc"
72
73
74
75#include "com01_c.inc"
76#include "com04_c.inc"
77#include "param_c.inc"
78#include "remesh_c.inc"
79#include "scr08_c.inc"
80#include "scr17_c.inc"
81#include "units_c.inc"
82
83
84
85 INTEGER NRT, NINT, NTY, NOINT,NSN,IGAP,
86 . INACTI,INTFRIC, NRT_IGE
87 INTEGER (4,*), IXS(NIXS,*), IXC(NIXC,*),
88 . NSV(*), (NIXTG,*), IXT(NIXT,*), IXP(NIXP,*),
89 . KNOD2ELS(*), KNOD2ELC(*), KNOD2ELTG(*), NOD2ELS(*), NOD2ELC(*),
90 . NOD2ELTG(*), IELES(*), , IELEC(*),
91 . SH3TREE(KSH3TREE,*), SH4TREE(KSH4TREE,*),IXR(NIXR,*) ,
92 . IPART(LIPART1,*),IPARTC(*),IPARTTG(*),NOD2EL1D(*),KNOD2EL1D(*),
93 . ITAB(*), IXS10(6,*), IXS16(*), IXS20(*),IDDLEVEL,IGEO(NPROPGI,*),
94 . IWORKSH(3,*),IT19,KXIG3D(NIXIG3D,*),IXIG3D(*),TAGPRT_FRIC(*),
95 . IPARTFRICS(*),IPARTFRICM(*),IPARTS(*),IREM_GAP
97 . stfac, gap,gapmin,gapinf, gapmax,percent_size, bgapsmx,
98 . gapinf_s, gapinf_m, drad, gapm_mx, gaps_mx, gaps_l_mx, gapm_l_mx
100 . x(3,*), stf(*), pm(npropm,*), geo(npropg,*), stfn(*),
101 . ms(*),wa(*),gap_s(*),gap_m(*),
102 . areas(*),thk(*),thk_part(*),
103 . gap_s_l(*),gap_m_l(*), fillsol(*),pm_stack(20,*)
104 INTEGER, DIMENSION(NUMELT), INTENT(IN) :: IPARTT
105 INTEGER, DIMENSION(NUMELP), INTENT(IN) :: IPARTP
106 INTEGER, DIMENSION(NUMELR), INTENT(IN) :: IPARTR
107 INTEGER, INTENT(IN) :: ID
108 CHARACTER(LEN=NCHARTITLE)::TITR
109 TYPE(INTBUF_FRIC_STRUCT_) INTBUF_FRIC_TAB(*)
110 TYPE (SURF_) :: IGRSURF
111 INTEGER, DIMENSION(NUMELS), INTENT(INOUT):: ELEM_LINKED_TO_SEGMENT
112 INTEGER, INTENT(IN) :: FLAG_ELEM_INTER25(NINTER25,NUMELS)
113
114
115
116 INTEGER NDX, I, J, INRT, NELS, MT, JJ, JJJ, NELC,
117 . MG, NUM, NPT, LL, L, NN, NELTG,N1,N2,N3,N4,IE,
118 . IP, NLEV, MYLEV, K, P, R, T,IGTYP,IPGMAT,IGMAT,
119 . ISUBSTACK,NELIG3D, COIN_IGE(8), IPID, PX, PY, PZ, IAD ,IPFMAX,IPL,
120 . IPFLMAX,IPG,NINV,ICONTR,NIN25
121 INTEGER JPERM(4)
122 LOGICAL TYPE18
124 . dxm, gapmx, gapmn,
area, vol, dx, gapm, ddx,
125 . gaptmp, gapscale,sx1,sy1,sz1,sx2,sy2,sz2,sx3,sy3,sz3,
126 . slsfac,xl,bulk
127 my_real,
dimension(:),
allocatable :: gap_s_l_tmp
128 INTEGER, DIMENSION(:),ALLOCATABLE :: TAGELEMS,INDEXE
129 LOGICAL :: PRINT_ERROR
130 INTEGER, DIMENSION(4) :: NODE_ID
131 INTEGER :: IELEM(2)
132 DATA jperm/2,3,4,1/
133
134
135
136
137
138
139
140 ALLOCATE( gap_s_l_tmp(numnod) )
141 type18 = .false.
142 IF(inacti == 7)type18=.true.
143 ipgmat = 700
144 IF(nty == 20)THEN
145 slsfac = one
146 ELSE
147 slsfac = stfac
148 ENDIF
149 IF(igap == 3)THEN
150 DO i=1,numnod
151 gap_s_l_tmp(i)=zero
152 ENDDO
153 DO i=1,nrt+nrt_ige
154 gap_m_l(i)=ep30
155 ENDDO
156 DO i=1,nsn
157 gap_s_l(i)=ep30
158 ENDDO
159 ENDIF
160 dxm=0.
161 ndx=0
162 gapmx=ep30
163 gapmn=ep30
164 gaps_mx=zero
165 gaps_l_mx=zero
166 gapm_mx=zero
167 gapm_l_mx=zero
168
169 IF(igap == 2)THEN
170 IF(iddlevel == 1) igap = 1
171 gapscale = gapmin
172 gapmin = zero
173 ELSEIF(igap == 3)THEN
174 gapscale=gapmin
175 gapmin=zero
176 ELSE
177 gapscale = one
178 ENDIF
179
180 IF(igap == 3)THEN
181 DO i=1,nrt+nrt_ige
182 xl = ep30
183 DO j=1,4
184 n1=irect(j,i)
185 n2=irect(jperm(j),i)
186 IF(n1 /= n2 .AND. n1 /= 0)xl=
min(xl,sqrt((x(1,n1)-x(1,n2))**2+(x(2,n1)-x(2,n2))**2+(x(3,n1)-x(3,n2))**2))
187 ENDDO
188 DO j=1,4
189 IF(gap_s_l_tmp(irect(j,i)) == zero) THEN
190 gap_s_l_tmp(irect(j,i))= percent_size*xl
191 ELSE
192 gap_s_l_tmp(irect(j,i))=
min(gap_s_l_tmp(irect(j,i)),percent_size*xl)
193 ENDIF
194 ENDDO
195
196 DO j=1,4
197 n1=irect(j,i)
198 DO k=knod2el1d(n1)+1,knod2el1d(n1+1)
199 IF (nod2el1d(k) <= numelt .AND. nod2el1d(k) /= zero) THEN
200 t=nod2el1d(k)
201 xl=
min(xl,sqrt((x(1,ixt(2,t))-x(1,ixt(3,t)))**2 + (x(2,ixt(2,t))-x(2,ixt(3,t)))**2 + (x(3,ixt(2,t))-x(3,ixt(3,t)))**2))
202 ELSEIF (nod2el1d(k) <= numelt+numelp .AND. nod2el1d(k) /= zero) THEN
203 p=nod2el1d(k) - numelt
204 xl=
min(xl,sqrt((x(1,ixp(2,p))-x(1,ixp(3,p)))**2 + (x(2,ixp(2,p))-x(2,ixp(3,p)))**2 + (x(3,ixp(2,p))-x(3,ixp(3,p)))**2))
205 ELSEIF (nod2el1d(k) <= numelt+numelp+numelr .AND. nod2el1d(k) /= zero) THEN
206 r=nod2el1d(k) - numelt - numelp
207 xl=
min(xl,sqrt((x(1,ixr(2,r))-x(1,ixr(3,r)))**2 + (x(2,ixr(2,r))-x(2,ixr(3,r)))**2 + (x(3,ixr(2,r))-x(3,ixr(3,r)))**2))
208 ENDIF
209 ENDDO
210 ENDDO
211 gap_m_l(i)=percent_size*xl
212 gapm_l_mx =
max(gapm_l_mx,gap_m_l(i))
213 DO j=1,4
214 gap_s_l_tmp(irect(j,i))=
min(gap_s_l_tmp(irect(j,i)),percent_size*xl)
215 ENDDO
216 ENDDO
217 ENDIF
218
219
220
221 IF(igap >= 1)THEN
222 DO i=1,numnod
223 wa(i)=zero
224 ENDDO
225 DO i=1,numelc
226 mg=ixc(6,i)
227 igtyp = igeo(11,mg)
228 ip = ipartc(i)
229 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
230 dx=half*thk_part(ip)
231 ELSEIF ( thk(i) /= zero .AND. iintthick == 0) THEN
232 dx=half*thk(i)
233 ELSEIF(igtyp == 17 .OR. igtyp ==51 .OR. igtyp == 52) THEN
234 dx=half*thk(i)
235 ELSE
236 dx=half*geo(1,mg)
237 ENDIF
238 wa(ixc(2,i))=
max(wa(ixc(2,i)),dx)
239 wa(ixc(3,i))=
max(wa(ixc(3,i)),dx)
240 wa(ixc(4,i))=
max(wa(ixc(4,i)),dx)
241 wa(ixc(5,i))=
max(wa(ixc(5,i)),dx)
242 ENDDO
243 DO i=1,numeltg
244 mg=ixtg(5,i)
245 igtyp = igeo(11,mg)
246 ip = iparttg(i)
247 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
248 dx=half*thk_part(ip)
249 ELSEIF ( thk(numelc+i) /= zero .AND. iintthick == 0) THEN
250 dx=half*thk(numelc+i)
251 ELSEIF(igtyp == 17 .OR. igtyp ==51 .OR. igtyp == 52) THEN
252 dx=half*thk(numelc+i)
253 ELSE
254 dx=half*geo(1,mg)
255 ENDIF
256 wa(ixtg(2,i))=
max(wa(ixtg(2,i)),dx)
257 wa(ixtg(3,i))=
max(wa(ixtg(3,i)),dx)
258 wa(ixtg(4,i))=
max(wa(ixtg(4,i)),dx)
259 ENDDO
260 DO i=1,numelt
261 mg=ixt(4,i)
262 ip = ipartt(i)
263 IF ( thk_part(ip) > zero ) THEN
264 dx=half*thk_part(ip)
265 ELSE
266 dx=half*sqrt(geo(1,mg))
267 END IF
268 wa(ixt(2,i))=
max(wa(ixt(2,i)),dx)
269 wa(ixt(3,i))=
max(wa(ixt(3,i)),dx)
270 ENDDO
271 DO i=1,numelp
272 mg=ixp(5,i)
273 ip = ipartp(i)
274 IF ( thk_part(ip) > zero ) THEN
275 dx=half*thk_part(ip)
276 ELSE
277 dx=half*sqrt(geo(1,mg))
278 END IF
279 wa(ixp(2,i))=
max(wa(ixp(2,i)),dx)
280 wa(ixp(3,i))=
max(wa(ixp(3,i)),dx)
281 ENDDO
282 DO i=1,numelr
283 ip = ipartr(i)
284 IF ( thk_part(ip) > zero ) THEN
285 mg=ixr(1,i)
286 igtyp = igeo(11,mg)
287 dx=half*thk_part(ip)
288 wa(ixr(2,i))=
max(wa(ixr(2,i)),dx)
289 wa(ixr(3,i))=
max(wa(ixr(3,i)),dx)
290 IF (igtyp==12) wa(ixr(4,i))=
max(wa(ixr(4,i)),dx)
291 END IF
292 ENDDO
293 IF(type18)THEN
294 gaps_mx = zero
295 gap_s(1:nsn) = zero
296 ELSE
297 DO i=1,nsn
298 gap_s(i)=gapscale * wa(nsv(i))
299 IF(igap == 3)THEN
300 IF(gap_s_l_tmp(nsv(i)) /= zero)gap_s_l(i)=
min(gap_s_l(i),gap_s_l_tmp(nsv(i)))
301 gaps_mx =
max(gaps_mx,gap_s(i))
302 gaps_l_mx =
max(gaps_l_mx,gap_s_l(i))
303 ELSE
304 gaps_mx=
max(gaps_mx,gap_s(i))
305 END IF
306 ENDDO
307 ENDIF
308 ENDIF
309
310
311 IF(intth > 0 ) THEN
312 IF(nadmesh == 0)THEN
313 DO i = 1,nsn
314 areas(i) = zero
315 DO j= knod2elc(nsv(i))+1,knod2elc(nsv(i)+1)
316 ie = nod2elc(j)
317 sx1 = x(1,ixc(4,ie)) - x(1,ixc(2,ie))
318 sy1 = x(2,ixc(4,ie)) - x(2,ixc(2,ie))
319 sz1 = x(3,ixc(4,ie)) - x(3,ixc(2,ie))
320 sx2 = x(1,ixc(5,ie)) - x(1,ixc(3,ie))
321 sy2 = x(2,ixc(5,ie)) - x(2,ixc(3,ie))
322 sz2 = x(3,ixc(5,ie)) - x(3,ixc(3,ie))
323 sx3 = sy1*sz2 - sz1*sy2
324 sy3 = sz1*sx2 - sx1*sz2
325 sz3 = sx1*sy2 - sy1*sx2
326 areas(i) = areas(i) + one_over_8*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
327
328 ielec(i) = ixc(1,ie)
329 END DO
330
331 DO j= knod2eltg(nsv(i))+1,knod2eltg(nsv(i)+1)
332 ie = nod2eltg(j)
333 sx1 = x(1,ixtg(3,ie)) - x(1,ixtg(2,ie))
334 sy1 = x(2,ixtg(3,ie)) - x(2,ixtg(2,ie))
335 sz1 = x(3,ixtg(3,ie)) - x(3,ixtg(2,ie))
336 sx2 = x(1,ixtg(4,ie)) - x(1,ixtg(2,ie))
337 sy2 = x(2,ixtg(4,ie)) - x(2,ixtg(2,ie))
338 sz2 = x(3,ixtg(4,ie)) - x(3,ixtg(2,ie))
339 sx3 = sy1*sz2 - sz1*sy2
340 sy3 = sz1*sx2 - sx1*sz2
341 sz3 = sx1*sy2 - sy1*sx2
342 areas(i) = areas(i) + one_over_6*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
343
344 ielec(i) = ixtg(1,ie)
345 END DO
346 END DO
347 ELSE
348 DO i = 1,nsn
349 areas(i) = zero
350 DO j=knod2elc(nsv(i))+1,knod2elc(nsv(i)+1)
351 ie = nod2elc(j)
352
353 ip = ipartc(ie)
354 nlev =ipart(10,ip)
355 mylev=sh4tree(3,ie)
356 IF(mylev < 0) mylev=-(mylev+1)
357
358 IF(mylev == nlev)THEN
359 sx1 = x(1,ixc(4,ie)) - x(1,ixc(2,ie))
360 sy1 = x(2,ixc(4,ie)) - x(2,ixc(2,ie))
361 sz1 = x(3,ixc(4,ie)) - x(3,ixc(2,ie))
362 sx2 = x(1,ixc(5,ie)) - x(1,ixc(3,ie))
363 sy2 = x(2,ixc(5,ie)) - x(2,ixc(3,ie))
364 sz2 = x(3,ixc(5,ie)) - x(3,ixc(3,ie))
365 sx3 = sy1*sz2 - sz1*sy2
366 sy3 = sz1*sx2 - sx1*sz2
367 sz3 = sx1*sy2 - sy1*sx2
368 areas(i) = areas(i) + one_over_8*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
369
370 ielec(i) = ixc(1,ie)
371 END IF
372
373 END DO
374
375 DO j= knod2eltg(nsv(i))+1,knod2eltg(nsv(i)+1)
376 ie = nod2eltg(j)
377 ip = iparttg(ie)
378 nlev =ipart(10,ip)
379 mylev=sh3tree(3,ie)
380 IF(mylev < 0) mylev=-(mylev+1)
381 IF(mylev == nlev)THEN
382 sx1 = x(1,ixtg(3,ie)) - x(1,ixtg(2,ie))
383 sy1 = x(2,ixtg(3,ie)) - x(2,ixtg(2,ie))
384 sz1 = x(3,ixtg(3,ie)) - x(3,ixtg(2,ie))
385 sx2 = x(1,ixtg(4,ie)) - x(1,ixtg(2,ie))
386 sy2 = x(2,ixtg(4,ie)) - x(2,ixtg(2,ie))
387 sz2 = x(3,ixtg(4,ie)) - x(3,ixtg(2,ie))
388 sx3 = sy1*sz2 - sz1*sy2
389 sy3 = sz1*sx2 - sx1*sz2
390 sz3 = sx1*sy2 - sy1*sx2
391 areas(i) = areas(i) + one_over_6*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
392
393 ielec(i) = ixtg(1,ie)
394 END IF
395
396 END DO
397 END DO
398 END IF
399 END IF
400
401
402
403
404 IF(numels > 0) THEN
405 CALL my_alloc(tagelems,numels)
406 tagelems = 0
407 CALL my_alloc(indexe,numels)
408 indexe = 0
409 ENDIF
410 ninv = 0
411 DO i=1,nrt
412 stf(i)=zero
413 IF(intth > 0 ) ieles(i) = 0
414 IF(slsfac < zero)stf(i)=slsfac
415 gapm=zero
416 inrt=i
417 CALL i4gmx3(x,irect,inrt,gapmx)
418
419 CALL inelts(x ,irect,ixs ,nint,nels ,
420 . inrt ,
area ,noint,0 ,igrsurf%ELTYP,
421 . igrsurf%ELEM)
422
423
424
425
426 IF(nels /= 0)THEN
427 mt=ixs(1,nels)
428 mg=ixs(nixs-1,nels)
429 icontr = igeo(97,mg)
430 IF(mt > 0)THEN
431 DO jj=1,8
432 jjj=ixs(jj+1,nels)
433 xc(jj)=x(1,jjj)
434 yc(jj)=x(2,jjj)
435 zc(jj)=x(3,jjj)
436 END DO
438 IF (icontr==1 ) THEN
439 bulk = pm(107,mt)
440 ELSE
441 bulk = pm(32,mt)
442 END IF
443 stf(i)=slsfac*fillsol(nels)*
area*
area*bulk/vol
444 ELSE
445 IF(nint >= 0) THEN
447 . msgtype=msgwarning,
448 . anmode=aninfo_blind_2,
450 . c1=titr,
451 . i2=ixs(nixs,nels),
452 . c2='SOLID',
453 . i3=i)
454 ENDIF
455 IF(nint < 0) THEN
457 . msgtype=msgwarning,
458 . anmode=aninfo_blind_2,
460 . c1=titr,
461 . i2=ixs(nixs,nels),
462 . c2='SOLID',
463 . i3=i)
464 ENDIF
465 ENDIF
466 IF(igap /= 0 .OR. (nty /=7 .AND. nty /= 20)) gap_m(i)=gapm
467
468 IF(intfric > 0) THEN
469 ip= iparts(nels)
470 ipg = tagprt_fric(ip)
471 IF(ipg > 0) THEN
473 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
474 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
475 ipartfricm(inrt) = ipl
476 ENDIF
477 ENDIF
478
479 cycle
480 ENDIF
481
482
483 CALL ineltc(nelc ,neltg ,inrt ,igrsurf%ELTYP, igrsurf%ELEM)
484
485
486
487 IF(neltg /= 0) THEN
488 mt=ixtg(1,neltg)
489 mg=ixtg(5,neltg)
490 igtyp = igeo(11,mg)
491 igmat = igeo(98,mg)
492 ip = iparttg(neltg)
493 IF (thk_part(ip) /= zero .AND. iintthick == 0) THEN
494 dx=thk_part(ip)*gapscale
495 ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0)THEN
496 dx=thk(numelc+neltg)*gapscale
497 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
498 dx=thk(numelc+neltg)*gapscale
499 ELSE
500 dx=geo(1,mg)*gapscale
501 ENDIF
502 gapm=half*dx
503 gapm_mx=
max(gapm_mx,gapm)
505 dxm=dxm+dx
506 ndx=ndx+1
507 IF(mt > 0)THEN
508 IF(igtyp == 11 .AND. igmat > 0) THEN
509 IF(slsfac < zero)THEN
510 stf(i)=-slsfac
511 ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0)THEN
512 stf(i)=slsfac*thk(numelc+neltg)*geo(ipgmat + 2 ,mg)
513 ELSE
514 stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
515 ENDIF
516 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
517 isubstack = iworksh(3,numelc + neltg)
518 IF(slsfac < zero)THEN
519 stf(i)=-slsfac
520 ELSE
521 stf(i)=slsfac*thk(numelc+neltg)*pm_stack(2 ,isubstack)
522 ENDIF
523 ELSE
524 IF(slsfac < zero)THEN
525 stf(i)=-slsfac
526 ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0)THEN
527 stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
528 ELSE
529 stf(i)=slsfac*geo(1,mg)*pm(20,mt)
530 ENDIF
531 ENDIF
532 ELSE
533 IF(nint >= 0) THEN
535 . msgtype=msgwarning,
536 . anmode=aninfo_blind_2,
538 . c1=titr,
539 . i2=ixtg(nixtg,neltg),
540 . c2='SHELL',
541 . i3=i)
542 END IF
543 IF(nint < 0) THEN
545 . msgtype=msgwarning,
546 . anmode=aninfo_blind_2,
548 . c1=titr,
549 . i2=ixtg(nixtg,neltg),
550 . c2='SHELL',
551 . i3=i)
552 END IF
553 END IF
554 IF(igap /= 0 .OR. (nty /= 7 .AND. nty /= 20)) gap_m(i)=gapm
555
556 IF(intfric > 0) THEN
557 ip= iparttg(neltg)
558 ipg = tagprt_fric(ip)
559 IF(ipg > 0) THEN
561 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
562 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
563 ipartfricm(inrt) = ipl
564 ENDIF
565 ENDIF
566
567 cycle
568 ENDIF
569
570
571
572
573 IF(nelc /= 0) THEN
574 mt=ixc(1,nelc)
575 mg=ixc(6,nelc)
576 ip = ipartc(nelc)
577 igtyp = igeo(11,mg)
578 igmat = igeo(98,mg)
579 IF (thk_part(ip) /= zero .AND. iintthick == 0) THEN
580 dx=thk_part(ip)*gapscale
581 ELSEIF (thk(nelc) /= zero .AND. iintthick == 0) THEN
582 dx=thk(nelc)*gapscale
583 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
584 dx=thk(nelc)*gapscale
585 ELSE
586 dx=geo(1,mg)*gapscale
587 ENDIF
588 gapm=half*dx
589 gapm_mx=
max(gapm_mx,gapm)
590 gapmn =
min(gapmn,dx)
591 dxm=dxm+dx
592 ndx=ndx+1
593
594 IF(mt > 0)THEN
595 IF(igtyp == 11 .AND. igmat > 0) THEN
596 IF(slsfac < zero)THEN
597 stf(i)=-slsfac
598 ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0) THEN
599 stf(i)=slsfac*thk(nelc)*geo(ipgmat + 2 ,mg)
600 ELSE
601 stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
602 ENDIF
603 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
604 isubstack = iworksh(3,nelc)
605 IF(slsfac < zero)THEN
606 stf(i)=-slsfac
607 ELSE
608 stf(i)=slsfac*thk(nelc)*pm_stack(2 ,isubstack )
609 ENDIF
610 ELSE
611 IF(slsfac < zero)THEN
612 stf(i)=-slsfac
613 ELSEIF (thk(nelc) /= zero .AND. iintthick == 0) THEN
614 stf(i)=slsfac*thk(nelc)*pm(20,mt)
615 ELSE
616 stf(i)=slsfac*geo(1,mg)*pm(20,mt)
617 ENDIF
618 ENDIF
619 ELSE
620 IF(nint >= 0) THEN
622 . msgtype=msgwarning,
623 . anmode=aninfo_blind_2,
625 . c1=titr,
626 . i2=ixc(nixc,nelc),
627 . c2='SHELL',
628 . i3=i)
629 END IF
630 IF(nint < 0) THEN
632 . msgtype=msgwarning,
633 . anmode=aninfo_blind_2,
635 . c1=titr,
636 . i2=ixc(nixc,nelc),
637 . c2='SHELL',
638 . i3=i)
639 END IF
640 END IF
641 IF(igap /=0 .OR. (nty /=7 .AND. nty /= 20)) gap_m(i)=gapm
642
643 IF(intfric > 0) THEN
644 ip= ipartc(nelc)
645 ipg = tagprt_fric(ip)
646 IF(ipg > 0) THEN
648 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
649 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
650 ipartfricm(inrt) = ipl
651 ENDIF
652 ENDIF
653
654 cycle
655 END IF
656
657
658
659
660 print_error = .false.
661 nin25 = 0
662 CALL insol3d(x ,irect ,ixs ,nint ,nels,inrt,
663 .
area ,noint ,knod2els ,nod2els ,0 ,
664 . ixs10 ,ixs16 ,ixs20 ,tagelems,indexe,ninv,ielem,
665 . elem_linked_to_segment ,print_error,
666 . nin25,nty, flag_elem_inter25 )
667 IF(print_error) THEN
668 node_id(1:4) = itab(irect(1:4,inrt))
669
671 . msgtype=msgwarning,
672 . anmode=aninfo_blind_1,
674 . i2=node_id(1),
675 . i3=node_id(2),
676 . i4=node_id(3),
677 . i5=node_id(4),
678 . c1=titr ,
679 . prmod=msg_print)
680 ENDIF
681
682
683
684
685 IF(nels /= 0) THEN
686 mt=ixs(1,nels)
687 IF(intth > 0 ) ieles(i) = nels
688 IF(mt > 0)THEN
689 DO jj=1,8
690 jjj=ixs(jj+1,nels)
691 xc(jj)=x(1,jjj)
692 yc(jj)=x(2,jjj)
693 zc(jj)=x(3,jjj)
694 ENDDO
696 stf(i)=slsfac*fillsol(nels)*
area*
area*pm(32,mt)/vol
697 ELSE
698 IF(nint >= 0) THEN
700 . msgtype=msgwarning,
701 . anmode=aninfo_blind_2,
703 . c1=titr,
704 . i2=ixs(nixs,nels),
705 . c2='SOLID',
706 . i3=i)
707 ENDIF
708 IF(nint < 0) THEN
710 . msgtype=msgwarning,
711 . anmode=aninfo_blind_2,
713 . c1=titr,
714 . i2=ixs(nixs,nels),
715 . c2='SOLID',
716 . i3=i)
717 ENDIF
718 ENDIF
719
720 IF(intfric > 0) THEN
721 ip= iparts(nels)
722 ipg = tagprt_fric(ip)
723 IF(ipg > 0) THEN
725 . ipg , intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
726 . intbuf_fric_tab(intfric)%TABPARTS_FRIC, ipl )
727 ipartfricm(inrt) = ipl
728 ENDIF
729 ENDIF
730
731 ENDIF
732
733
734
735
736 CALL incoq3(irect ,ixc ,ixtg ,nint ,nelc ,
737 . neltg ,inrt ,geo ,pm ,knod2elc ,
738 . knod2eltg ,nod2elc ,nod2eltg ,thk ,nty ,
739 . igeo ,pm_stack ,iworksh )
740
741
742
743 IF(neltg /= 0) THEN
744 mt=ixtg(1,neltg)
745 mg=ixtg(5,neltg)
746 igtyp = igeo(11,mg)
747 igmat = igeo(98,mg)
748 ip = iparttg(neltg)
749 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
750 dx=thk_part(ip)*gapscale
751 ELSEIF (thk(numelc+neltg) /= zero .AND. iintthick == 0)THEN
752 dx=thk(numelc+neltg)*gapscale
753 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
754 dx=thk(numelc+neltg)*gapscale
755 ELSE
756 dx=geo(1,mg)*gapscale
757 ENDIF
758 gapm=half*dx
759 gapm_mx=
max(gapm_mx,gapm)
760 gapmn =
min(gapmn,dx)
761 dxm=dxm+dx
762 ndx=ndx+1
763 IF(mt > 0)THEN
764 IF(igtyp == 11 .AND. igmat > 0) THEN
765 IF(slsfac < zero)THEN
766 stf(i)=-slsfac
767 ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0) THEN
768 stf(i)=slsfac*thk(numelc+neltg)*geo(ipgmat + 2 ,mg)
769 ELSE
770 stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
771 ENDIF
772 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
773 isubstack = iworksh(3,numelc+neltg)
774 IF(slsfac < zero)THEN
775 stf(i)=-slsfac
776 ELSE
777 stf(i)=slsfac*thk(numelc+neltg)*pm_stack(2 ,isubstack)
778 ENDIF
779 ELSE
780 IF(slsfac < zero)THEN
781 stf(i)=-slsfac
782 ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0) THEN
783 stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
784 ELSE
785 stf(i)=slsfac*geo(1,mg)*pm(20,mt)
786 ENDIF
787 ENDIF
788 ELSE
789 IF(nint >= 0) THEN
791 . msgtype=msgwarning,
792 . anmode=aninfo_blind_2,
794 . c1=titr,
795 . i2=ixtg(nixtg,neltg),
796 . c2='SHELL',
797 . i3=i)
798 ENDIF
799 IF(nint < 0) THEN
801 . msgtype=msgwarning,
802 . anmode=aninfo_blind_2,
804 . c1=titr,
805 . i2=ixtg(nixtg,neltg),
806 . c2='SHELL',
807 . i3=i)
808 ENDIF
809 ENDIF
810
811 IF(intfric > 0) THEN
812 ip= iparttg(neltg)
813 ipg = tagprt_fric(ip)
814 IF(ip > 0) THEN
816 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
817 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
818 ipartfricm(inrt) = ipl
819 ENDIF
820 ENDIF
821
822
823
824
825
826 ELSEIF(nelc /= 0) THEN
827 mt=ixc(1,nelc)
828 mg=ixc(6,nelc)
829 igtyp = igeo(11,mg)
830 igmat = igeo(98,mg)
831 ip = ipartc(nelc)
832 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
833 dx=thk_part(ip)*gapscale
834 ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0) THEN
835 dx=thk(nelc)*gapscale
836 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
837 dx=thk(nelc)*gapscale
838 ELSE
839 dx=geo(1,mg)*gapscale
840 ENDIF
841 gapm=half*dx
842 gapm_mx=
max(gapm_mx,gapm)
843 gapmn =
min(gapmn,dx)
844 dxm=dxm+dx
845 ndx=ndx+1
846 IF(mt > 0)THEN
847 IF(igtyp == 11 .AND. igmat > 0) THEN
848 IF(slsfac < zero)THEN
849 stf(i)=-slsfac
850 ELSEIF (thk(nelc) /= zero .AND. iintthick == 0) THEN
851 stf(i)=slsfac*thk(nelc)*geo(ipgmat + 2 ,mg)
852 ELSE
853 stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
854 ENDIF
855 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
856 isubstack = iworksh(3,nelc)
857 IF(slsfac < zero)THEN
858 stf(i)=-slsfac
859 ELSE
860 stf(i)=slsfac*thk(nelc)*pm_stack(2 ,isubstack)
861 ENDIF
862 ELSE
863 IF(slsfac < zero)THEN
864 stf(i)=-slsfac
865 ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0) THEN
866 stf(i)=slsfac*thk(nelc)*pm(20,mt)
867 ELSE
868 stf(i)=slsfac*geo(1,mg)*pm(20,mt)
869 ENDIF
870 ENDIF
871 ELSE
872 IF(nint >= 0) THEN
874 . msgtype=msgwarning,
875 . anmode=aninfo_blind_2,
877 . c1=titr,
878 . i2=ixc(nixc,nelc),
879 . c2='SHELL',
880 . i3=i)
881 ENDIF
882 IF(nint < 0) THEN
884 . msgtype=msgwarning,
885 . anmode=aninfo_blind_2,
887 . c1=titr,
888 . i2=ixc(nixc,nelc),
889 . c2='SHELL',
890 . i3=i)
891 ENDIF
892 ENDIF
893
894 IF(intfric > 0) THEN
895 ip= ipartc(nelc)
896 ipg = tagprt_fric(ip)
897 IF(ipg > 0) THEN
899 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
900 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
901 ipartfricm(inrt) = ipl
902 ENDIF
903 ENDIF
904
905 ENDIF
906
907 IF(igap /= 0 .OR. (nty /= 7 .AND. nty /= 20)) gap_m(i)=gapm
908
909
910
911
912
913 IF(nels+nelc+neltg == 0)THEN
914
915
916 IF(nint > 0) THEN
918 . msgtype=msgerror,
919 . anmode=aninfo_blind_2,
921 . c1=titr,
922 . i2=i)
923 ENDIF
924 IF(nint < 0) THEN
926 . msgtype=msgerror,
927 . anmode=aninfo_blind_2,
929 . c1=titr,
930 . i2=i)
931 ENDIF
932
933 ENDIF
934
935 enddo
936
937 IF(numels > 0) DEALLOCATE(tagelems,indexe)
938
940 . msgtype=msgwarning,
941 . anmode=aninfo_blind_1,
943 . c1=titr,
944 . prmod=msg_print)
946 . msgtype=msgwarning,
947 . anmode=aninfo_blind_1,
949 . c1=titr,
950 . prmod=msg_print)
951 IF(ninv > 0 .AND.nint>0)
953 . msgtype=msgwarning,
954 . anmode=aninfo_blind_1,
956 . c1=titr,
957 . i2=ninv)
958
959 IF(ninv > 0 .AND.nint< 0)
961 . msgtype=msgwarning,
962 . anmode=aninfo_blind_1,
964 . c1=titr,
965 . i2=ninv)
966
967
968
969 DO i=nrt+1,nrt+nrt_ige
970 stf(i)=zero
971 IF(intth > 0) ieles(i) = 0
972 IF(slsfac < zero)stf(i)=slsfac
973 gapm =zero
974 inrt=i
975 CALL i4gmx3(x,irect,inrt,gapmx)
976
977
978
979 CALL ineltigeo(x ,irect ,ixs ,nint ,nelig3d ,
980 . inrt ,
area ,noint ,0 ,igrsurf%ELTYP_IGE,
981 . ixig3d ,kxig3d ,igeo ,igrsurf%ELEM_IGE)
982 IF(nelig3d /= 0)THEN
983 mt=kxig3d(1,nelig3d)
984 IF(mt > 0)THEN
985 ipid = kxig3d(2,nelig3d)
986 px = igeo(41,ipid)-1
987 py = igeo(42,ipid)-1
988 pz = igeo(43,ipid)-1
989 coin_ige(1) = (px+1)*py+1
990 coin_ige(2) = (px+1)*(py+1)
991 coin_ige(3) = px+1
992 coin_ige(4) = 1
993 coin_ige(5) = (px+1)*(py+1)*pz+(px+1)*py+1
994 coin_ige(6) = (px+1)*(py+1)*(pz+1)
995 coin_ige(7) = (px+1)*(py+1)*pz+px+1
996 coin_ige(8) = (px+1)*(py+1)*pz+1
997 DO jj=1,8
998 xc(jj)=x(1,ixig3d(kxig3d(4,nelig3d)+coin_ige(jj)-1))
999 yc(jj)=x(2,ixig3d(kxig3d(4,nelig3d)+coin_ige(jj)-1))
1000 zc(jj)=x(3,ixig3d(kxig3d(4,nelig3d)+coin_ige(jj)-1))
1001 END DO
1003 stf(i)=slsfac*
area*
area*pm(32,mt)/vol
1004 stf(i)=stf(i)*((px+1)*(py+1)+(py+1)*(pz+1)+(pz+1)*(px+1))/3
1005 ELSE
1006 IF(nint >= 0) THEN
1008 . msgtype=msgwarning,
1009 . anmode=aninfo_blind_2,
1011 . c1=titr,
1012 . i2=kxig3d(5,nelig3d),
1013 . c2='ISOGEOMETRIC SOLID',
1014 . i3=i)
1015 ENDIF
1016 IF(nint < 0) THEN
1018 . msgtype=msgwarning,
1019 . anmode=aninfo_blind_2,
1021 . c1=titr,
1022 . i2=kxig3d(5,nelig3d),
1023 . c2='ISOGEOMETRIC SOLID',
1024 . i3=i)
1025 ENDIF
1026 ENDIF
1027 ELSEIF(nelig3d == 0)THEN
1028
1029
1030 IF(nint > 0) THEN
1032 . msgtype=msgerror,
1033 . anmode=aninfo_blind_2,
1035 . c1=titr,
1036 . i2=i)
1037 ENDIF
1038 IF(nint < 0) THEN
1040 . msgtype=msgerror,
1041 . anmode=aninfo_blind_2,
1043 . c1=titr,
1044 . i2=i)
1045 ENDIF
1046
1047 ENDIF
1048 enddo
1049
1050
1051
1052
1053 gapmx=sqrt(gapmx)
1054 IF(igap == 0)THEN
1055
1056 IF(gap <= zero)THEN
1057 IF(ndx /= 0)THEN
1058 gap = dxm/ndx
1059 gap =
min(half*gapmx,gap)
1060 ELSE
1061 gap = em01 * gapmx
1062 ENDIF
1063 IF (it19 <= 0 .AND. .NOT.type18) WRITE(iout,1300)gap
1064 ENDIF
1065 gapmin = gap
1066
1067 IF (gapmin <= 0) THEN
1069 . msgtype=msgerror,
1070 . anmode=aninfo,
1072 . c1=titr)
1073 ENDIF
1074 IF ((inacti /= 7).AND.(gap > 0.5*gapmx) .AND. (irem_gap /= 2)) THEN
1075 gaptmp = half*gapmx
1077 . msgtype=msgwarning,
1078 . anmode=aninfo_blind_2,
1080 . c1=titr,
1081 . r1=gap,
1082 . r2=gaptmp)
1083 ENDIF
1084 ELSE
1085
1086
1087
1088 IF(gap <= zero)THEN
1089 IF(ndx /= 0)THEN
1090 gapmin = gapmn
1091 gapmin =
min(half*gapmx,gapmin)
1092 ELSE
1093 gapmin = em01 * gapmx
1094 ENDIF
1095 IF (gapmin <= 0) THEN
1097 . msgtype=msgerror,
1098 . anmode=aninfo,
1100 . c1=titr)
1101 ENDIF
1102 IF (it19 <= 0 .AND. .NOT.type18) WRITE(iout,1300)gapmin
1103 ELSE
1104 gapmin = gap
1105 ENDIF
1106
1107 IF(igap == 3) THEN
1108 gap =
max(
min(gaps_mx+gapm_mx,gaps_l_mx+gapm_l_mx) ,gapmin)
1109 ELSE
1110 gap =
max(gaps_mx+gapm_mx,gapmin)
1111 ENDIF
1113 IF ((igap /= 3).AND.(irem_gap /= 2)) THEN
1114 IF(inacti /= 7.AND.gap > half*gapmx .AND. iddlevel == 1)THEN
1115 gaptmp = 0.5*gapmx
1117 . msgtype=msgwarning,
1118 . anmode=aninfo_blind_2,
1120 . c1=titr,
1121 . r1=gap)
1122 ENDIF
1123 ENDIF
1124 ENDIF
1125
1126 IF(intth /= 0)THEN
1127 IF(drad == zero)THEN
1128
1130 ELSEIF(drad < gap)THEN
1131
1132 drad=gap
1133 END IF
1134 WRITE(iout,2001)drad
1135
1136
1137 IF(drad > gapmx)THEN
1139 . msgtype=msgwarning,
1140 . anmode=aninfo_blind_2,
1142 . c1=titr,
1143 . r1=drad ,
1144 . r2=gapmx,
1146 END IF
1147 END IF
1148
1149
1150 IF(intfric > 0) THEN
1151 IF(numels /= 0)THEN
1152 DO i = 1,nsn
1153 ipfmax = 0
1154 ipflmax = 0
1155 DO j= knod2els(nsv(i))+1,knod2els(nsv(i)+1)
1156 ie = nod2els(j)
1157 ip = iparts(ie)
1158 ipg = tagprt_fric(ip)
1159 IF(ipg > 0 .AND. ip > ipfmax) THEN
1161 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
1162 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
1163 IF(ipl /= 0) THEN
1164 ipfmax = ip
1165 ipflmax = ipl
1166 ENDIF
1167 ENDIF
1168 ENDDO
1169 IF(ipfmax /= 0) THEN
1170 ipartfrics(i) = ipflmax
1171 ENDIF
1172
1173 ENDDO
1174 ENDIF
1175
1176 IF(numelc /= 0 .OR. numeltg /= 0) THEN
1177 DO i = 1,nsn
1178 ipfmax = 0
1179 ipflmax = 0
1180 DO j= knod2elc(nsv(i))+1,knod2elc(nsv(i)+1)
1181 ie = nod2elc(j)
1182 ip = ipartc(ie)
1183 ipg = tagprt_fric(ip)
1184 IF(ipg > 0 .AND. ip > ipfmax) THEN
1186 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
1187 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
1188 IF(ipl /= 0) THEN
1189 ipfmax = ip
1190 ipflmax = ipl
1191 ENDIF
1192 ENDIF
1193 ENDDO
1194
1195 DO j= knod2eltg(nsv(i))+1,knod2eltg(nsv(i)+1)
1196 ie = nod2eltg(j)
1197 ip = iparttg(ie)
1198 ipg = tagprt_fric(ip)
1199 IF(ipg > 0.AND.ip > ipfmax) THEN
1201 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
1202 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
1203
1204 IF(ipl /= 0) THEN
1205 ipfmax = ip
1206 ipflmax = ipl
1207 ENDIF
1208 ENDIF
1209 ENDDO
1210 IF(ipfmax /= 0) THEN
1211 ipartfrics(i) = ipflmax
1212 ENDIF
1213
1214 ENDDO
1215 ENDIF
1216 ENDIF
1217
1218
1219
1220
1221
1222 DO l=1,nsn
1223 stfn(l) = 1.
1224 ENDDO
1225
1226
1227
1228 bgapsmx = zero
1229 IF (igap == 0) THEN
1230 gapinf=gap
1231 ELSE
1232 gapinf = ep30
1233 gapinf_s = ep30
1234 gapinf_m = ep30
1235 DO i = 1, nsn
1236 gapinf_s =
min(gapinf_s,gap_s(i))
1237 bgapsmx =
max(bgapsmx,gap_s(i))
1238 ENDDO
1239 DO i = 1, nrt+nrt_ige
1240 gapinf_m =
min(gapinf_m,gap_m(i))
1241 ENDDO
1242 gapinf=gapinf_s+gapinf_m
1243 gapinf=
max(gapinf,gapmin)
1244 ENDIF
1245 DEALLOCATE( gap_s_l_tmp )
1246 RETURN
1247 1300 FORMAT(2x,'GAP MIN = ',1pg20.13)
1248 2001 FORMAT(2x,'Maximum distance for radiation computation = ',
1249 . 1pg20.13)
subroutine area(d1, x, x2, y, y2, eint, stif0)
subroutine i4gmx3(x, irect, i, gapmax)
subroutine friction_parts_search(ip, npartsfric, partsfric, ipl)
subroutine incoq3(irect, ixc, ixtg, nint, nel, neltg, is, geo, pm, knod2elc, knod2eltg, nod2elc, nod2eltg, thk, nty, igeo, pm_stack, iworksh)
subroutine inelts(x, irect, ixs, nint, nel, i, area, noint, ir, surf_eltyp, surf_elem)
subroutine ineltigeo(xe, irect, ixs, nint, nelig3d, i, area, noint, ir, surf_eltyp_ige, ixig3d, kxig3d, igeo, surf_elem_ige)
subroutine ineltc(nelc, neltg, is, surf_eltyp, surf_elem)
subroutine insol3d(x, irect, ixs, nint, nel, i, area, noint, knod2els, nod2els, ir, ixs10, ixs16, ixs20, tagelems, indexe, ninv, ielem_m, elem_linked_to_segment, print_error, nin25, nty, flag_elem_inter25)
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)