56
57
58
59
60
61
62
63
64
65
66
67
68
71 USE intbuf_fric_mod
73
74
75
76#include "implicit_f.inc"
77
78
79
80#include "units_c.inc"
81#include "param_c.inc"
82#include "com01_c.inc"
83#include "com04_c.inc"
84#include "scr08_c.inc"
85#include "scr23_c.inc"
86#include "r2r_c.inc"
87
88
89
90 INTEGER,INTENT(IN) :: NRT, NINTR, NTY, NOINT,IGAP,INTTH,INTFRIC,IREM_GAP
91 INTEGER,INTENT(IN) :: IRECT(2,*), IXS(NIXS,NUMELS), IXC(NIXC,NUMELC),
92 . IXTG(NIXTG,NUMELTG),IXT(NIXT,NUMELT),IXP(NIXP,NUMELP),IXR(NIXR,NUMELR),
93 . IPARTC(NUMELC), IPARTTG(NUMELTG),NOD2EL1D(*),KNOD2EL1D(*),ITAB(NUMNOD),
94 . IXS10(6,*),KXX(NIXX,*),IXX(*),IGEO(NPROPGI,NUMGEO),
95 . KNOD2ELS(*), KNOD2ELC(*), KNOD2ELTG(*),
96 . NOD2ELS(*), NOD2ELC(*), NOD2ELTG(*),IWORKSH(3,*),
97 . TAGPRT_FRIC(*),IPARTS(*),
98 . IPARTT(*) ,IPARTP(*) ,IPARTX(*) ,IPARTR(*)
99 INTEGER,INTENT(INOUT) :: IPARTFRIC(*),IELEC(*)
100 INTEGER,INTENT(IN) :: ID,IT19
101 my_real,
INTENT(IN) :: x(3,numnod), pm(npropm,nummat), geo(npropg,numgeo),
102 . ms(*),thk(*),thk_part(*),
103 . lelx(*), fillsol(*),pm_stack(20,*)
104 my_real,
INTENT(IN) :: slsfac,gap0,percent_size
105 my_real,
INTENT(INOUT) :: gap_l(*),stf(*),gap_sm(*),
area(*),drad,gapinf,gapmax,gapmin,bgapsmx
106 CHARACTER(LEN=NCHARTITLE) ,INTENT(IN) :: TITR
107 TYPE(INTBUF_FRIC_STRUCT_),INTENT(IN) :: INTBUF_FRIC_TAB(*)
108
109
110
111 INTEGER NDX, I, , NELS, MT, JJ, JJJ, NELC, J,
112 . MG, NELTG,NELT,NELP,NELR,
113 . IGTYP, IP,N1,N2,K,T,P,R,NELX,IPGMAT,IGMAT,IE,
114 . JJ1,JJ2,IEC,K1,K2,IPL,IPC,
115 . IE1(50,2),IE2(50,2),ISUBSTACK,IPG,N3,N4,N5,N6,N7,N8,ICONTR
116 my_real dxm, gapmx, gapmn, areass, vol, dx,gap1,gaps1,gaptmp,xl,
117 . sx1,sy1,sz1,sx2,sy2,sz2,sx3,sy3,sz3,
118 . x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4,x5,y5,z5,x6,y6,z6,x7,y7,z7,x8,y8,z8,
119 . xx1(4),xx2(4), xx3(4) ,face(6),
120 . n_1, n_2, n_3,dx1,xl2,bulk
121
123
124
125
126 dxm=zero
127 ndx=0
128 gaps1=zero
129 gapmx=ep30
130 gapmn=ep30
131 ipgmat = 700
132 IF(igap == 3)THEN
133 DO i=1,nrt
134 gap_l(i)=ep30
135 ENDDO
136 ENDIF
137
138 DO i=1,nrt
139 stf(i)=zero
140 gap_sm(i)=zero
141 inrt=i
142 CALL i11gmx3(x,irect,inrt,gapmx,xl2)
143
144
145
146 CALL i11sol(x,irect,ixs,nintr,nels,inrt,areass,noint,knod2els,nod2els,ixs10)
147 IF(nels /= 0) THEN
148 mt=ixs(1,nels)
149 mg=ixs(nixs-1,nels)
150 icontr = igeo(97,mg)
151 IF(mt > 0)THEN
152 DO jj=1,8
153 jjj=ixs(jj+1,nels)
154 xc(jj)=x(1,jjj)
155 yc(jj)=x(2,jjj)
156 zc(jj)=x(3,jjj)
157 ENDDO
159 IF (icontr==1 ) THEN
160 bulk = pm(107,mt)
161 ELSE
162 bulk = pm(32,mt)
163 END IF
164 IF(xl2 > 0.0)THEN
165 stf(i)=slsfac*fillsol(nels)*vol*bulk/xl2
166 ELSE
167 stf(i)=zero
168 ENDIF
169 ELSE
170 IF(nintr >= 0) THEN
171 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
173 . c1=titr,
174 . i2=ixs(nixs,nels),
175 . c2='SOLID',
176 . i3=i)
177 ENDIF
178 IF(nintr < 0) THEN
179 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
181 . c1=titr,
182 . i2=ixs(nixs,nels),
183 . c2='SOLID',
184 . i3=i)
185 ENDIF
186 ENDIF
187
188
189 IF(intfric > 0) THEN
190 ip= iparts(nels)
191 ipg = tagprt_fric(ip)
192 IF(ipg > 0) THEN
194 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
195 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
196 ipartfric(inrt) = ipl
197 ENDIF
198 ENDIF
199
200 ENDIF
201 CALL i11coq(irect ,ixc ,ixtg ,nintr ,nelc ,
202 . neltg ,inrt ,geo ,pm ,thk ,
203 . igeo ,knod2elc,knod2eltg,nod2elc,nod2eltg,
204 . pm_stack,iworksh)
205 IF(neltg /= 0) THEN
206 mt=ixtg(1,neltg)
207 mg=ixtg(5,neltg)
208 igtyp = igeo(11,mg)
209 ip = iparttg(neltg)
210 igmat = igeo(98,mg)
211 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
212 dx=thk_part(ip)
213 ELSEIF (thk(numelc+neltg) /= zero .AND. iintthick ==0)THEN
214 dx=thk(numelc+neltg)
215 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
216 dx=thk(numelc+neltg)
217 ELSE
218 dx=geo(1,mg)
219 ENDIF
220 gap_sm(i)=half*dx
221 gaps1=
max(gaps1,gap_sm(i))
222 gapmn =
min(gapmn,dx)
223 dxm=dxm+dx
224 ndx=ndx+1
225 IF(mt > 0)THEN
226 IF(igtyp == 11 .AND. igmat > 0) THEN
227 stf(i)=slsfac*dx*geo(ipgmat + 2 ,mg)
228 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0)) THEN
229 isubstack = iworksh(3,numelc+neltg)
230 stf(i)=slsfac*dx*pm_stack( 2 ,isubstack)
231 ELSE
232 stf(i)=slsfac*dx*pm(20,mt)
233 ENDIF
234 ELSE
235 IF(nintr >= 0) THEN
236 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
238 . c1=titr,
239 . i2=ixs(nixs,nels),
240 . c2='SOLID',
241 . i3=i)
242 ENDIF
243 IF(nintr < 0) THEN
244 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
246 . c1=titr,
247 . i2=ixs(nixs,nels),
248 . c2='SOLID',
249 . i3=i)
250 ENDIF
251
252 ENDIF
253
254
255 IF(intfric > 0) THEN
256 ip= iparttg(neltg)
257 ipg = tagprt_fric(ip)
258 IF(ipg > 0) THEN
260 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
261 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
262 ipartfric(inrt) = ipl
263 ENDIF
264 ENDIF
265
266
267 ELSEIF(nelc /= 0) THEN
268 mt=ixc(1,nelc)
269 mg=ixc(6,nelc)
270 igtyp = igeo(11,mg)
271 ip = ipartc(nelc)
272 igmat = igeo(98,mg)
273 IF (thk_part(ip) /= zero .AND. iintthick == 0) THEN
274 dx=thk_part(ip)
275 ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0 ) THEN
276 dx=thk(nelc)
277 ELSEIF(igtyp ==17 .OR. igtyp == 51 .OR. igtyp == 52)THEN
278 dx=thk(nelc)
279 ELSE
280 dx=geo(1,mg)
281 ENDIF
282 gap_sm(i)=half*dx
283 gaps1=
max(gaps1,gap_sm(i))
284 gapmn =
min(gapmn,dx)
285 dxm=dxm+dx
286 ndx=ndx+1
287 IF(mt > 0)THEN
288 IF(igtyp == 11 .AND. igmat > 0) THEN
289 stf(i)=slsfac*dx*geo(ipgmat + 2 ,mg)
290 ELSEIF(igtyp == 52 .OR. ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0)) THEN
291 isubstack = iworksh(3,nelc)
292 stf(i)=slsfac*dx*pm_stack( 2 ,isubstack)
293 ELSE
294 stf(i)=slsfac*dx*pm(20,mt)
295 ENDIF
296 ELSE
297 IF(nintr >= 0) THEN
298 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
300 . c1=titr,
301 . i2=ixc(nixc,nelc),
302 . c2='SHELL',
303 . i3=i)
304 ENDIF
305 IF(nintr < 0) THEN
306 CALL ancmsg(msgid=96,msgtype=msgwarning, anmode=aninfo_blind_2,
308 . c1=titr,
309 . i2=ixc(nixc,nelc),
310 . c2='SHELL',
311 . i3=i)
312 ENDIF
313 ENDIF
314
315
316 IF(intfric > 0) THEN
317 ip= ipartc(nelc)
318
319 ipg = tagprt_fric(ip)
320 IF(ipg > 0) THEN
322 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
323 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
324 ipartfric(inrt) = ipl
325 ENDIF
326 ENDIF
327
328 ENDIF
329 CALL i11fil(irect,ixt ,ixp,ixr,nintr,nelt ,
330 . nelp,nelr,nelx,inrt,nod2el1d,
331 . knod2el1d,kxx,ixx)
332 IF(nelt /= 0) THEN
333 mt=ixt(1,nelt)
334 mg=ixt(4,nelt)
335 ip = ipartt(nelt)
336 IF (thk_part(ip) > zero ) THEN
337 dx1=thk_part(ip)
338 ELSE
339 dx1=sqrt(geo(1,mg))
340 END IF
341 dx=sqrt(geo(1,mg))
342 gap_sm(i)=
max(gap_sm(i),half*dx1)
343 gaps1=
max(gaps1,gap_sm(i))
344 gapmn =
min(gapmn,dx1)
345 dxm=dxm+dx
346 ndx=ndx+1
347 IF(mt > 0)THEN
348 stf(i)=slsfac*dx*pm(20,mt)
349 IF (nsubdom>0) stf(i)=slsfac*dx*pm_r2r(mt)
350 ELSE
351 IF(nintr >= 0) THEN
352 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
354 . c1=titr,
355 . i2=ixt(nixt,nelt),
356 . c2='TRUSS',
357 . i3=i)
358 ENDIF
359 IF(nintr < 0) THEN
360 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
362 . c1=titr,
363 . i2=ixt(nixt,nelt),
364 . c2='TRUSS',
365 . i3=i)
366 ENDIF
367 ENDIF
368
369 IF(intfric > 0) THEN
370 ip= ipartt(nelt)
371 ipg = tagprt_fric(ip)
372 IF(ipg > 0) THEN
374 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
375 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
376 ipartfric(inrt) = ipl
377 ENDIF
378 ENDIF
379
380 ELSEIF(nelp /= 0) THEN
381 mt=ixp(1,nelp)
382 mg=ixp(5,nelp)
383 ip = ipartp(nelp)
384 IF (thk_part(ip) > zero ) THEN
385 dx1=thk_part(ip)
386 ELSE
387 dx1=sqrt(geo(1,mg))
388 END IF
389 dx=sqrt(geo(1,mg))
390 gap_sm(i)=
max(gap_sm(i),half*dx1)
391 gaps1=
max(gaps1,gap_sm(i))
392 gapmn =
min(gapmn,dx1)
393 dxm=dxm+dx
394 ndx=ndx+1
395 IF(mt > 0)THEN
396 stf(i)=slsfac*dx*pm(20,mt)
397 IF (nsubdom>0) stf(i)=slsfac*dx*pm_r2r(mt)
398 ELSE
399 IF(nintr >= 0) THEN
402 . c1=titr,
403 . i2=ixp(nixp,nelp),
404 . c2='BEAM',
405 . i3=i)
406 ENDIF
407 IF(nintr<0) THEN
408 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
410 . c1=titr,
411 . i2=ixp(nixp,nelp),
412 . c2='beam',
413 . I3=I)
414 ENDIF
415 ENDIF
416
417 IF(INTFRIC > 0) THEN
418 IP= IPARTP(NELP)
419 IPG = TAGPRT_FRIC(IP)
420 IF(IPG > 0) THEN
421 CALL FRICTION_PARTS_SEARCH (
422 . IPG,INTBUF_FRIC_TAB(INTFRIC)%S_TABPARTS_FRIC,
423 . INTBUF_FRIC_TAB(INTFRIC)%TABPARTS_FRIC,IPL )
424 IPARTFRIC(INRT) = IPL
425 ENDIF
426 ENDIF
427
428 ELSEIF(NELR /= 0) THEN
429 MG=IXR(1,NELR)
430 MT = IXR(5,NELR)
431 IP = IPARTR(NELR)
432 IF(IP > 0)THEN
433 IF (THK_PART(IP) > ZERO ) THEN
434 DX1=THK_PART(IP)
435 GAP_SM(I)=MAX(GAP_SM(I),HALF*DX1)
436 GAPS1=MAX(GAPS1,GAP_SM(I))
437 GAPMN = MIN(GAPMN,DX1)
438 END IF
439 ENDIF
440 IF(MG > 0)THEN
441 IGTYP=NINT(GEO(12,MG))
442.OR. IF(IGTYP==4IGTYP==12)THEN
443 STF(I)=SLSFAC*GEO(2,MG)
444.OR. ELSEIF(IGTYP==8IGTYP==13)THEN
445 STF(I)=SLSFAC*MAX(GEO(3,MG),GEO(10,MG),GEO(15,MG))
446 ELSEIF(IGTYP == 23)THEN
447 STF(I)=SLSFAC*MAX(PM(191,MT),PM(192,MT),PM(193,MT))
448 ELSEIF(IGTYP==25)THEN
449 STF(I)=SLSFAC*GEO(10,MG)
450 ELSEIF(IGTYP>=29)THEN
451 STF(I)=SLSFAC*GEO(3,MG)
452 ELSE
453 WRITE(6,'(a)') 'INTERNAL ERROR 987'
455 ENDIF
456 ELSE
457 IF(nintr >= 0) THEN
458 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
460 . c1=titr,
461 . i2=ixr(nixr,nelr),
462 . c2='SPRING',
463 . i3=i)
464 ENDIF
465 IF(nintr < 0) THEN
466 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
468 . c1=titr,
469 . i2=ixr(nixr,nelr),
470 . c2='SPRING',
471 . i3=i)
472 ENDIF
473 ENDIF
474
475 IF(intfric > 0) THEN
476 ip= ipartr(nelr)
477 ipg = tagprt_fric(ip)
478 IF(ipg > 0) THEN
480 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
481 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
482 ipartfric(inrt) = ipl
483 ENDIF
484 ENDIF
485
486 ELSEIF(nelx /= 0) THEN
487 mg=kxx(2,nelx)
488 IF(mg>0)THEN
489 stf(i)=slsfac*get_u_geo(4,mg)*(kxx(3,nelx)-1)/lelx(nelx)
490 ELSE
491 IF(nintr >= 0) THEN
492 CALL ancmsg(msgid=95,msgtype=msgwarning,anmode=aninfo_blind_2,
494 . c1=titr,
495 . i2=kxx(nixx,nelx),
496 . c2='XELEM',
497 . i3=i)
498 ENDIF
499 IF(nintr < 0) THEN
500 CALL ancmsg(msgid=96,msgtype=msgwarning,anmode=aninfo_blind_2,
502 . c1=titr,
503 . i2=kxx(nixx,nelx),
504 . c2='XELEM',
505 . i3=i)
506 ENDIF
507 ENDIF
508
509 IF(intfric > 0) THEN
510 ip= ipartx(nelx)
511 ipg = tagprt_fric(ip)
512 IF(ipg > 0) THEN
514 . ipg,intbuf_fric_tab(intfric)%S_TABPARTS_FRIC,
515 . intbuf_fric_tab(intfric)%TABPARTS_FRIC,ipl )
516 ipartfric(inrt) = ipl
517 ENDIF
518 ENDIF
519
520 ENDIF
521
522 IF(nels+nelc+neltg+nelt+nelp+nelr+numelx==0.)THEN
523 IF(nintr > 0) THEN
524 CALL ancmsg(msgid=481,msgtype=msgerror,anmode=aninfo_blind_2,
526 . c1=titr,
527 . i2=i)
528 ENDIF
529 IF(nintr < 0) THEN
530 CALL ancmsg(msgid=482,msgtype=msgerror,anmode=aninfo_blind_2,
532 . c1=titr,
533 . i2=i)
534 ENDIF
535 ENDIF
536 enddo
537
538
539
540
541 IF(igap == 3)THEN
542 DO i=1,nrt
543 xl = ep30
544 n1=irect(1,i)
545 n2=irect(2,i)
546 IF(n1 /= n2 .AND. n1 /= 0)
547 . xl=
min(xl,sqrt((x(1,n1)-x(1,n2))**2+(x(2,n1)-x(2,n2))**2+
548 . (x(3,n1)-x(3,n2))**2))
549
550 DO j=1,2
551 n1=irect(j,i)
552 DO k=knod2el1d(n1)+1,knod2el1d(n1+1)
553 IF (nod2el1d(k) <= numelt .AND. nod2el1d(k) /= zero) THEN
554 t = nod2el1d(k)
555 xl=
min(xl,sqrt((x(1,ixt(2,t))-x(1,ixt(3,t)))**2+
556 . (x(2,ixt(2,t))-x(2,ixt(3,t)))**2+
557 . (x(3,ixt(2,t))-x(3,ixt(3,t)))**2))
558 ELSEIF (nod2el1d(k) <= numelt+numelp .AND. nod2el1d(k) /= zero) THEN
559 p = nod2el1d(k) - numelt
560 xl=
min(xl,sqrt((x(1,ixp(2,p))-x(1,ixp(3,p)))**2+
561 . (x(2,ixp(2,p))-x(2,ixp(3,p)))**2+
562 . (x(3,ixp(2,p))-x(3,ixp(3,p)))**2))
563 ELSEIF (nod2el1d(k) <= numelt+numelp+numelr .AND. nod2el1d(k) /= zero) THEN
564 r = nod2el1d(k) - numelt - numelp
565 xl=
min(xl,sqrt((x(1,ixr(2,r))-x(1,ixr(3,r)))**2+
566 . (x(2,ixr(2,r))-x(2,ixr(3,r)))**2+
567 . (x(3,ixr(2,r))-x(3,ixr(3,r)))**2))
568 ENDIF
569 ENDDO
570 ENDDO
571 DO j=1,2
572 gap_l(i) =
min(gap_l(i),percent_size*xl)
573 ENDDO
574 ENDDO
575 ENDIF
576
577
578
579 gapmx=sqrt(gapmx)
580 IF(igap==0)THEN
581
582
583
584 IF(gap0 > zero)THEN
585 gap1 = gap0
586 ELSE
587 IF(ndx/=0)THEN
588 gap1 =
min(half*gapmx,dxm/ndx)
589 ELSE
590 gap1 = em01* gapmx
591 ENDIF
592 IF ((nintr<0).AND.(it19==0)) WRITE(iout,1300)half*(gapmin+gap1)
593 ENDIF
594
595 IF(nintr < 0) gap1 = half*(gapmin+gap1)
596 gapmin = gap1
597 gapmax = gap1
598
599 IF ((gap1 > half*gapmx) .AND. (irem_gap /= 2)) THEN
600 gaptmp = half*gapmx
601 CALL ancmsg(msgid=94,msgtype=msgwarning,anmode=aninfo_blind_2,
603 . c1=titr,
604 . r1=gap1,
605 . r2=gaptmp)
606 ENDIF
607 ELSE
608
609
610
611 bgapsmx = zero
612 IF(gap0 > zero)THEN
613 gap1 = gap0
614 ELSE
615 IF(ndx /= 0)THEN
616 gap1 =
min(half*gapmx,gapmn)
617 ELSE
618 gap1 = em01 * gapmx
619 ENDIF
620 IF ((nintr<0) .AND. (it19==0)) WRITE(iout,1300)half*(gapmin+gap1)
621 ENDIF
622
623 IF(nintr>0)THEN
624 gapmin = gap1
625 gapmax = gaps1
626 ELSE
627 gapmin = half*(gapmin+gap1)
628 gapmax =
max(gapmax+gaps1,gapmin)
629 bgapsmx =
max(bgapsmx,gaps1)
630 ENDIF
631
632 IF ((gapmax>half*gapmx) .AND. (igap/=3) .AND. (irem_gap/=2))THEN
633 gaptmp = half*gapmx
634 CALL ancmsg(msgid=94,msgtype=msgwarning,anmode=aninfo_blind_2,
636 . c1=titr,
637 . r1=gapmax,
638 . r2=gaptmp)
639 ENDIF
640 ENDIF
641
642
643
644 IF(slsfac < zero)THEN
645 DO i=1,nrt
646 stf(i)=-slsfac
647 ENDDO
648 ENDIF
649
650
651
652
653 IF(igap == 0) THEN
654 gapinf=gapmax
655 ELSEIF(igap==1 .OR. igap==2) THEN
656 DO i = 1, nrt
657 gapinf =
min(gapinf,gap_sm(i))
658 ENDDO
659 ELSEIF(igap==3) THEN
660 DO i = 1, nrt
661 gapinf =
min(gapinf,
min(gap_sm(i),gap_l(i)))
662 ENDDO
663 ENDIF
664
665 IF(intth /= 0)THEN
666 IF(drad == zero)THEN
667
669 ELSEIF(drad < gap1)THEN
670
671 drad=gap1
672 END IF
673 WRITE(iout,2001)drad
674
675
676 IF(drad > gapmx)THEN
677 CALL ancmsg(msgid=918, msgtype=msgwarning, anmode
679 . c1=titr,
680 . r1=drad ,
681 . r2=gapmx,
683 END IF
684 END IF
685
686
687 IF(intth > 0 ) THEN
688 IF(nelc /=0 .OR. neltg /= 0) THEN
689 DO i=1,nrt
691 jj1 = 0
692 DO j=knod2elc(irect(1,i))+1,knod2elc(irect(1,i)+1)
693 jj1 = jj1 +1
694 ie1(jj1,1) = nod2elc(j)
695 ie1(jj1,2) = 4
696 ENDDO
697 DO j= knod2eltg(irect(1,i))+1,knod2eltg(irect(1,i)+1)
698 jj1 = jj1 +1
699 ie1(jj1,1) = nod2eltg(j)
700 ie1(jj1,2) = 3
701 ENDDO
702 jj2 = 0
703 DO j=knod2elc(irect(2,i))+1,knod2elc(irect(2,i)+1)
704 jj2 = jj2 +1
705 ie2(jj2,1) = nod2elc(j)
706 ie2(jj2,2) = 4
707 ENDDO
708 DO j= knod2eltg(irect(2,i))+1,knod2eltg(irect(2,i)+1)
709 jj2 = jj2 +1
710 ie2(jj2,1) = nod2eltg(j)
711 ie2(jj2,2) = 3
712 ENDDO
713 iec = 0
714 DO j=1,jj1
715 DO k=1,jj2
716 IF(ie1(j,1) == ie2(k,1)) THEN
717 ie = ie1(j,1)
718 iec = iec +1
719 IF(ie1(j,2)==4) THEN
720 sx1 = x(1,ixc(4,ie)) - x(1,ixc(2,ie))
721 sy1 = x(2,ixc(4,ie)) - x(2,ixc(2,ie))
722 sz1 = x(3,ixc(4,ie)) - x(3,ixc(2,ie))
723 sx2 = x(1,ixc(5,ie)) - x(1,ixc(3,ie))
724 sy2 = x(2,ixc(5,ie)) - x(2,ixc(3,ie))
725 sz2 = x(3,ixc(5,ie)) - x(3,ixc(3,ie))
726 sx3 = sy1*sz2 - sz1*sy2
727 sy3 = sz1*sx2 - sx1*sz2
728 sz3 = sx1*sy2 - sy1*sx2
729 area(i) =
area(i)+ one_over_8*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
730 ielec(i) = ixc(1,ie)
731 mg = ixc(6,ie)
732 ENDIF
733 IF(ie1(j,2) == 3) THEN
734 sx1 = x(1,ixtg(3,ie)) - x(1,ixtg(2,ie))
735 sy1 = x(2,ixtg(3,ie)) - x(2,ixtg(2,ie))
736 sz1 = x(3,ixtg(3,ie)) - x(3,ixtg(2,ie))
737 sx2 = x(1,ixtg(4,ie)) - x(1,ixtg(2,ie))
738 sy2 = x(2,ixtg(4,ie)) - x(2,ixtg(2,ie))
739 sz2 = x(3,ixtg(4,ie)) - x(3,ixtg(2,ie))
740 sx3 = sy1*sz2 - sz1*sy2
741 sy3 = sz1*sx2 - sx1*sz2
742 sz3 = sx1*sy2 - sy1*sx2
743 area(i) =
area(i) + one_over_6*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
744 ielec(i) = ixtg(1,ie)
745 mg = ixtg(5,ie)
746 ENDIF
747 ENDIF
748 ENDDO
749 ENDDO
751 IF(iec == 1) THEN
752 dx=geo(1,mg)
754 ENDIF
755 ENDDO
756
757 ELSEIF(nels /= 0) THEN
758 DO i=1,nrt
760 jj1 = 0
761 DO j=knod2els(irect(1,i))+1,knod2els(irect(1,i)+1)
762 jj1 = jj1 +1
763 ie1(jj1,1) = nod2els(j)
764 ENDDO
765 jj2 = 0
766 DO j=knod2els(irect(2,i))+1,knod2els(irect(2,i)+1)
767 jj2 = jj2 +1
768 ie2(jj2,1) = nod2els(j)
769 ENDDO
770 iec = 0
771 DO j=1,jj1
772 DO k=1,jj2
773 IF(ie1(j,1) == ie2(k,1)) THEN
774 ie = ie1(j,1)
775 iec= iec +1
776 ielec(i) = ixs(1,ie)
777
778 n1=ixs(2,ie)
779 n2=ixs(3,ie)
780 n3=ixs(4,ie)
781 n4=ixs(5,ie)
782 n5=ixs(6,ie)
783 n6=ixs(7,ie)
784 n7=ixs(8,ie)
785 n8=ixs(9,ie)
786
787 x1=x(1,n1)
788 y1=x(2,n1)
789 z1=x(3,n1)
790 x2=x(1,n2)
791 y2=x(2,n2)
792 z2=x(3,n2)
793 x3=x(1,n3)
794 y3=x(2,n3)
795 z3=x(3,n3)
796 x4=x(1,n4)
797 y4=x(2,n4)
798 z4=x(3,n4)
799 x5=x(1,n5)
800 y5=x(2,n5)
801 z5=x(3,n5)
802 x6=x(1,n6)
803 y6=x(2,n6)
804 z6=x(3,n6)
805 x7=x(1,n7)
806 y7=x(2,n7)
807 z7=x(3,n7)
808 x8=x(1,n8)
809 y8=x(2,n8)
810 z8=x(3,n8)
811
812
813 xx1(1)=x1
814 xx2(1)=y1
815 xx3(1)=z1
816 xx1(2)=x2
817 xx2(2)=y2
818 xx3(2)=z2
819 xx1(3)=x3
820 xx2(3)=y3
821 xx3(3)=z3
822 xx1(4)=x4
823 xx2(4)=y4
824 xx3(4)=z4
825 CALL norma1(n_1,n_2,n_3,face(1),xx1,xx2,xx3)
826 IF( n4/=n3
827 . .AND.n3/=n2
828 . .AND.n2/=n1
829 . .AND.n1/=n4)THEN
830 face(1) = fourth*face(1)
831 ELSE
832 face(1) = third*face(1)
833 ENDIF
834
835 xx1(1)=x5
836 xx2(1)=y5
837 xx3(1)=z5
838 xx1(2)=x6
839 xx2(2)=y6
840 xx3(2)=z6
841 xx1(3)=x7
842 xx2(3)=y7
843 xx3(3)=z7
844 xx1(4)=x8
845 xx2(4)=y8
846 xx3(4)=z8
847 CALL norma1(n_1,n_2,n_3,face(2),xx1,xx2,xx3)
848 IF( n8/=n7
849 . .AND.n7/=n6
850 . .AND.n6/=n5
851 . .AND.n5/=n8)THEN
852 face(2) = fourth*face(2)
853 ELSE
854 face(2) = third*face(2)
855 ENDIF
856
857 xx1(1)=x2
858 xx2(1)=y2
859 xx3(1)=z2
860 xx1(2)=x3
861 xx2(2)=y3
862 xx3(2)=z3
863 xx1(3)=x7
864 xx2(3)=y7
865 xx3(3)=z7
866 xx1(4)=x6
867 xx2(4)=y6
868 xx3(4)=z6
869 CALL norma1(n_1,n_2,n_3,face(3),xx1,xx2,xx3)
870 IF( n6/=n7
871 . .AND.n7/=n3
872 . .AND.n3/=n2
873 . .AND.n2/=n6)THEN
874 face(3) = fourth*face(3)
875 ELSE
876 face(3) = third*face(3)
877 ENDIF
878
879 xx1(1)=x1
880 xx2(1)=y1
881 xx3(1)=z1
882 xx1(2)=x4
883 xx2(2)=y4
884 xx3(2)=z4
885 xx1(3)=x8
886 xx2(3)=y8
887 xx3(3)=z8
888 xx1(4)=x5
889 xx2(4)=y5
890 xx3(4)=z5
891 CALL norma1(n_1,n_2,n_3,face(4),xx1,xx2,xx3)
892 IF( n5/=n8
893 . .AND.n8/=n4
894 . .AND.n4/=n1
895 . .AND.n1/=n5)THEN
896 face(4) = fourth*face(4)
897 ELSE
898 face(4) = third*face(4)
899 ENDIF
900
901 xx1(1)=x1
902 xx2(1)=y1
903 xx3(1)=z1
904 xx1(2)=x2
905 xx2(2)=y2
906 xx3(2)=z2
907 xx1(3)=x6
908 xx2(3)=y6
909 xx3(3)=z6
910 xx1(4)=x5
911 xx2(4)=y5
912 xx3(4)=z5
913 CALL norma1(n_1,n_2,n_3,face(5),xx1,xx2,xx3)
914 IF( n5/=n6
915 . .AND.n6/=n2
916 . .AND.n2/=n1
917 . .AND.n1/=n5)THEN
918 face(5) = fourth*face(5)
919 ELSE
920 face(5) = third*face(5)
921 ENDIF
922
923 xx1(1)=x4
924 xx2(1)=y4
925 xx3(1)=z4
926 xx1(2)=x3
927 xx2(2)=y3
928 xx3(2)=z3
929 xx1(3)=x7
930 xx2(3)=y7
931 xx3(3)=z7
932 xx1(4)=x8
933 xx2(4)=y8
934 xx3(4)=z8
935 CALL norma1(n_1,n_2,n_3,face(6),xx1,xx2,xx3)
936 IF( n8/=n7
937 . .AND.n7/=n3
938 . .AND.n3/=n4
939 . .AND.n4/=n8)THEN
940
941 ELSE
942 face(6) = third*face(6)
943 ENDIF
944
945 DO k1=1,8
946 n1 = ixs(k1+1,ie)
947 IF (n1 == irect(1,iTHEN
948 DO k2=1,8
949 n2 = ixs(k2+1,ie)
950 IF (k1 == 1 .AND. k2 == 2) THEN
951 area(i) =
area(i) + (face(1)+face(5))
952 ELSEIF (k1 == 1 .AND. k2THEN
954 ELSEIF (k1 == 1 .AND. k2 == 5) THEN
955 area(i) =
area(i) + (face(4)+face(5))
956 ELSEIFTHEN
957 area(i) =
area(i) + (face(1)+face(3))
958 ELSEIF (k1 == 2 .AND. k2 == 6) THEN
959 area(i) =
area(i) + (face(5)+face(3))
960 ELSEIF (k1 == 3 .AND. k2 == 4) THEN
961 area(i) =
area(i) + (face(1)+face(6))
962 ELSEIF (k1 == 3 .AND. k2 == 7) THEN
963 area(i) =
area(i) + (face(3)+face(6))
964 ELSEIF (k1 == 4 .AND. k2 == 1) THEN
965 area(i) =
area(i) + (face(1)+face(4))
966 ELSEIF (k1 == 4 .AND. k2 == 8) THEN
968 ELSEIF (k1 == 5 .AND. k2 == 6) THEN
969 area(i) =
area(i) + (face(2)+face(5))
970 ELSEIF (k1 == 6 .AND. k2 == 7) THEN
971 area(i) =
area(i) + (face(2)+face(5))
972 ELSEIF (k1 == 7 .AND. k2 == 8) THEN
973 area(i) =
area(i) + (face(2)+face(6))
974 ELSEIF (k1 == 8.AND. k2 == 5) THEN
975 area(i) =
area(i) + (face(2)+face(4))
976 ENDIF
977 ENDDO
978 ENDIF
979 ENDDO
980 ENDIF
981 ENDDO
982 ENDDO
983 ENDDO
984 ELSEIF(nelp /= 0) THEN
985 DO i=1,nrt
987 jj1 = 0
988 DO j=knod2el1d(irect(1,i))+1,knod2el1d(irect(1,i)+1)
989 ie = nod2el1d(j)
990 jj1 = jj1+1
991 mg = ixp(5,ie)
992 dx = sqrt(geo(1,mg))
994 ielec(i) = ixp(1,ie)
995 ENDDO
996 DO j=knod2el1d(irect(2,i))+1,knod2el1d(irect(2,i)+1)
997 ie = nod2el1d(j)
998 jj1 = jj1+1
999 mg = ixp(5,ie)
1000 dx = sqrt(geo(1,mg))
1002 ENDDO
1004 ENDDO
1005 ELSEIF(nelt /= 0) THEN
1006 DO i=1,nrt
1008 jj1 = 0
1009 DO j=knod2el1d(irect(1,i))+1,knod2el1d(irect(1,i)+1)
1010 ie = nod2el1d(j)
1011 jj1 = jj1+1
1012 mg=ixt(4,ie)
1013 dx=sqrt(geo(1,mg))
1015 ielec(i) = ixt(1,ie)
1016 ENDDO
1017 DO j=knod2el1d(irect(2,i))+1,knod2el1d(irect(2,i)+1)
1018 ie = nod2el1d(j)
1019 jj1 = jj1+1
1020 mg = ixt(4,ie)
1021 dx = sqrt(geo(1,mg))
1023 ENDDO
1025 ENDDO
1026 ENDIF
1027
1028 ENDIF
1029
1030
1031 RETURN
1032
1033
1034 1300 FORMAT(2x,'COMPUTED GAP = ',1pg20.13)
1035 2001 FORMAT(2x,'Maximum distance for radiation computation = ',1pg20.13)
1036
subroutine area(d1, x, x2, y, y2, eint, stif0)
subroutine i11coq(irect, ixc, ixtg, nint, nel, neltg, is, geo, pm, thk, igeo, knod2elc, knod2eltg, nod2elc, nod2eltg, pm_stack, iworksh)
subroutine i11fil(irect, ixt, ixp, ixr, nint, nelt, nelp, nelr, nelx, is, nod2el1d, knod2el1d, kxx, ixx)
subroutine i11gmx3(x, irect, i, gapmax, xl2)
subroutine i11sol(x, irect, ixs, nint, nel, i, area, noint, knod2els, nod2els, ixs10)
subroutine friction_parts_search(ip, npartsfric, partsfric, ipl)
integer, parameter nchartitle
subroutine norma1(n1, n2, n3, area, xx1, xx2, xx3)
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)