OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sponof2.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| sponof2 ../engine/source/elements/sph/sponof2.f
25!||--- called by ------------------------------------------------------
26!|| sphprep ../engine/source/elements/sph/sphprep.F
27!||--- calls -----------------------------------------------------
28!|| initbuf ../engine/share/resol/initbuf.F
29!|| sph_nodseg ../engine/source/elements/sph/sph_nodseg.F
30!|| sphreqs ../engine/source/elements/sph/sphreq.F
31!|| spmd_sphgetimp ../engine/source/mpi/elements/spmd_sph.F
32!|| spmd_sphgetvois_off ../engine/source/mpi/elements/spmd_sph.F
33!||--- uses -----------------------------------------------------
34!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
35!|| groupdef_mod ../common_source/modules/groupdef_mod.F
36!|| initbuf_mod ../engine/share/resol/initbuf.F
37!|| sphbox ../engine/share/modules/sphbox.F
38!||====================================================================
39 SUBROUTINE sponof2(X ,V ,D ,MS ,SPBUF ,
40 2 ITAB ,KXSP ,IXSP ,NOD2SP ,NPC ,
41 3 PLD ,IPARG ,ELBUF_TAB,ISPHIO ,VSPHIO ,
42 4 PM ,IPART ,IPARTSP ,IGRSURF ,
43 5 LPRTSPH ,LONFSPH ,IWA ,MWA ,WA ,
44 6 VNORMAL ,SPHVELN ,XDP,IBUFSSG_IO ,OFF_SPH_R2R,
45 7 WFEXT)
46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE initbuf_mod
50 USE elbufdef_mod
51 USE sphbox
52 USE groupdef_mod
53C----6------------------------------------------
54C I m p l i c i t T y p e s
55C-----------------------------------------------
56#include "implicit_f.inc"
57#include "comlock.inc"
58C-----------------------------------------------
59C C o m m o n B l o c k s
60C-----------------------------------------------
61#include "vect01_c.inc"
62#include "com01_c.inc"
63#include "com04_c.inc"
64#include "com06_c.inc"
65#include "com08_c.inc"
66#include "param_c.inc"
67#include "sphcom.inc"
68#include "scr05_c.inc"
69#include "tabsiz_c.inc"
70#include "scr17_c.inc"
71#include "rad2r_c.inc"
72C-----------------------------------------------------------------
73C D u m m y A r g u m e n t s
74C-----------------------------------------------
75 INTEGER KXSP(NISP,*),IXSP(KVOISPH,*),NOD2SP(*),ITAB(*),NPC(*),
76 . IPARG(NPARG,*),ISPHIO(NISPHIO,*),IPART(LIPART1,*),
77 . LPRTSPH(2,0:NPART),LONFSPH(*),
78 . IWA(*),MWA(*),IBUFSSG_IO(SIBUFSSG_IO),IPARTSP(*),OFF_SPH_R2R(*)
79C REAL
80 my_real
81 . x(3,*) ,v(3,*) ,d(3,*) ,ms(*) ,spbuf(nspbuf,*) ,
82 . pld(*) ,vsphio(*) ,pm(npropm,*),wa(*),
83 . vnormal(3,*), sphveln(2,*)
84 double precision
85 . xdp(3,*)
86 TYPE (ELBUF_STRUCT_), TARGET, DIMENSION (NGROUP) :: ELBUF_TAB
87 TYPE (SURF_) , DIMENSION(NSURF) :: IGRSURF
88 DOUBLE PRECISION,INTENT(INOUT) :: WFEXT
89C-----------------------------------------------
90C L o c a l V a r i a b l e s
91C-----------------------------------------------
92 INTEGER I, ITYPE, IVAD,
93 . ii,ipt,jj(6),npf,
94 . ifvitsx,ifvitsy,ifvitsz,ifdens,ifpres,ifener,
95 . j,n,inod,iactive,
96 . ilaston,nlaston,i1stoff,n1stoff,
97 . isu,nseg,in1,in2,in3,in4,
98 . imat,iprop,iprt,nel,kad,ng,k,
99 . kp,np2sort,nband,
100 . knsx,kisx,knsy,kisy,knsz,kisz,
101 . knpx,kipx,knpy,kipy,knpz,kipz,kvnorm,kf,
102 . ippv,m,jnod,knod,jk,nvois1,nvois2,nvoiss1,nvoiss2,
103 . impose,iun,kk,i1,iad2,
104 . nsegb,impose2,inod0,nn,off,tflag,np2sort_old,koff,ij
105C REAL
106 my_real
107 . wfextt,volo,
108 . x1,x2,x3,x4,xp,y1,y2,y3,y4,yp,z1,z2,z3,z4,zp,
109 . xi,yi,zi,xj,yj,zj,di,dd,dmin,dmax,dps,
110 . xg,yg,zg,xx(12),xx_old(12),vmax,vi,dps_old,xi_old,
111 . yi_old,zi_old,dt_old,vmaxs
112
113 INTEGER, DIMENSION(:), ALLOCATABLE :: OFF_SPH, TAG_SPH
114 TYPE(G_BUFEL_) ,POINTER :: GBUF
115
116 LOGICAL :: lBOOL
117C-----------------------------------------------
118 my_real
119 . get_u_geo
120 EXTERNAL get_u_geo
121 DATA iun/1/
122C-----------------------------------------------
123 ALLOCATE(tag_sph(nsphr))
124 ALLOCATE(off_sph(numsph))
125 off_sph(1:numsph) = zero
126 tag_sph(1:nsphr) = zero
127
128 wfextt=zero
129C-----------------------------------------------
130 DO i=1,nsphio
131 itype=isphio(1,i)
132 vmax = em20
133 IF(itype/=1)THEN
134C outlets.
135 iprt =isphio(2,i)
136 ivad=isphio(4,i)
137 DO iactive=lprtsph(2,iprt-1)+1,lprtsph(1,iprt)
138 iwa(iactive)=0
139 IF (itype==4) THEN
140C determination of max velocity
141 n=lonfsph(iactive)
142 inod=kxsp(3,n)
143 vi=sqrt(v(1,inod)**2+v(2,inod)**2+v(3,inod)**2)
144 vmax = max(vmax,vi)
145 ENDIF
146 ENDDO
147 dt_old = vsphio(ivad+14)
148 IF (itype==4) vsphio(ivad+3)=vmax*max(dt_old,dt12)
149 ENDIF
150 ENDDO
151C-----------------------------------------------
152 DO i=1,nsphio
153 itype=isphio(1,i)
154C-------
155 iprt =isphio(2,i)
156 IF(lprtsph(1,iprt)>lprtsph(2,iprt-1))THEN
157C-------
158 IF ((itype>1).AND.(isphio(12,i)==0)) THEN
159C-------
160C general outlets & silent boundary.
161C-------
162 ivad=isphio(4,i)
163C-------
164C control up to DMAX
165 dmax=-vsphio(ivad+3)
166C-------
167C non interaction outside DI:
168 di=zero
169 DO iactive=lprtsph(2,iprt-1)+1,lprtsph(1,iprt)
170 n =lonfsph(iactive)
171 inod=kxsp(3,n)
172 di=max(di,spbuf(1,n))
173 ENDDO
174 di=two*di
175C-------
176 xbmin =ep20
177 ybmin =ep20
178 zbmin =ep20
179 xbmax=-ep20
180 ybmax=-ep20
181 zbmax=-ep20
182C------
183 isu =isphio(3,i)
184 nsegb=igrsurf(isu)%NSEG
185 nseg=isphio(10,i)
186 iad2 =isphio(11,i)
187
188 DO j=0,nseg-1
189 in1=ibufssg_io(iad2+nibsph*j )
190 in2=ibufssg_io(iad2+nibsph*j+1)
191 in3=ibufssg_io(iad2+nibsph*j+2)
192 in4=ibufssg_io(iad2+nibsph*j+3)
193
194c IN1=IBUFSSG(IAD+NISX*J )
195c IN2=IBUFSSG(IAD+NISX*J+1)
196c IN3=IBUFSSG(IAD+NISX*J+2)
197c IN4=IBUFSSG(IAD+NISX*J+3)
198
199 xp=x(1,in1)
200 yp=x(2,in1)
201 zp=x(3,in1)
202
203 xbmin=min(xbmin,xp)
204 ybmin=min(ybmin,yp)
205 zbmin=min(zbmin,zp)
206 xbmax=max(xbmax,xp)
207 ybmax=max(ybmax,yp)
208 zbmax=max(zbmax,zp)
209
210 xp=x(1,in2)
211 yp=x(2,in2)
212 zp=x(3,in2)
213 xbmin=min(xbmin,xp)
214 ybmin=min(ybmin,yp)
215 zbmin=min(zbmin,zp)
216 xbmax=max(xbmax,xp)
217 ybmax=max(ybmax,yp)
218 zbmax=max(zbmax,zp)
219
220 xp=x(1,in3)
221 yp=x(2,in3)
222 zp=x(3,in3)
223
224 xbmin=min(xbmin,xp)
225 ybmin=min(ybmin,yp)
226 zbmin=min(zbmin,zp)
227 xbmax=max(xbmax,xp)
228 ybmax=max(ybmax,yp)
229 zbmax=max(zbmax,zp)
230
231 IF(in4/=in3)THEN
232 xp=x(1,in4)
233 yp=x(2,in4)
234 zp=x(3,in4)
235 xbmin=min(xbmin,xp)
236 ybmin=min(ybmin,yp)
237 zbmin=min(zbmin,zp)
238 xbmax=max(xbmax,xp)
239 ybmax=max(ybmax,yp)
240 zbmax=max(zbmax,zp)
241 ENDIF
242 ENDDO
243
244 xbmin=xbmin+dmax
245 ybmin=ybmin+dmax
246 zbmin=zbmin+dmax
247 xbmax=xbmax-dmax
248 ybmax=ybmax-dmax
249 zbmax=zbmax-dmax
250C------
251 np2sort=0
252 DO iactive=lprtsph(2,iprt-1)+1,lprtsph(1,iprt)
253 n=lonfsph(iactive)
254 inod=kxsp(3,n)
255 xi=x(1,inod)
256 yi=x(2,inod)
257 zi=x(3,inod)
258 IF( xi>xbmin.AND.yi>ybmin.AND.zi>zbmin
259 . .AND.xi<xbmax.AND.yi<ybmax.AND.zi<zbmax)THEN
260 np2sort=np2sort+1
261 mwa(np2sort)=iactive
262 ENDIF
263 ENDDO
264C------
265C bucket.
266c tri bucket
267 dbucs=-dmax
268 DO j=1,nseg
269 kk =iad2+nibsph*(j-1)
270c KKB =IAD+NISX*(J-1)
271 in1=ibufssg_io(kk )
272 in2=ibufssg_io(kk+1)
273 in3=ibufssg_io(kk+2)
274 in4=ibufssg_io(kk+3)
275c IN1=IBUFSSG(KKB )
276c IN2=IBUFSSG(KKB+1)
277c IN3=IBUFSSG(KKB+2)
278c IN4=IBUFSSG(KKB+3)
279 x1=x(1,in1)
280 y1=x(2,in1)
281 z1=x(3,in1)
282 x2=x(1,in2)
283 y2=x(2,in2)
284 z2=x(3,in2)
285 x3=x(1,in3)
286 y3=x(2,in3)
287 z3=x(3,in3)
288
289 dbucs=max(dbucs,abs(x1-x2))
290 dbucs=max(dbucs,abs(y1-y2))
291 dbucs=max(dbucs,abs(z1-z2))
292 dbucs=max(dbucs,abs(x2-x3))
293 dbucs=max(dbucs,abs(y2-y3))
294 dbucs=max(dbucs,abs(z2-z3))
295 dbucs=max(dbucs,abs(x3-x1))
296 dbucs=max(dbucs,abs(y3-y1))
297 dbucs=max(dbucs,abs(z3-z1))
298 IF(in4/=in3)THEN
299 x4=x(1,in4)
300 y4=x(2,in4)
301 z4=x(3,in4)
302 ELSE
303 x4=x3
304 y4=y3
305 z4=z3
306 ENDIF
307 dbucs=max(dbucs,abs(x1-x4))
308 dbucs=max(dbucs,abs(y1-y4))
309 dbucs=max(dbucs,abs(z1-z4))
310 dbucs=max(dbucs,abs(x2-x4))
311 dbucs=max(dbucs,abs(y2-y4))
312 dbucs=max(dbucs,abs(z2-z4))
313 dbucs=max(dbucs,abs(x3-x4))
314 dbucs=max(dbucs,abs(y3-y4))
315 dbucs=max(dbucs,abs(z3-z4))
316 ENDDO
317C if more than 1 box, each box size >= DBUCS
318 nbox =max(iun,int((xbmax-xbmin)/dbucs))
319 nboy =max(iun,int((ybmax-ybmin)/dbucs))
320 nboz =max(iun,int((zbmax-zbmin)/dbucs))
321 nband=max(nbox,nboy,nboz)+1
322 knsx=numsph+1
323 kisx=knsx+nband+1
324 knsy=kisx+4*nseg
325 kisy=knsy+nband+1
326 knsz=kisy+4*nseg
327 kisz=knsz+nband+1
328 knpx=kisz+4*nseg
329 kipx=knpx+nband+1
330 knpy=kipx+3*np2sort
331 kipy=knpy+nband+1
332 knpz=kipy+3*np2sort
333 kipz=knpz+nband+1
334 kvnorm=kipz+3*np2sort
335C KF =KVNORM+3*NUMSPH
336
337C-----------------------------------------------
338C Build buckets from segments & particles.
339C-----------------------------------------------
340 CALL sphreqs(nseg ,ibufssg_io(iad2) ,x ,np2sort ,mwa ,
341 2 lonfsph ,kxsp ,wa(knsx) ,wa(kisx) ,wa(knsy) ,
342 3 wa(kisy) ,wa(knsz) ,wa(kisz) ,wa(knpx) ,wa(kipx) ,
343 4 wa(knpy) ,wa(kipy) ,wa(knpz) ,wa(kipz) ,wa ,
344 5 wa(kvnorm),vsphio(ivad+14),vsphio(ivad+13),v,spbuf ,
345 6 isphio(1,i))
346C
347 vsphio(ivad+14) = dt12
348C------
349 ELSEIF (isphio(12,i)>0) THEN
350C------
351C control section / new outlet - silent boundary
352C------
353 tflag = 0
354 ivad = isphio(4,i)
355 dmax=-vsphio(ivad+3)
356C-------
357 xbmin =ep20
358 ybmin =ep20
359 zbmin =ep20
360 xbmax=-ep20
361 ybmax=-ep20
362 zbmax=-ep20
363 vmaxs=-ep20
364C-------
365 dt_old = vsphio(ivad+14)
366 vsphio(ivad+14) = dt12
367C-------
368 DO ii=1,4
369 IF (ii<4) THEN
370 IF (isphio(12,i)==1) THEN
371 off = 3*ii
372 xi=vsphio(ivad+off+1)
373 yi=vsphio(ivad+off+2)
374 zi=vsphio(ivad+off+3)
375 vi=zero
376 ELSE
377 xi=x(1,isphio(12+ii,i))
378 yi=x(2,isphio(12+ii,i))
379 zi=x(3,isphio(12+ii,i))
380 vi=sqrt(v(1,isphio(12+ii,i))**2+v(2,isphio(12+ii,i))**2+v(3,isphio(12+ii,i))**2)
381 ENDIF
382 xx(3*(ii-1)+1) = xi
383 xx(3*(ii-1)+2) = yi
384 xx(3*(ii-1)+3) = zi
385 ELSE
386C calcul nod4
387 xi = xx(4)+xx(7)-xx(1)
388 yi = xx(5)+xx(8)-xx(2)
389 zi = xx(6)+xx(9)-xx(3)
390C permutation nod3 nod4
391 xx(10) = xx(7)
392 xx(11) = xx(8)
393 xx(12) = xx(9)
394 xx(7) = xi
395 xx(8) = yi
396 xx(9) = zi
397 ENDIF
398 xbmin=min(xbmin,xi)
399 ybmin=min(ybmin,yi)
400 zbmin=min(zbmin,zi)
401 xbmax=max(xbmax,xi)
402 ybmax=max(ybmax,yi)
403 zbmax=max(zbmax,zi)
404 vmaxs=max(vmaxs,vi)
405 END DO
406C-------
407C computation of old coordinates of the nodes of the segment
408 IF (isphio(12,i)==1) THEN
409 DO ii=1,12
410 xx_old(ii) = xx(ii)
411 END DO
412 ELSE
413 DO ii=1,3
414 xx_old(3*(ii-1)+1) = x(1,isphio(12+ii,i))-dt_old*v(1,isphio(12+ii,i))
415 xx_old(3*(ii-1)+2) = x(2,isphio(12+ii,i))-dt_old*v(2,isphio(12+ii,i))
416 xx_old(3*(ii-1)+3) = x(3,isphio(12+ii,i))-dt_old*v(3,isphio(12+ii,i))
417 END DO
418 xx_old(10) = xx_old(7)
419 xx_old(11) = xx_old(8)
420 xx_old(12) = xx_old(9)
421 xx_old(7) = xx_old(4)+xx_old(7)-xx_old(1)
422 xx_old(8) = xx_old(5)+xx_old(8)-xx_old(2)
423 xx_old(9) = xx_old(6)+xx_old(9)-xx_old(3)
424 ENDIF
425C-------
426 IF (isphio(1,i)==4) dmax = dmax + vmaxs*max(dt_old,dt12)
427 xbmin=xbmin+dmax
428 ybmin=ybmin+dmax
429 zbmin=zbmin+dmax
430 xbmax=xbmax-dmax
431 ybmax=ybmax-dmax
432 zbmax=zbmax-dmax
433C------
434 kvnorm=numsph+1
435C------
436 np2sort=0
437 DO iactive=lprtsph(2,iprt-1)+1,lprtsph(1,iprt)
438 n=lonfsph(iactive)
439 inod=kxsp(3,n)
440 xi=x(1,inod)
441 yi=x(2,inod)
442 zi=x(3,inod)
443 IF( xi>xbmin.AND.yi>ybmin.AND.zi>zbmin
444 . .AND.xi<xbmax.AND.yi<ybmax.AND.zi<zbmax)THEN
445 np2sort=np2sort+1
446 mwa(np2sort)=iactive
447 xi_old = xi - v(1,inod)*dt_old
448 yi_old = yi - v(2,inod)*dt_old
449 zi_old = zi - v(3,inod)*dt_old
450 CALL sph_nodseg(xi_old,yi_old,zi_old,xx_old,tflag,np2sort,lonfsph,mwa,wa,wa(kvnorm),1)
451 dps_old=wa(np2sort)
452 CALL sph_nodseg(xi,yi,zi,xx,tflag,np2sort,lonfsph,mwa,wa,wa(kvnorm),1)
453 dps=wa(np2sort)
454 IF (dps _old*dps<zero) THEN
455C pur outlet/silent boundaries on permute le sens de comptage
456 IF (isphio(1,i)/=4) dps_old = -dps_old
457C la particule a traverse la surface dans le sens +
458 IF (dps_old>zero) vsphio(ivad+13) = vsphio(ivad+13)-spbuf(12,n)
459C la particule a traverse la surface dans le sens -
460 IF (dps_old<zero) vsphio(ivad+13) = vsphio(ivad+13)+spbuf(12,n)
461 ENDIF
462 ENDIF
463 ENDDO
464C------
465 ENDIF
466C------
467 IF ((itype==2).OR.(itype==3)) THEN
468C------
469 DO kp=1,np2sort
470 dps=wa(kp)
471 iactive=mwa(kp)
472 IF(dps>dmax.AND.dps<zero)THEN
473C particule imposee :
474C conditions multiples => la derniere condition rencontree impose ..
475 n=lonfsph(iactive)
476 iwa(iactive)=i
477 ng =mod(kxsp(2,n),ngroup+1)
478 kxsp(2,n)=ng+i*(ngroup+1)
479 IF(tt==zero)THEN
480 inod=kxsp(3,n)
481 kk=kvnorm+3*(n-1)
482 sphveln(1,n)= wa(kk )*v(1,inod)
483 . +wa(kk+1)*v(2,inod)
484 . +wa(kk+2)*v(3,inod)
485 ng=mod(kxsp(2,n),ngroup+1)
486 CALL initbuf(iparg ,ng ,
487 2 mtn ,nel ,nft ,kad ,ity ,
488 3 npt ,jale ,ismstr ,jeul ,jtur ,
489 4 jthe ,jlag ,jmult ,jhbe ,jivf ,
490 5 nvaux ,jpor ,jcvt ,jclose ,jplasol ,
491 6 irep ,iint ,igtyp ,israt ,isrot ,
492 7 icsen ,isorth ,isorthg ,ifailure,jsms)
493!
494 DO ij=1,6
495 jj(ij) = nel*(ij-1)
496 ENDDO
497!
498 gbuf => elbuf_tab(ng)%GBUF
499 k=n-nft
500!
501 sphveln(2,n)=-third*
502 . (gbuf%SIG(jj(1)+k)+gbuf%SIG(jj(2)+k)+gbuf%SIG(jj(3)+k))
503 ENDIF
504 kk=kvnorm+3*(n-1)
505 vnormal(1,n)=wa(kk )
506 vnormal(2,n)=wa(kk+1)
507 vnormal(3,n)=wa(kk+2)
508 ELSEIF(dps>=zero.AND.iwa(iactive)==0)THEN
509 n=lonfsph(iactive)
510 inod=kxsp(3,n)
511 kk=kvnorm+3*(n-1)
512 sphveln(1,n)= wa(kk )*v(1,inod)
513 . +wa(kk+1)*v(2,inod)
514 . +wa(kk+2)*v(3,inod)
515 ng=mod(kxsp(2,n),ngroup+1)
516 CALL initbuf(iparg ,ng ,
517 2 mtn ,nel ,nft ,kad ,ity ,
518 3 npt ,jale ,ismstr ,jeul ,jtur ,
519 4 jthe ,jlag ,jmult ,jhbe ,jivf ,
520 5 nvaux ,jpor ,jcvt ,jclose ,jplasol ,
521 6 irep ,iint ,igtyp ,israt ,isrot ,
522 7 icsen ,isorth ,isorthg ,ifailure,jsms)
523!
524 DO ij=1,6
525 jj(ij) = nel*(ij-1)
526 ENDDO
527!
528 gbuf => elbuf_tab(ng)%GBUF
529 k=n-nft
530!
531 sphveln(2,n)=-third*
532 . (gbuf%SIG(jj(1)+k)+gbuf%SIG(jj(2)+k)+gbuf%SIG(jj(3)+k))
533 ENDIF
534 ENDDO
535 ENDIF
536 ENDIF
537 ENDDO !fin DO I=1,NSPHIO
538C-----------------------------------------------
539 IF(nspmd>1)THEN
540C send IMPOSE value
541 CALL spmd_sphgetimp(kxsp)
542 END IF
543C requete ppv <=> non zero interaction / impose
544c ppv: plus proche voisin
545C-----------------------------------------------
546 DO i=1,nsphio
547 itype=isphio(1,i)
548 IF ((itype==2).OR.(itype==3)) THEN
549C-------
550C general outlets & silent boundary.
551 iprt =isphio(2,i)
552 DO iactive=lprtsph(2,iprt-1)+1,lprtsph(1,iprt)
553 IF(iwa(iactive)/=0)THEN
554 n =lonfsph(iactive)
555 inod=kxsp(3,n)
556 xi=x(1,inod)
557 yi=x(2,inod)
558 zi=x(3,inod)
559 ippv=0
560 dmin=ep20
561
562 DO j=1,kxsp(4,n)
563 jnod=ixsp(j,n)
564
565 IF(jnod>0)THEN ! particule locale
566 m=nod2sp(jnod)
567 IF(kxsp(2,m)>=0)THEN
568 impose=kxsp(2,m)/(ngroup+1)
569 lbool=.false.
570 IF(impose == 0)THEN
571 lbool=.true.
572 ELSE
573 IF(isphio(1,impose)==1)lbool=.true.
574 ENDIF
575 IF(lbool)THEN
576 xj =x(1,jnod)
577 yj =x(2,jnod)
578 zj =x(3,jnod)
579 dd =(xi-xj)*(xi-xj)+(yi-yj)*(yi-yj)+(zi-zj)*(zi-zj)
580 IF(dd<dmin)THEN
581 ippv=jnod
582 dmin=dd
583 ENDIF
584 ENDIF
585 ENDIF
586 ELSE ! particule remote
587 nn = -jnod
588 impose = nint(xsphr(12,nn))
589 IF(impose>0)THEN
590 impose2=isphio(1,impose)
591 ELSE
592 impose2=0
593 ENDIF
594 IF(impose2==0.OR.impose2==1)THEN
595 xj = xsphr(3,nn)
596 yj = xsphr(4,nn)
597 zj = xsphr(5,nn)
598 dd =(xi-xj)*(xi-xj)+(yi-yj)*(yi-yj)+(zi-zj)*(zi-zj)
599 IF(dd<dmin)THEN
600 ippv=jnod
601 dmin=dd
602 ENDIF
603 ENDIF
604 ENDIF
605 ENDDO
606C-------
607 IF(ippv==0)THEN
608C no interaction between the outgoing particle
609C and the particles inside the domain.
610 ng=mod(abs(kxsp(2,n)),ngroup+1)
611 kxsp(2,n)=-ng
612c elmt delete OFF_SPH set to 1
613 off_sph(n)=1
614 IF (r2r_siu/=0) off_sph_r2r(inod) = -1
615 ENDIF
616 ENDIF
617 ENDDO
618 END IF
619 ENDDO
620C-----------------------------------------------
621C supprime particules n'inter-agissant plus avec l'interieur du domaine.
622C-----------------------------------------------
623 DO i1=1,nsphio
624 itype=isphio(1,i1)
625 IF ((itype==2).OR.(itype==3)) THEN
626C-------
627C general outlets & silent boundary.
628 iprt =isphio(2,i1)
629 iactive=lprtsph(2,iprt-1)+1
630 DO WHILE(iactive<=lprtsph(1,iprt))
631 IF(iwa(iactive)/=0)THEN
632 n =lonfsph(iactive)
633 IF(kxsp(2,n)<0)THEN
634 inod=kxsp(3,n)
635 IF(tt/=zero)wfextt=wfextt-half*ms(inod)
636 . *(v(1,inod)*v(1,inod)+v(2,inod)*v(2,inod)+v(3,inod)*v(3,inod))
637 v(1,inod)=zero
638 v(2,inod)=zero
639 v(3,inod)=zero
640 i =iwa(iactive)
641 IF (isphio(12,i)==0) THEN
642 isu=isphio(3,i)
643 iad2=isphio(11,i)
644c IN1=IBUFSSG_IO(IAD2 )
645c IN2=IBUFSSG_IO(IAD2+1)
646c IN1B=IBUFSSG(IAD )
647c IN2B=IBUFSSG(IAD+1)
648 x1 =x(1,in1)
649 y1 =x(2,in1)
650 z1 =x(3,in1)
651 x2 =x(1,in2)
652 y2 =x(2,in2)
653 z2 =x(3,in2)
654 ELSEIF (isphio(12,i)==1) THEN
655 ivad = isphio(4,i)
656 x1 =vsphio(ivad+4)
657 y1 =vsphio(ivad+5)
658 z1 =vsphio(ivad+6)
659 x2 =vsphio(ivad+7)
660 y2 =vsphio(ivad+8)
661 z2 =vsphio(ivad+9)
662 ELSE
663 x1 =x(1,isphio(13,i))
664 y1 =x(2,isphio(14,i))
665 z1 =x(3,isphio(15,i))
666 x2 =x(1,isphio(16,i))
667 y2 =x(2,isphio(17,i))
668 z2 =x(3,isphio(18,i))
669 ENDIF
670 xg =half*(x1+x2)
671 yg =half*(y1+y2)
672 zg =half*(z1+z2)
673 IF(itab(inod) >= 1000000000)THEN
674 d(1,inod)=zero
675 d(2,inod)=zero
676 d(3,inod)=zero
677 ELSE
678 d(1,inod)=xg-x(1,inod)+d(1,inod)
679 d(2,inod)=yg-x(2,inod)+d(2,inod)
680 d(3,inod)=zg-x(3,inod)+d(3,inod)
681 ENDIF
682 x(1,inod)=xg
683 x(2,inod)=yg
684 x(3,inod)=zg
685 IF (iresp==1) THEN
686 xdp(1,inod)=xg
687 xdp(2,inod)=yg
688 xdp(3,inod)=zg
689 END IF
690 ng=mod(abs(kxsp(2,n)),ngroup+1)
691 CALL initbuf(iparg ,ng ,
692 2 mtn ,nel ,nft ,kad ,ity ,
693 3 npt ,jale ,ismstr ,jeul ,jtur ,
694 4 jthe ,jlag ,jmult ,jhbe ,jivf ,
695 5 nvaux ,jpor ,jcvt ,jclose ,jplasol ,
696 6 irep ,iint ,igtyp ,israt ,isrot ,
697 7 icsen ,isorth ,isorthg ,ifailure,jsms)
698!
699 DO ij=1,6
700 jj(ij) = nel*(ij-1)
701 ENDDO
702!
703 gbuf => elbuf_tab(ng)%GBUF
704 k=n-nft
705 volo=gbuf%VOL(k)
706 wfextt=wfextt-volo*gbuf%EINT(k)
707 gbuf%EINT(k)=zero
708!
709 gbuf%SIG(jj(1)+k)=zero
710 gbuf%SIG(jj(2)+k)=zero
711 gbuf%SIG(jj(3)+k)=zero
712 gbuf%SIG(jj(4)+k)=zero
713 gbuf%SIG(jj(5)+k)=zero
714 gbuf%SIG(jj(6)+k)=zero
715 gbuf%OFF(k)=zero
716C
717 kxsp(2,n)=-ng
718 iparg(10,ng)=iparg(10,ng)-1
719 IF(iparg(10,ng)==0)iparg(8,ng)=1
720C ISPHBUC=1
721C chainage:
722 ilaston =lprtsph(1,iprt)
723 nlaston =lonfsph(ilaston)
724 lonfsph(iactive) =nlaston
725 iwa(iactive) =iwa(ilaston)
726 lonfsph(ilaston) =n
727 ilaston =ilaston-1
728 lprtsph(1,iprt) =ilaston
729
730C mise a jour des voisinages:
731c 1) ON NE TRAITE QUE LES VOISINS LOCAUX
732
733 DO j=1,kxsp(5,n)
734
735 inod0=ixsp(j,n)
736
737 IF(inod0<0)THEN
738c will be treated in SPMD_SPHGETVOIS_OFF
739c TAG_SPH set to 1
740 tag_sph(-inod0)=1
741 ELSE
742 m =nod2sp(inod0)
743 nvois1=0
744 DO k=1,kxsp(4,m)
745 knod=ixsp(k,m)
746 IF(knod/=inod)THEN
747 nvois1=nvois1+1
748 ixsp(nvois1,m)=knod
749 ENDIF
750 ENDDO
751
752 nvois2=nvois1
753 DO k=kxsp(4,m)+1,kxsp(5,m)
754 knod=ixsp(k,m)
755 IF(knod/=inod)THEN
756 nvois2=nvois2+1
757 ixsp(nvois2,m)=knod
758 ENDIF
759 ENDDO
760
761 nvoiss1=0
762 DO k=kxsp(5,m)+1,kxsp(5,m)+kxsp(6,m)
763 jk =ixsp(k,m)
764 IF(jk<0)THEN
765 nvoiss1=nvoiss1+1
766 ixsp(nvois2+nvoiss1,m)=jk
767 ELSE
768 knod=kxsp(3,jk/(nspcond+1))
769 IF(knod/=inod)THEN
770 nvoiss1=nvoiss1+1
771 ixsp(nvois2+nvoiss1,m)=jk
772 ENDIF
773 ENDIF
774 ENDDO
775
776 nvoiss2=nvoiss1
777 DO k=kxsp(5,m)+kxsp(6,m)+1,kxsp(5,m)+kxsp(7,m)
778 jk =ixsp(k,m)
779 IF(jk<0)THEN
780 nvoiss2=nvoiss2+1
781 ixsp(nvois2+nvoiss2,m)=jk
782 ELSE
783 knod=kxsp(3,jk/(nspcond+1))
784 IF(knod/=inod)THEN
785 nvoiss2=nvoiss2+1
786 ixsp(nvois2+nvoiss2,m)=jk
787 ENDIF
788 ENDIF
789 ENDDO
790 kxsp(4,m)= nvois1
791 kxsp(5,m)= nvois2
792 kxsp(6,m)=nvoiss1
793 kxsp(7,m)=nvoiss2
794 ENDIF ! fin IF(INOD0<0)THEN
795
796 ENDDO !fin do J=1,KXSP(5,N)
797
798 ELSE
799 iactive=iactive+1
800 ENDIF
801 ELSE
802 iactive=iactive+1
803 ENDIF
804 ENDDO
805 END IF
806 ENDDO
807
808c 2) treatment of remote neighbours only
809 IF(nspmd>1)
810 . CALL spmd_sphgetvois_off(off_sph, tag_sph,
811 . kxsp, ixsp)
812C-------------------------------------------
813!$OMP ATOMIC
814 wfext=wfext+wfextt
815C-------------------------------------------
816 DEALLOCATE(off_sph,tag_sph)
817
818 RETURN
819 END
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine initbuf(iparg, ng, mtn, llt, nft, iad, ity, npt, jale, ismstr, jeul, jtur, jthe, jlag, jmult, jhbe, jivf, mid, jpor, jcvt, jclose, jpla, irep, iint, igtyp, israt, isrot, icsen, isorth, isorthg, ifailure, jsms)
Definition initbuf.F:261
integer nsphr
Definition sphbox.F:83
subroutine sph_nodseg(xi, yi, zi, xx, tflag, np, lonfsph, ixp, dps, wnormal, flag)
Definition sph_nodseg.F:30
subroutine sphreqs(nseg, iseg, x, ncel, ixp, lonfsph, kxsp, nsx, isx, nsy, isy, nsz, isz, npx, ipx, npy, ipy, npz, ipz, dps, wnormal, dt_old, mass_cum, v, spbuf, itype)
Definition sphreq.F:40
subroutine spmd_sphgetimp(kxsp)
Definition spmd_sph.F:1832
subroutine spmd_sphgetvois_off(off_sph, tag_sph, kxsp, ixsp)
Definition spmd_sph.F:1627
subroutine sponof2(x, v, d, ms, spbuf, itab, kxsp, ixsp, nod2sp, npc, pld, iparg, elbuf_tab, isphio, vsphio, pm, ipart, ipartsp, igrsurf, lprtsph, lonfsph, iwa, mwa, wa, vnormal, sphveln, xdp, ibufssg_io, off_sph_r2r, wfext)
Definition sponof2.F:46