33 SUBROUTINE volpre(IVOLU ,RVOLU ,NJET ,IBAGJET ,RBAGJET ,
34 2 NSENSOR ,SENSOR_TAB,X ,V ,A ,
35 3 NORMAL ,NPC ,TF ,NN ,SURF_NODES,
36 4 FEXT ,H3D_DATA,SURF_ELTYP,SURF_ELEM,WFEXT,PYTHON)
47#include "implicit_f.inc"
60 INTEGER ,
INTENT(IN) :: NSENSOR
61 INTEGER IVOLU(*),NJET,IBAGJET(NIBJET,*),NPC(*),NN,SURF_NODES(NN,4),SURF_ELTYP(NN),SURF_ELEM(NN)
62 my_real RVOLU(*),RBAGJET(NRBJET,*), X(3,*), V(3,*), A(3,*), NORMAL(3,NN),TF(*),FEXT(3,*)
64 TYPE () ,
DIMENSION(NSENSOR) ,
INTENT(IN) :: SENSOR_TAB
65 DOUBLE PRECISION,
INTENT(INOUT) :: WFEXT
66 TYPE(python_),
intent(inout) :: PYTHON
70 INTEGER I,II,NJ1,NJ2,NJ3,IPT,IPA,IPZ,
73 my_real pres,fx,fy,fz,p,pext,dp,q,pj,pj0,dydx,xx,yy,zz,
74 . theta,nx1,nx2,ny1,ny2,nz1,nz2,an1,an,
75 . xm,ym,zm,nx3,ny3,nz3,
76 . ps1,ps2,ps3,an2,an3,det,mu,facr,
77 . api,dist,tstart,ts,tmpo,fpt,fpa,fpz,scalt,scala,scald
90 IF(surf_eltyp(i) == 3)
THEN
118 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
119 fext(1,n1) = fext(1,n1)+fx
120 fext(2,n1) = fext(2,n1)+fy
121 fext(3,n1) = fext(3,n1)+fz
122 fext(1,n2) = fext(1,n2)+fx
123 fext(2,n2) = fext(2,n2)+fy
124 fext(3,n2) = fext(3,n2)+fz
125 fext(1,n3) = fext(1,n3)+fx
126 fext(2,n3) = fext(2,n3)+fy
127 fext(3,n3) = fext(3,n3)+fz
128 fext(1,n4) = fext(1,n4)+fx
129 fext(2,n4) = fext(2,n4)+fy
130 fext(3,n4) = fext(3,n4)+fz
132 ELSEIF(surf_eltyp(i) == 7)
THEN
134 ii=surf_elem(i) + numelc
155 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
156 fext(1,n1) = fext(1,n1)+fx
157 fext(2,n1) = fext(2,n1)+fy
158 fext(3,n1) = fext(3,n1)+fz
159 fext(1,n2) = fext(1,n2)+fx
160 fext(2,n2) = fext(2,n2)+fy
161 fext(3,n2) = fext(3,n2)+fz
162 fext(1,n3) = fext(1,n3)+fx
163 fext(2,n3) = fext(2,n3)+fy
164 fext(3,n3) = fext(3,n3)+fz
194 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
195 fext(1,n1) = fext(1,n1)+fx
196 fext(2,n1) = fext(2,n1)+fy
197 fext(3,n1) = fext(3,n1)+fz
198 fext(1,n2) = fext(1,n2)+fx
199 fext(2,n2) = fext(2,n2)+fy
200 fext(3,n2) = fext(3,n2)+fz
201 fext(1,n3) = fext(1,n3)+fx
202 fext(2,n3) = fext(2,n3)+fy
203 fext(3,n3) = fext(3,n3)+fz
204 fext(1,n4) = fext(1,n4)+fx
205 fext(2,n4) = fext(2,n4)+fy
206 fext(3,n4) = fext(3,n4)+fz
212 IF(ityp /= 4 .AND. ityp /= 5 .AND. ityp /= 7 .AND. ityp /= 9)
RETURN
231 tstart=sensor_tab(isens)%TSTART
234 IF(nj1 /= 0 .AND. tt >= tstart)
THEN
236 pj0=fpt*finter_mixed(python,nfunct,ipt,ts*scalt,npc
238 IF(surf_eltyp(i)==3)
THEN
248 xx = fourth*(x(1,n1)+x(1,n2)+x(1,n3)+x(1,n4))
249 yy = fourth*(x(2,n1)+x(2,n2)+x(2,n3)+x(2,n4))
250 zz = fourth*(x(3,n1)+x(3,n2)+x(3,n3)+x(3,n4))
252 xm=half*(x(1,nj1)+x(1,nj3))
253 ym=half*(x(2,nj1)+x(2,nj3))
254 zm=half*(x(3,nj1)+x(3,nj3))
270 ps1 = nx1*nx2+ny1*ny2+nz1*nz2
271 ps2 = nx2*nx3+ny2*ny3+nz2*nz3
272 ps3 = nx1*nx3+ny1*ny3+nz1*nz3
273 an2 = nx2*nx2+ny2*ny2+nz2*nz2
274 an3 = nx3*nx3+ny3*ny3+nz3*nz3
275 det = ps2*ps2-an2*an3
277 mu = (ps2*ps1-ps3*an2)/sign(
max(em30,abs(det)),det)
279 facr =
min(one,
max(-one,mu))
284 an1 = (nx1**2+ny1**2+nz1**2)
285 an =
max(em30,sqrt((normal(1,i)**2+normal(2,i)**2+normal(3,i)**2)*an1))
286 pj=pj*
max(zero,(nx1*normal(1,i)+ny1*normal(2,i)+nz1*normal(3,i))/an)
287 an =
max(em30,sqrt((nx2**2+ny2**2+nz2**2)*an1))
288 tmpo = (nx1*nx2+ny1*ny2+nz1*nz2)/an
289 tmpo = sign(
min(one,abs(tmpo)),tmpo)
291 pj=pj*fpa*finter_mixed(python,nfunct,ipa,theta*scala,npc,tf)
293 IF(ipz/=0)pj=pj*fpz*finter_mixed(python,nfunct,ipz,dist*scald,npc,tf)
315 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
316 fext(1,n1) = fext(1,n1)+fx
318 fext(3,n1) = fext(3,n1)+fz
319 fext(1,n2) = fext(1,n2)+fx
320 fext(2,n2) = fext(2,n2)+fy
321 fext(3,n2) = fext(3,n2)+fz
322 fext(1,n3) = fext(1,n3)+fx
323 fext(2,n3) = fext(2,n3)+fy
324 fext(3,n3) = fext(3,n3)+fz
325 fext(1,n4) = fext(1,n4)+fx
326 fext(2,n4) = fext(2,n4)+fy
327 fext(3,n4) = fext(3,n4)+fz
330 wfext=wfext+dt1*(fx*(v(1,n1)+v(1,n2)+v(1,n3)+v(1,n4))
331 + +fy*(v(2,n1)+v(2,n2)+v(2,n3)+v(2,n4))
332 + +fz*(v(3,n1)+v(3,n2)+v(3,n3)+v(3,n4)))
333 ELSEIF(surf_eltyp(i)==7)
THEN
340 ii=surf_elem(i) + numelc
342 xx = (x(1,n1)+x(1,n2)+x(1,n3)) * third
343 yy = (x(2,n1)+x(2,n2)+x(2,n3)) * third
344 zz = (x(3,n1)+x(3,n2)+x(3,n3)) * third
346 xm=0.5*(x(1,nj1)+x(1,nj3))
347 ym=0.5*(x(2,nj1)+x(2,nj3))
348 zm=0.5*(x(3,nj1)+x(3,nj3))
363 ps1 = nx1*nx2+ny1*ny2+nz1*nz2
364 ps2 = nx2*nx3+ny2*ny3+nz2*nz3
365 ps3 = nx1*nx3+ny1*ny3+nz1*nz3
366 an2 = nx2*nx2+ny2*ny2+nz2*nz2
367 an3 = nx3*nx3+ny3*ny3+nz3*nz3
368 det = ps2*ps2-an2*an3
370 mu = (ps2*ps1-ps3*an2)/sign(
max(em30,abs(det)),det)
372 facr =
min(one,
max(-one,mu))
377 an1 = (nx1**2+ny1**2+nz1**2)
378 an =
max(em30,sqrt((normal(1,i)**2+normal(2,i)**2+normal(3,i)**2)*an1))
379 pj=pj*
max(zero,(nx1*normal(1,i)+ny1*normal(2,i)+nz1*normal(3,i))/an)
380 an =
max(em30,sqrt((nx2**2+ny2**2+nz2**2)*an1))
381 tmpo = (nx1*nx2+ny1*ny2+nz1*nz2)/an
382 tmpo = sign(
min(one,abs(tmpo)),tmpo)
384 pj=pj*fpa*finter_mixed(python,nfunct,ipa,theta*scala,npc,tf)
386 IF(ipz/=0)pj=pj*fpz*finter_mixed(python,nfunct,ipz,dist
404 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+
405 . anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
406 fext(1,n1) = fext(1,n1)+fx
407 fext(2,n1) = fext(2,n1)+fy
408 fext(3,n1) = fext(3,n1)+fz
409 fext(1,n2) = fext(1,n2)+fx
410 fext(2,n2) = fext(2,n2)+fy
411 fext(3,n2) = fext(3,n2)+fz
412 fext(1,n3) = fext(1,n3)+fx
413 fext(2,n3) = fext(2,n3)+fy
414 fext(3,n3) = fext(3,n3)+fz
417 wfext=wfext+dt1*(fx*(v(1,n1)+v(1,n2)+v(1,n3))
418 + +fy*(v(2,n1)+v(2,n2)+v(2,n3))
419 + +fz*(v(3,n1)+v(3,n2)+v(3,n3)))
434 xx = fourth*(x(1,n1)+x(1,n2)+x(1,n3)+x(1,n4))
435 yy = fourth*(x(2,n1)+x(2,n2)+x(2,n3)+x(2,n4))
436 zz = fourth*(x(3,n1)+x(3,n2)+x(3,n3)+x(3,n4))
438 xm=half*(x(1,nj1)+x(1,nj3))
439 ym=half*(x(2,nj1)+x(2,nj3))
440 zm=half*(x(3,nj1)+x(3,nj3))
455 ps1 = nx1*nx2+ny1*ny2+nz1*nz2
456 ps2 = nx2*nx3+ny2*ny3+nz2*nz3
457 ps3 = nx1*nx3+ny1*ny3+nz1*nz3
458 an2 = nx2*nx2+ny2*ny2+nz2*nz2
459 an3 = nx3*nx3+ny3*ny3+nz3*nz3
460 det = ps2*ps2-an2*an3
462 mu = (ps2*ps1-ps3*an2)/sign(
max(em30,abs(det)),det)
464 facr =
min(one,
max(-one,mu))
469 an1 = (nx1**2+ny1**2+nz1**2)
470 an =
max(em30,sqrt((normal(1,i)**2+normal(2,i)**2+normal(3,i)**2)*an1))
471 pj=pj*
max(zero,(nx1*normal(1,i)+ny1*normal(2,i)+nz1*normal(3,i))/an)
472 an =
max(em30,sqrt((nx2**2+ny2**2+nz2**2)*an1))
473 tmpo = (nx1*nx2+ny1*ny2+nz1*nz2)/an
474 tmpo = sign(
min(one,abs(tmpo)),tmpo)
476 pj=pj*fpa*finter_mixed(python,nfunct,ipa,theta*scala,npc,tf)
478 IF(ipz/=0)pj=pj*fpz*finter_mixed(python,nfunct,ipz,dist*scald,npc,tf)
500 IF(anim_v(5)+outp_v(5)+h3d_data%N_VECT_FINT+
501 . anim_v(6)+outp_v(6)+h3d_data%N_VECT_FEXT >0)
THEN
502 fext(1,n1) = fext(1,n1)+fx
503 fext(2,n1) = fext(2,n1)+fy
504 fext(3,n1) = fext(3,n1)+fz
505 fext(1,n2) = fext(1,n2)+fx
506 fext(2,n2) = fext(2,n2)+fy
507 fext(3,n2) = fext(3,n2)+fz
508 fext(1,n3) = fext(1,n3)+fx
509 fext(2,n3) = fext(2,n3)+fy
510 fext(3,n3) = fext(3,n3)+fz
511 fext(1,n4) = fext(1,n4)+fx
512 fext(2,n4) = fext(2,n4)+fy
513 fext(3,n4) = fext(3,n4)+fz
516 wfext=wfext+dt1*(fx*(v(1,n1)+v(1,n2)+v(1,n3)+v(1,n4))
517 + +fy*(v(2,n1)+v(2,n2)+v(2,n3)+v(2,n4))
518 + +fz*(v(3,n1)+v(3,n2)+v(3,n3)+v(3,n4)))