38 . X ,SPBUF ,KXSP ,IXSP ,IPARTSP ,
39 . IPARG ,NGROUNC ,IGROUNC ,ELBUF_TAB,ITASK ,
40 . NOD2SP ,SOL2SPH ,SPH2SOL ,IXS ,MS ,
41 . PM ,IADS ,ADDCNE ,FSKYD ,DMSPH ,
49 use element_mod ,
only : nixs
53#include "implicit_f.inc"
67#include "vect01_c.inc"
72 . IPARTSP(*), IPARG(NPARG,*), NGROUNC,
73 . IGROUNC(*), ITASK, IXSP(KVOISPH,*), NOD2SP(*),
74 . SOL2SPH(2,*), SPH2SOL(*), IXS(NIXS,*),
75 . IADS(8,*), ADDCNE(*), ICONTACT(*)
77 . x(3,*), spbuf(nspbuf,*), ms(*), pm(npropm,*), fskyd(*),
79 TYPE (elbuf_struct_),
TARGET,
DIMENSION(NGROUP) :: elbuf_tab
83 INTEGER I, N, KP, NG, MG, NP, KFT, IG, NELEM,
84 . nel, offset, inod, imat,
85 . n1, n2, n3, n4, n5, n6, n7, n8,
86 . k1, k2, k3, k4, k5, k6, k7, k8,
89 . dm, rho0, ehourt, ek, vi2, vxi, vyi, vzi,
90 . vx1, vx2, vx3, vx4, vx5, vx6, vx7, vx8,
91 . vy1, vy2, vy3, vy4, vy5, vy6, vy7, vy8,
92 . vz1, vz2, vz3, vz4, vz5, vz6, vz7, vz8
95 TYPE(g_bufel_) ,
POINTER :: GBUF, GBUFSP
96 TYPE(L_BUFEL_) ,
POINTER :: LBUF
97 TYPE(BUF_MAT_) ,
POINTER :: MBUF
102 IF(iparg(8,ng)==1)
GOTO 50
104 DO nelem = 1,iparg(2,ng),nvsiz
107 nft =iparg(3,ng) + offset
110 ipartsph=iparg(69,ng)
112 llt=
min(nvsiz,nel-nelem+1)
113 IF(ity==1.AND.ipartsph/=0)
THEN
115 gbuf => elbuf_tab(ng)%GBUF
116 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(1,1,1)
117 mbuf => elbuf_tab(ng)%BUFLY(1)%MAT(1,1,1)
120 IF(gbuf%OFF(i)/=zero)
THEN
122 DO kp=1,sol2sph(2,n)-sol2sph(1,n)
125 IF(icontact(inod)/=0)
THEN
128 gbuf%OFF(i)=four_over_5
132 .
' -- PARTICLE INTO CONTACT => DELETE SOLID ELEMENT AT NEXT CYCLE',
135 .
' -- PARTICLE INTO CONTACT => DELETE SOLID ELEMENT AT NEXT CYCLE',
137#include "lockoff.inc"
159 IF(iparg(8,ng)==1)
GOTO 100
161 DO nelem = 1,iparg(2,ng),nvsiz
164 nft =iparg(3,ng) + offset
167 ipartsph=iparg(69,ng)
169 llt=
min(nvsiz,nel-nelem+1)
170 IF(ity==1.AND.ipartsph/=0)
THEN
172 gbuf => elbuf_tab(ng)%GBUF
173 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(1,1,1)
174 mbuf => elbuf_tab(ng)%BUFLY(1)%MAT(1,1,1)
177 IF(gbuf%OFF(i)==zero)
THEN
186 DO kp=1,sol2sph(2,n)-sol2sph(1,n)
188 mg =mod(-kxsp(2,np),ngroup+1)
190 gbufsp => elbuf_tab(mg)%GBUF
192 gbufsp%OFF(np-kft)=one
196 vi2= v(1,inod)*v(1,inod)
197 . +v(2,inod)*v(2,inod)
198 . +v(3,inod)*v(3,inod)
199 ek=ek+half*ms(inod)*vi2
211 dm=one_over_8*gbuf%VOL(i)*rho0
213 dmsph(n1)=dmsph(n1)+dm
214 dmsph(n2)=dmsph(n2)+dm
215 dmsph(n3)=dmsph(n3)+dm
216 dmsph(n4)=dmsph(n4)+dm
217 dmsph(n5)=dmsph(n5)+dm
218 dmsph(n6)=dmsph(n6)+dm
219 dmsph(n7)=dmsph(n7)+dm
220 dmsph(n8)=dmsph(n8)+dm
254 vxi=vx1+vx2+vx3+vx4+vx5+vx6+vx7+vx8
255 vyi=vy1+vy2+vy3+vy4+vy5+vy6+vy7+vy8
256 vzi=vz1+vz2+vz3+vz4+vz5+vz6+vz7+vz8
257 vi2=vx1*vx1+vx2*vx2+vx3*vx3+vx4*vx4
258 1 +vx5*vx5+vx6*vx6+vx7*vx7+vx8*vx8
259 2 +vy1*vy1+vy2*vy2+vy3*vy3+vy4*vy4
260 3 +vy5*vy5+vy6*vy6+vy7*vy7+vy8*vy8
261 4 +vz1*vz1+vz2*vz2+vz3*vz3+vz4*vz4
262 5 +vz5*vz5+vz6*vz6+vz7*vz7+vz8*vz8
265 ehourt=ehourt+half*dm*vi2-ek
280 nodft = 1+itask*numnod/ nthread
281 nodlt = (itask+1)*numnod/nthread
283 fskyd(addcne(n):addcne(n+1)-1)=zero
291 IF(iparg(8,ng)==1)
GOTO 200
293 DO nelem = 1,iparg(2,ng),nvsiz
296 nft =iparg(3,ng) + offset
299 ipartsph=iparg(69,ng)
301 llt=
min(nvsiz,nel-nelem+1)
302 IF(ity==1.AND.ipartsph/=0)
THEN
304 gbuf => elbuf_tab(ng)%GBUF
305 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(1,1,1)
306 mbuf => elbuf_tab(ng)%BUFLY(1)%MAT(1,1,1)
309 IF(gbuf%OFF(i)==zero)
THEN
318 DO kp=1,sol2sph(2,n)-sol2sph(1,n)
320 mg =mod(-kxsp(2,np),ngroup+1)
322 gbufsp => elbuf_tab(mg)%GBUF
323 kxsp(2,np) =abs(kxsp(2,np))
324 gbufsp%OFF(np-kft)=one
328 vi2= v(1,inod)*v(1,inod)
329 . +v(2,inod)*v(2,inod)
330 . +v(3,inod)*v(3,inod)
331 ek=ek+half*ms(inod)*vi2
335 dm=one_over_8*gbuf%VOL(i)*rho0
386 vxi=vx1+vx2+vx3+vx4+vx5+vx6+vx7+vx8
387 vyi=vy1+vy2+vy3+vy4+vy5+vy6+vy7+vy8
388 vzi=vz1+vz2+vz3+vz4+vz5+vz6+vz7+vz8
389 vi2=vx1*vx1+vx2*vx2+vx3*vx3+vx4*vx4
390 1 +vx5*vx5+vx6*vx6+vx7*vx7+vx8*vx8
391 2 +vy1*vy1+vy2*vy2+vy3*vy3+vy4*vy4
392 3 +vy5*vy5+vy6*vy6+vy7*vy7+vy8*vy8
393 4 +vz1*vz1+vz2*vz2+vz3*vz3+vz4*vz4
394 5 +vz5*vz5+vz6*vz6+vz7*vz7+vz8*vz8
397 ehourt=ehourt+half*dm*vi2-ek
413#include "lockoff.inc"
431 . X ,SPBUF ,KXSP ,IXSP ,IPARTSP ,
432 . IPARG ,NGROUNC ,IGROUNC ,ELBUF_TAB,ITASK ,
433 . NOD2SP ,SOL2SPH ,SPH2SOL ,IXS ,MS ,
434 . PM ,IADS ,ADDCNE ,FSKYD ,DMSPH ,
442 use element_mod ,
only : nixs
446#include "implicit_f.inc"
447#include "comlock.inc"
451#include "com01_c.inc"
452#include "com04_c.inc"
453#include "param_c.inc"
454#include "scr17_c.inc"
457#include "units_c.inc"
458#include "vect01_c.inc"
463 . IPARTSP(*), IPARG(NPARG,*), NGROUNC,
464 . IGROUNC(*), ITASK, IXSP(KVOISPH,*), NOD2SP(*),
465 . SOL2SPH(2,*), SPH2SOL(*), IXS(NIXS,*),
466 . IADS(8,*), ADDCNE(*), ICONTACT(*), IPART(LIPART1,*)
468 . X(3,*), SPBUF(NSPBUF,*), MS(*), PM(NPROPM,*), FSKYD(*),
470 TYPE (ELBUF_STRUCT_),
TARGET,
DIMENSION(NGROUP) :: ELBUF_TAB
474 INTEGER I, N, KP, NG, J, NP, IG, NELEM,
475 . NEL, OFFSET, NVOIS, M, INOD, JNOD, NN
478 TYPE(g_bufel_) ,
POINTER :: GBUF
479 TYPE(L_BUFEL_) ,
POINTER :: LBUF
480 TYPE(buf_mat_) ,
POINTER :: MBUF
485 IF(iparg(8,ng)==1)
GOTO 50
487 DO nelem = 1,iparg(2,ng),nvsiz
490 nft =iparg(3,ng) + offset
493 ipartsph=iparg(69,ng)
495 llt=
min(nvsiz,nel-nelem+1)
496 IF(ity==1.AND.ipartsph/=0)
THEN
498 gbuf => elbuf_tab(ng)%GBUF
499 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(1,1,1)
500 mbuf => elbuf_tab(ng)%BUFLY(1)%MAT(1,1,1)
502 IF ((itsol2sph==1).OR.(nsubs==0))
THEN
506 IF(gbuf%OFF(i)/=zero)
THEN
508 DO kp=1,sol2sph(2,n)-sol2sph(1,n)
516 IF(ipartsp(m)/=ipartsp(np))
THEN
519 gbuf%OFF(i)=four_over_5
522 WRITE(iout,5000) ixs(nixs,n)
523 WRITE(istdo,5000) ixs(nixs,n)
524#include "lockoff.inc"
529 IF(nint(xsphr(14,nn))/=ipartsp(np))
THEN
532 gbuf%OFF(i)=four_over_5
535 WRITE(iout,5000) ixs(nixs,n)
536 WRITE(istdo,5000) ixs(nixs,n)
537#include "lockoff.inc"
547 ELSEIF (itsol2sph==2)
THEN
551 IF(gbuf%OFF(i)/=zero)
THEN
553 DO kp=1,sol2sph(2,n)-sol2sph(1,n)
561 IF((ipart(3,ipartsp(m))/=ipart(3,ipartsp(np))).OR.
562 . (((ipart(3,ipartsp(m))+ipart(3,ipartsp(np)))==2*nsubs).
563 . and.(ipartsp(m)/=ipartsp(np))))
THEN
566 gbuf%OFF(i)=four_over_5
569 WRITE(iout,6000) ixs(nixs,n)
570 WRITE(istdo,6000) ixs(nixs,n)
571#include
"lockoff.inc"
576 IF((ipart(3,nint(xsphr(14,nn)))/=ipart(3,ipartsp(np))).OR.
577 . (((ipart(3,ipartsp(np))+ipart(3,nint(xsphr(14,nn)))==2*nsubs).
578 . and.(nint(xsphr(14,nn))/=ipartsp(np)))))
THEN
584 WRITE(iout,6000) ixs(nixs,n)
585 WRITE(istdo,6000) ixs(nixs,n)
586#include "lockoff.inc"
607 &
' -- PARTICLE INTERACTING W/OTHER SPH PART',
608 .
' => DELETE SOLID ELEMENT AT NEXT CYCLE',i10)
610 &
' -- PARTICLE INTERACTING W/OTHER SPH PART OR SUBSET',
611 .
' => DELETE SOLID ELEMENT AT NEXT CYCLE',i10)