OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fvrezone.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!|| fvrezone0 ../engine/source/airbag/fvrezone.F
25!||--- called by ------------------------------------------------------
26!|| resol ../engine/source/engine/resol.F
27!||--- calls -----------------------------------------------------
28!|| fvrezone1 ../engine/source/airbag/fvrezone.F
29!||--- uses -----------------------------------------------------
30!|| fvbag_mod ../engine/share/modules/fvbag_mod.F
31!||====================================================================
32 SUBROUTINE fvrezone0(MONVOL, X)
33C-----------------------------------------------
34C M o d u l e s
35C-----------------------------------------------
36 USE fvbag_mod
37C-----------------------------------------------
38C I m p l i c i t T y p e s
39C-----------------------------------------------
40#include "implicit_f.inc"
41C-----------------------------------------------
42C C o m m o n B l o c k s
43C-----------------------------------------------
44#include "com01_c.inc"
45#include "com04_c.inc"
46#include "param_c.inc"
47C-----------------------------------------------
48C D u m m y A r g u m e n t s
49C-----------------------------------------------
50 INTEGER MONVOL(*)
51 my_real
52 . x(3,*)
53C-----------------------------------------------
54C L o c a l V a r i a b l e s
55C-----------------------------------------------
56 INTEGER K1, K2, KIBJET, KIBHOL, KIBALE, IADPOLH, IADPOLHO, N,
57 . ITYP, NNS, NTG, NBRIC, KI1, KI2, NNFV1, NTRFV1, NPOLY1,
58 . NPOLH1, NNFV0, NTRFV0, NPOLY0, NPOLH0, IFV, IMESH, NSTEP,
59 . ILVOUT, NNA
60C
61 k1=1
62 k2=1+nimv*nvolu
63 kibjet=k2+licbag
64 kibhol=kibjet+libagjet
65 kibale=kibhol+libaghol
66 ifv=0
67 DO n=1,nvolu
68 ityp=monvol(k1-1+2)
69 IF (ityp==6.OR.ityp==8) THEN
70 ifv=monvol(k1-1+45)
71 imesh=monvol(k1-1+56)
72 IF (imesh/=0) THEN
73C
74 nns=monvol(k1-1+32)
75 ntg=monvol(k1-1+33)
76 nbric=monvol(k1-1+35)
77 ki1=kibale+monvol(k1-1+31)
78 ki2=ki1+nns
79C
80 nnfv1=monvol(k1-1+46)
81 ntrfv1=monvol(k1-1+47)
82 npoly1=monvol(k1-1+48)
83 npolh1=monvol(k1-1+49)
84C
85 nnfv0=monvol(k1-1+50)
86 ntrfv0=monvol(k1-1+51)
87 npoly0=monvol(k1-1+52)
88 npolh0=monvol(k1-1+53)
89 nstep=monvol(k1-1+58)
90 ilvout=monvol(k1-1+44)
91 nna=monvol(k1-1+64)
92C
93 CALL fvrezone1(
94 1monvol(ki1), monvol(ki2), x, ntg, fvdata(ifv)%IFVNOD,
95 2fvdata(ifv)%RFVNOD, fvdata(ifv)%IFVTRI, fvdata(ifv)%IFVPOLY,
96 3fvdata(ifv)%IFVTADR, fvdata(ifv)%IFVPOLH, fvdata(ifv)%IFVPADR,
97 4fvdata(ifv)%MPOLH, fvdata(ifv)%QPOLH, fvdata(ifv)%EPOLH,
98 5fvdata(ifv)%PPOLH, fvdata(ifv)%RPOLH, fvdata(ifv)%GPOLH,
99 6fvdata_old(ifv)%IFVNOD, fvdata_old(ifv)%RFVNOD,
100 7fvdata_old(ifv)%IFVTRI, fvdata_old(ifv)%IFVPOLY,
101 8fvdata_old(ifv)%IFVTADR, fvdata_old(ifv)%IFVPOLH,
102 9fvdata_old(ifv)%IFVPADR, fvdata_old(ifv)%MPOLH,
103 afvdata_old(ifv)%QPOLH, fvdata_old(ifv)%EPOLH,
104 bfvdata_old(ifv)%PPOLH, fvdata_old(ifv)%RPOLH,
105 cfvdata_old(ifv)%GPOLH, nnfv1, ntrfv1, npolh1, nnfv0,
106 dntrfv0, npolh0, npoly1, npoly0, nstep,
107 emonvol(k1), fvdata(ifv)%CPAPOLH, fvdata(ifv)%CPBPOLH,
108 ffvdata(ifv)%CPCPOLH, fvdata(ifv)%RMWPOLH,
109 gfvdata_old(ifv)%CPAPOLH, fvdata_old(ifv)%CPBPOLH,
110 hfvdata_old(ifv)%CPCPOLH, fvdata_old(ifv)%RMWPOLH, ilvout, nns,nna,
111 iifv )
112 ENDIF
113 ENDIF
114 k1=k1+nimv
115 ENDDO
116C
117 RETURN
118 END
119C
120!||====================================================================
121!|| fvrezone1 ../engine/source/airbag/fvrezone.F
122!||--- called by ------------------------------------------------------
123!|| fvrezone0 ../engine/source/airbag/fvrezone.F
124!||--- calls -----------------------------------------------------
125!|| pinpolh ../engine/source/airbag/fvrezone.F
126!|| progbar_c ../common_source/sortie/progbar_c.c
127!|| spmd_fvb_gath ../engine/source/mpi/airbags/spmd_fvb_gath.F
128!||--- uses -----------------------------------------------------
129!|| fvbag_mod ../engine/share/modules/fvbag_mod.F
130!||====================================================================
131 SUBROUTINE fvrezone1(
132 1 IBUF , ELEM , X , NEL , IFVNOD1,
133 2 RFVNOD1 , IFVTRI1 , IFVPOLY1,
134 3 IFVTADR1, IFVPOLH1, IFVPADR1,
135 4 MPOLH1 , QPOLH1 , EPOLH1 ,
136 5 PPOLH1 , RPOLH1 , GPOLH1 ,
137 6 IFVNOD0 , RFVNOD0 ,
138 7 IFVTRI0 , IFVPOLY0,
139 8 IFVTADR0, IFVPOLH0,
140 9 IFVPADR0, MPOLH0 ,
141 A QPOLH0 , EPOLH0 ,
142 B PPOLH0 , RPOLH0 ,
143 C GPOLH0 , NNS1 , NNTR1 , NPOLH1, NNS0 ,
144 D NNTR0 , NPOLH0 , NPOLY1 , NPOLY0, NSTEP ,
145 E ID , CPAPOLH1, CPBPOLH1,
146 F CPCPOLH1, RMWPOLH1,
147 G CPAPOLH0, CPBPOLH0,
148 H CPCPOLH0, RMWPOLH0, ILVOUT , NNF , NNA ,
149 I IFV )
150C-----------------------------------------------
151C M o d u l e s
152C-----------------------------------------------
153 USE fvbag_mod
154C-----------------------------------------------
155C I m p l i c i t T y p e s
156C-----------------------------------------------
157#include "implicit_f.inc"
158C-----------------------------------------------
159C C o m m o n B l o c k s
160C-----------------------------------------------
161#include "units_c.inc"
162#include "com01_c.inc"
163#include "task_c.inc"
164C-----------------------------------------------
165C D u m m y A r g u m e n t s
166C-----------------------------------------------
167 INTEGER IBUF(*), ELEM(3,*), NEL, IFVNOD1(3,*), IFVTRI1(6,*),
168 . IFVPOLY1(*), IFVTADR1(*), IFVPOLH1(*), IFVPADR1(*),
169 . IFVNOD0(3,*), IFVTRI0(6,*), IFVPOLY0(*), IFVTADR0(*),
170 . IFVPOLH0(*), IFVPADR0(*), NNS1, NNTR1, NPOLH1, NNS0,
171 . nntr0, npolh0, npoly1, npoly0, nstep, id, ilvout,
172 . nnf, nna, ifv
173 my_real
174 . x(3,*), rfvnod1(2,*), mpolh1(*), qpolh1(3,*), epolh1(*),
175 . ppolh1(*), rpolh1(*), gpolh1(*), rfvnod0(2,*),
176 . mpolh0(*), qpolh0(3,*), epolh0(*), ppolh0(*),
177 . rpolh0(*), gpolh0(*), cpapolh1(*), cpbpolh1(*),
178 . cpcpolh1(*), rmwpolh1(*), cpapolh0(*), cpbpolh0(*),
179 . cpcpolh0(*), rmwpolh0(*)
180C-----------------------------------------------
181C L o c a l V a r i a b l e s
182C-----------------------------------------------
183 INTEGER I, IEL, N1, N2, N3, NN1, NN2, NN3, J, JJ, K, KK, NNT,
184 . L, NN, PNTR0(NNTR0), PNTR1(NNTR1), NNT1, NNT0, LL, NNTI,
185 . inner, inn(3), nti, icut, nb, nnb, nbi, m, ii, i1, i2,
186 . nnsa
187 my_real
188 . ksi, eta, x1, y1, z1, x2, y2, z2, x3, y3, z3,
189 . px0(3,nns0), px1(3,nns1), x12, y12, z12, x13, y13, z13,
190 . nrx, nry, nrz, area2, tarea0(nntr0), norm0(3,nntr0),
191 . tarea1(nntr1), norm1(3,nntr1), volu1(npolh1), area,
192 . nx, ny, nz, xmin, xmax, ymin, ymax, zmin, zmax, xx, yy,
193 . zz, bbox0(6,npolh0), bbox1(6,npolh1), crit, xmin1, xmax1,
194 . ymin1, ymax1, zmin1, zmax1, xmin0, xmax0, ymin0,
195 . ymax0, zmin0, zmax0, volu0(npolh0), dxb, dyb, dzb, volb,
196 . vol, volg, volt, rr, mass0, qx0, qy0, qz0, ener0, mass1,
197 . qx1, qy1, qz1, ener1, gama, fac, xxx(3,nnf), xxxa(3,nna)
198 INTEGER, ALLOCATABLE :: PTRI1(:,:), PTRI0(:,:), TCUT(:), ITAGT(:),
199 . BB(:,:), INB(:), INB_TMP(:), LISTB(:)
200 my_real
201 . , ALLOCATABLE :: xb(:,:), xxxsa(:,:)
202 CHARACTER NAME*18
203C
204 my_real
205 . dlamch
206 EXTERNAL dlamch
207C
208C 0 = Anciens polyhedres
209C 1 = Nouveaux polyhedres
210C
211 nnsa=fvspmd(ifv)%NNSA
212 ALLOCATE(xxxsa(3,nnsa))
213 IF (nspmd == 1) THEN
214C traitement necessaire pour rester p/on
215 DO i=1,fvspmd(ifv)%NN_L
216 i1=fvspmd(ifv)%IBUF_L(1,i)
217 i2=fvspmd(ifv)%IBUF_L(2,i)
218 xxx(1,i1)=x(1,i2)
219 xxx(2,i1)=x(2,i2)
220 xxx(3,i1)=x(3,i2)
221 ENDDO
222 DO i=1,fvspmd(ifv)%NNA_L
223 i1=fvspmd(ifv)%IBUFA_L(1,i)
224 i2=fvspmd(ifv)%IBUFA_L(2,i)
225 xxxa(1,i1)=x(1,i2)
226 xxxa(2,i1)=x(2,i2)
227 xxxa(3,i1)=x(3,i2)
228 ENDDO
229 DO i=1,fvspmd(ifv)%NNSA_L
230 i1=fvspmd(ifv)%IBUFSA_L(1,i)
231 i2=fvspmd(ifv)%IBUFSA_L(2,i)
232 xxxsa(1,i1)=x(1,i2)
233 xxxsa(2,i1)=x(2,i2)
234 xxxsa(3,i1)=x(3,i2)
235 ENDDO
236c DO I=1,NNF
237c II=IBUF(I)
238c XXX(1,I)=X(1,II)
239c XXX(2,I)=X(2,II)
240c XXX(3,I)=X(3,II)
241c ENDDO
242 ELSE
243 CALL spmd_fvb_gath(ifv, x, xxx, xxxa, xxxsa,
244 . 2 )
245 IF (ispmd/=fvspmd(ifv)%PMAIN-1) RETURN
246 ENDIF
247C Critere pour test d'interiorite
248 crit=ten
249C---------------------------------------------------
250C Preliminaires
251C---------------------------------------------------
252C Coordonnees des noeuds
253 DO i=1,nns0
254 IF (ifvnod0(1,i)==1) THEN
255 iel=ifvnod0(2,i)
256 ksi=rfvnod0(1,i)
257 eta=rfvnod0(2,i)
258C
259 n1=elem(1,iel)
260 n2=elem(2,iel)
261 n3=elem(3,iel)
262 x1=xxx(1,n1)
263 x2=xxx(1,n2)
264 x3=xxx(1,n3)
265 y1=xxx(2,n1)
266 y2=xxx(2,n2)
267 y3=xxx(2,n3)
268 z1=xxx(3,n1)
269 z2=xxx(3,n2)
270 z3=xxx(3,n3)
271 px0(1,i)=(one-ksi-eta)*x1+ksi*x2+eta*x3
272 px0(2,i)=(one-ksi-eta)*y1+ksi*y2+eta*y3
273 px0(3,i)=(one-ksi-eta)*z1+ksi*z2+eta*z3
274 ELSEIF (ifvnod0(1,i)==2) THEN
275 ii=ifvnod0(2,i)
276c IF (NSPMD == 1) THEN
277c PX0(1,I)=X(1,II)
278c PX0(2,I)=X(2,II)
279c PX0(3,I)=X(3,II)
280c ELSE
281 px0(1,i)=xxxsa(1,ii)
282 px0(2,i)=xxxsa(2,ii)
283 px0(3,i)=xxxsa(3,ii)
284c ENDIF
285 ENDIF
286 ENDDO
287 DO i=1,nns0
288 IF (ifvnod0(1,i)==3) THEN
289 i1=ifvnod0(2,i)
290 i2=ifvnod0(3,i)
291 fac=rfvnod0(1,i)
292 px0(1,i)=fac*px0(1,i1)+(one-fac)*px0(1,i2)
293 px0(2,i)=fac*px0(2,i1)+(one-fac)*px0(2,i2)
294 px0(3,i)=fac*px0(3,i1)+(one-fac)*px0(3,i2)
295 ENDIF
296 ENDDO
297C
298 DO i=1,nns1
299 IF (ifvnod1(1,i)==1) THEN
300 iel=ifvnod1(2,i)
301 ksi=rfvnod1(1,i)
302 eta=rfvnod1(2,i)
303C
304 n1=elem(1,iel)
305 n2=elem(2,iel)
306 n3=elem(3,iel)
307 x1=xxx(1,n1)
308 x2=xxx(1,n2)
309 x3=xxx(1,n3)
310 y1=xxx(2,n1)
311 y2=xxx(2,n2)
312 y3=xxx(2,n3)
313 z1=xxx(3,n1)
314 z2=xxx(3,n2)
315 z3=xxx(3,n3)
316 px1(1,i)=(one-ksi-eta)*x1+ksi*x2+eta*x3
317 px1(2,i)=(one-ksi-eta)*y1+ksi*y2+eta*y3
318 px1(3,i)=(one-ksi-eta)*z1+ksi*z2+eta*z3
319 ELSEIF (ifvnod1(1,i)==2) THEN
320 ii=ifvnod1(2,i)
321c IF (NSPMD == 1) THEN
322c PX1(1,I)=X(1,II)
323c PX1(2,I)=X(2,II)
324c PX1(3,I)=X(3,II)
325c ELSE
326 px1(1,i)=xxxsa(1,ii)
327 px1(2,i)=xxxsa(2,ii)
328 px1(3,i)=xxxsa(3,ii)
329c ENDIF
330 ENDIF
331 ENDDO
332 DO i=1,nns1
333 IF (ifvnod1(1,i)==3) THEN
334 i1=ifvnod1(2,i)
335 i2=ifvnod1(3,i)
336 fac=rfvnod1(1,i)
337 px1(1,i)=fac*px1(1,i1)+(one-fac)*px1(1,i2)
338 px1(2,i)=fac*px1(2,i1)+(one-fac)*px1(2,i2)
339 px1(3,i)=fac*px1(3,i1)+(one-fac)*px1(3,i2)
340 ENDIF
341 ENDDO
342 DEALLOCATE(xxxsa)
343C Normale et aire des triangles
344 DO i=1,nntr0
345 n1=ifvtri0(1,i)
346 n2=ifvtri0(2,i)
347 n3=ifvtri0(3,i)
348 x1=px0(1,n1)
349 y1=px0(2,n1)
350 z1=px0(3,n1)
351 x2=px0(1,n2)
352 y2=px0(2,n2)
353 z2=px0(3,n2)
354 x3=px0(1,n3)
355 y3=px0(2,n3)
356 z3=px0(3,n3)
357 x12=x2-x1
358 y12=y2-y1
359 z12=z2-z1
360 x13=x3-x1
361 y13=y3-y1
362 z13=z3-z1
363 nrx=y12*z13-z12*y13
364 nry=z12*x13-x12*z13
365 nrz=x12*y13-y12*x13
366 area2=sqrt(nrx**2+nry**2+nrz**2)
367 tarea0(i)=half*area2
368 IF (area2>0) THEN
369 norm0(1,i)=nrx/area2
370 norm0(2,i)=nry/area2
371 norm0(3,i)=nrz/area2
372 ELSE
373 norm0(1,i)=zero
374 norm0(2,i)=zero
375 norm0(3,i)=zero
376 ENDIF
377 ENDDO
378C
379 DO i=1,nntr1
380 n1=ifvtri1(1,i)
381 n2=ifvtri1(2,i)
382 n3=ifvtri1(3,i)
383 x1=px1(1,n1)
384 y1=px1(2,n1)
385 z1=px1(3,n1)
386 x2=px1(1,n2)
387 y2=px1(2,n2)
388 z2=px1(3,n2)
389 x3=px1(1,n3)
390 y3=px1(2,n3)
391 z3=px1(3,n3)
392 x12=x2-x1
393 y12=y2-y1
394 z12=z2-z1
395 x13=x3-x1
396 y13=y3-y1
397 z13=z3-z1
398 nrx=y12*z13-z12*y13
399 nry=z12*x13-x12*z13
400 nrz=x12*y13-y12*x13
401 area2=sqrt(nrx**2+nry**2+nrz**2)
402 tarea1(i)=half*area2
403 IF (area2>1) THEN
404 norm1(1,i)=nrx/area2
405 norm1(2,i)=nry/area2
406 norm1(3,i)=nrz/area2
407 ELSE
408 norm1(1,i)=zero
409 norm1(2,i)=zero
410 norm1(3,i)=zero
411 ENDIF
412 ENDDO
413C Volume des polyhedres
414 DO i=1,npolh0
415 volu0(i)=zero
416C Boucle sur les polygones du polyhedre
417 DO j=ifvpadr0(i),ifvpadr0(i+1)-1
418 jj=ifvpolh0(j)
419C Boucle sur les triangles du polygone
420 DO k=ifvtadr0(jj), ifvtadr0(jj+1)-1
421 kk=ifvpoly0(k)
422 area=tarea0(kk)
423 iel=ifvtri0(4,kk)
424 IF (iel>0) THEN
425 nx=norm0(1,kk)
426 ny=norm0(2,kk)
427 nz=norm0(3,kk)
428 ELSE
429 IF (ifvtri0(5,kk)==i) THEN
430 nx=norm0(1,kk)
431 ny=norm0(2,kk)
432 nz=norm0(3,kk)
433 ELSEIF (ifvtri0(6,kk)==i) THEN
434 nx=-norm0(1,kk)
435 ny=-norm0(2,kk)
436 nz=-norm0(3,kk)
437 ENDIF
438 ENDIF
439 n1=ifvtri0(1,kk)
440 x1=px0(1,n1)
441 y1=px0(2,n1)
442 z1=px0(3,n1)
443 volu0(i)=volu0(i)+third*area*(x1*nx+y1*ny+z1*nz)
444 ENDDO
445 ENDDO
446 ENDDO
447C
448 DO i=1,npolh1
449 volu1(i)=zero
450C Boucle sur les polygones du polyhedre
451 DO j=ifvpadr1(i),ifvpadr1(i+1)-1
452 jj=ifvpolh1(j)
453C Boucle sur les triangles du polygone
454 DO k=ifvtadr1(jj), ifvtadr1(jj+1)-1
455 kk=ifvpoly1(k)
456 area=tarea1(kk)
457 iel=ifvtri1(4,kk)
458 IF (iel>0) THEN
459 nx=norm1(1,kk)
460 ny=norm1(2,kk)
461 nz=norm1(3,kk)
462 ELSE
463 IF (ifvtri1(5,kk)==i) THEN
464 nx=norm1(1,kk)
465 ny=norm1(2,kk)
466 nz=norm1(3,kk)
467 ELSEIF (ifvtri1(6,kk)==i) THEN
468 nx=-norm1(1,kk)
469 ny=-norm1(2,kk)
470 nz=-norm1(3,kk)
471 ENDIF
472 ENDIF
473 n1=ifvtri1(1,kk)
474 x1=px1(1,n1)
475 y1=px1(2,n1)
476 z1=px1(3,n1)
477 volu1(i)=volu1(i)+third*area*(x1*nx+y1*ny+z1*nz)
478 ENDDO
479 ENDDO
480 ENDDO
481C Bounding box pour les polyhedres
482 DO i=1,npolh0
483 xmin=ep30
484 xmax=-ep30
485 ymin=ep30
486 ymax=-ep30
487 zmin=ep30
488 zmax=-ep30
489 nnt=0
490 DO j=ifvpadr0(i),ifvpadr0(i+1)-1
491 jj=ifvpolh0(j)
492 DO k=ifvtadr0(jj),ifvtadr0(jj+1)-1
493 nnt=nnt+1
494 kk=ifvpoly0(k)
495 DO l=1,3
496 nn=ifvtri0(l,kk)
497 xx=px0(1,nn)
498 yy=px0(2,nn)
499 zz=px0(3,nn)
500 xmin=min(xmin,xx)
501 xmax=max(xmax,xx)
502 ymin=min(ymin,yy)
503 ymax=max(ymax,yy)
504 zmin=min(zmin,zz)
505 zmax=max(zmax,zz)
506 ENDDO
507 ENDDO
508 ENDDO
509 bbox0(1,i)=xmin
510 bbox0(2,i)=xmax
511 bbox0(3,i)=ymin
512 bbox0(4,i)=ymax
513 bbox0(5,i)=zmin
514 bbox0(6,i)=zmax
515 pntr0(i)=nnt
516 ENDDO
517C
518 DO i=1,npolh1
519 xmin=ep30
520 xmax=-ep30
521 ymin=ep30
522 ymax=-ep30
523 zmin=ep30
524 zmax=-ep30
525 nnt=0
526 DO j=ifvpadr1(i),ifvpadr1(i+1)-1
527 jj=ifvpolh1(j)
528 DO k=ifvtadr1(jj),ifvtadr1(jj+1)-1
529 nnt=nnt+1
530 kk=ifvpoly1(k)
531 DO l=1,3
532 nn=ifvtri1(l,kk)
533 xx=px1(1,nn)
534 yy=px1(2,nn)
535 zz=px1(3,nn)
536 xmin=min(xmin,xx)
537 xmax=max(xmax,xx)
538 ymin=min(ymin,yy)
539 ymax=max(ymax,yy)
540 zmin=min(zmin,zz)
541 zmax=max(zmax,zz)
542 ENDDO
543 ENDDO
544 ENDDO
545 bbox1(1,i)=xmin
546 bbox1(2,i)=xmax
547 bbox1(3,i)=ymin
548 bbox1(4,i)=ymax
549 bbox1(5,i)=zmin
550 bbox1(6,i)=zmax
551 pntr1(i)=nnt
552 ENDDO
553C---------------------------------------------------
554C Rezone
555C---------------------------------------------------
556 nb=nstep**3
557 ALLOCATE(bb(8,nb))
558 nn=0
559 DO i=1,nstep
560 DO j=1,nstep
561 DO k=1,nstep
562 nn=nn+1
563 bb(1,nn)=(i-1)*(nstep+1)**2+(j-1)*(nstep+1)+k
564 bb(2,nn)=(i-1)*(nstep+1)**2+(j-1)*(nstep+1)+k+1
565 bb(3,nn)=(i-1)*(nstep+1)**2+j*(nstep+1)+k+1
566 bb(4,nn)=(i-1)*(nstep+1)**2+j*(nstep+1)+k
567 bb(5,nn)=i*(nstep+1)**2+(j-1)*(nstep+1)+k
568 bb(6,nn)=i*(nstep+1)**2+(j-1)*(nstep+1)+k+1
569 bb(7,nn)=i*(nstep+1)**2+j*(nstep+1)+k+1
570 bb(8,nn)=i*(nstep+1)**2+j*(nstep+1)+k
571 ENDDO
572 ENDDO
573 ENDDO
574C
575 IF (ilvout/=0) WRITE(istdo,'(A25,I8,A14)')
576 .' ** MONITORED VOLUME ID: ',id,' - REZONING **'
577 nnb=(nstep+1)**3
578 DO i=1,npolh1
579 mpolh1(i)=zero
580 qpolh1(1,i)=zero
581 qpolh1(2,i)=zero
582 qpolh1(3,i)=zero
583 epolh1(i)=zero
584 gpolh1(i)=zero
585 cpapolh1(i)=zero
586 cpbpolh1(i)=zero
587 cpcpolh1(i)=zero
588 rmwpolh1(i)=zero
589C
590 xmin1=bbox1(1,i)
591 xmax1=bbox1(2,i)
592 ymin1=bbox1(3,i)
593 ymax1=bbox1(4,i)
594 zmin1=bbox1(5,i)
595 zmax1=bbox1(6,i)
596C
597 nnt1=pntr1(i)
598 ALLOCATE(ptri1(3,nnt1))
599 nnt=0
600 DO j=ifvpadr1(i),ifvpadr1(i+1)-1
601 jj=ifvpolh1(j)
602 DO k=ifvtadr1(jj),ifvtadr1(jj+1)-1
603 nnt=nnt+1
604 kk=ifvpoly1(k)
605 IF (ifvtri1(4,kk)>0) THEN
606 ptri1(1,nnt)=ifvtri1(1,kk)
607 ptri1(2,nnt)=ifvtri1(2,kk)
608 ptri1(3,nnt)=ifvtri1(3,kk)
609 ELSEIF (ifvtri1(5,kk)==i) THEN
610 ptri1(1,nnt)=ifvtri1(1,kk)
611 ptri1(2,nnt)=ifvtri1(2,kk)
612 ptri1(3,nnt)=ifvtri1(3,kk)
613 ELSEIF (ifvtri1(6,kk)==i) THEN
614 ptri1(1,nnt)=ifvtri1(1,kk)
615 ptri1(2,nnt)=ifvtri1(3,kk)
616 ptri1(3,nnt)=ifvtri1(2,kk)
617 ENDIF
618 ENDDO
619 ENDDO
620C
621 ALLOCATE(xb(3,nnb), inb(nnb))
622 dxb=(xmax1-xmin1)/nstep
623 dyb=(ymax1-ymin1)/nstep
624 dzb=(zmax1-zmin1)/nstep
625 volb=dxb*dyb*dzb
626 nn=0
627 DO j=1,nstep+1
628 zz=zmin1+(j-1)*dzb
629 DO k=1,nstep+1
630 yy=ymin1+(k-1)*dyb
631 DO l=1,nstep+1
632 xx=xmin1+(l-1)*dxb
633 nn=nn+1
634 xb(1,nn)=xx
635 xb(2,nn)=yy
636 xb(3,nn)=zz
637 ENDDO
638 ENDDO
639 ENDDO
640 CALL pinpolh(nnt1, ptri1, xb, nnb, px1,
641 . inb, crit , xmin1, xmax1, ymin1,
642 . ymax1, zmin1, zmax1)
643C
644 volt=zero
645 jj=0
646 DO j=1,npolh0
647 xmin0=bbox0(1,j)
648 xmax0=bbox0(2,j)
649 ymin0=bbox0(3,j)
650 ymax0=bbox0(4,j)
651 zmin0=bbox0(5,j)
652 zmax0=bbox0(6,j)
653C
654 IF (xmax1<xmin0.OR.ymax1<ymin0.OR.zmax1<zmin0.OR.
655 . xmin1>xmax0.OR.ymin1>ymax0.OR.zmin1>zmax0)
656 . cycle
657C
658 nnt0=pntr0(j)
659 ALLOCATE(ptri0(3,nnt0))
660 nnt=0
661 DO k=ifvpadr0(j),ifvpadr0(j+1)-1
662 kk=ifvpolh0(k)
663 DO l=ifvtadr0(kk),ifvtadr0(kk+1)-1
664 nnt=nnt+1
665 ll=ifvpoly0(l)
666 IF (ifvtri0(4,ll)>0) THEN
667 ptri0(1,nnt)=ifvtri0(1,ll)
668 ptri0(2,nnt)=ifvtri0(2,ll)
669 ptri0(3,nnt)=ifvtri0(3,ll)
670 ELSEIF (ifvtri0(5,ll)==j) THEN
671 ptri0(1,nnt)=ifvtri0(1,ll)
672 ptri0(2,nnt)=ifvtri0(2,ll)
673 ptri0(3,nnt)=ifvtri0(3,ll)
674 ELSEIF (ifvtri0(6,ll)==j) THEN
675 ptri0(1,nnt)=ifvtri0(1,ll)
676 ptri0(2,nnt)=ifvtri0(3,ll)
677 ptri0(3,nnt)=ifvtri0(2,ll)
678 ENDIF
679 ENDDO
680 ENDDO
681C
682 ALLOCATE(inb_tmp(nnb))
683 CALL pinpolh(nnt0, ptri0, xb, nnb, px0,
684 . inb_tmp, crit, xmin0, xmax0, ymin0,
685 . ymax0, zmin0, zmax0)
686 DO k=1,nnb
687 inb_tmp(k)=inb_tmp(k)*inb(k)
688 ENDDO
689 ALLOCATE(listb(nb))
690 nbi=0
691 vol=zero
692 DO k=1,nb
693 nn=0
694 DO l=1,8
695 ll=bb(l,k)
696 nn=nn+inb_tmp(ll)
697 ENDDO
698 IF (nn>0) THEN
699 nbi=nbi+1
700 listb(nbi)=k
701 ENDIF
702 vol=vol+nn*volb/eight
703 ENDDO
704C
705 rr=vol/volu0(j)
706 mpolh1(i)=mpolh1(i)+rr*mpolh0(j)
707 qpolh1(1,i)=qpolh1(1,i)+rr*qpolh0(1,j)
708 qpolh1(2,i)=qpolh1(2,i)+rr*qpolh0(2,j)
709 qpolh1(3,i)=qpolh1(3,i)+rr*qpolh0(3,j)
710 epolh1(i)=epolh1(i)+rr*epolh0(j)
711 gpolh1(i)=gpolh1(i)+rr*gpolh0(j)
712 cpapolh1(i)=cpapolh1(i)+rr*cpapolh0(j)
713 cpbpolh1(i)=cpbpolh1(i)+rr*cpbpolh0(j)
714 cpcpolh1(i)=cpcpolh1(i)+rr*cpcpolh0(j)
715 rmwpolh1(i)=rmwpolh1(i)+rr*rmwpolh0(j)
716C
717 volt=volt+vol
718C
719 DEALLOCATE(ptri0, inb_tmp)
720 ENDDO
721C
722 DEALLOCATE(ptri1, xb, inb)
723C
724 IF (ilvout<0) CALL progbar_c(i,npolh1)
725 ENDDO
726 DEALLOCATE(bb)
727C
728 mass0=zero
729 qx0=zero
730 qy0=zero
731 qz0=zero
732 ener0=zero
733 mass1=zero
734 qx1=zero
735 qy1=zero
736 qz1=zero
737 ener1=zero
738 DO i=1,npolh0
739 mass0=mass0+mpolh0(i)
740 qx0=qx0+qpolh0(1,i)
741 qy0=qy0+qpolh0(2,i)
742 qz0=qz0+qpolh0(3,i)
743 ener0=ener0+epolh0(i)
744 ENDDO
745 DO i=1,npolh1
746 mass1=mass1+mpolh1(i)
747 qx1=qx1+qpolh1(1,i)
748 qy1=qy1+qpolh1(2,i)
749 qz1=qz1+qpolh1(3,i)
750 ener1=ener1+epolh1(i)
751 ENDDO
752
753 WRITE(istdo,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
754 . ' INITIAL MASS: ',mass0,' REZONED MASS: ',mass1,
755 . ' ERR: ',min(abs((mass1-mass0)/mass0*hundred),99.99d0),'%'
756 WRITE(istdo,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
757 . ' INITIAL QX : ',qx0,' REZONED QX : ',qx1,
758 . ' ERR: ',min(abs((qx1-qx0)/qx0*hundred),99.99d0),'%'
759 WRITE(istdo,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
760 . ' INITIAL QY : ',qy0,' REZONED QY : ',qy1,
761 . ' ERR: ',min(abs((qy1-qy0)/qy0*hundred),99.99d0),'%'
762 WRITE(istdo,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
763 . ' INITIAL QZ : ',qz0,' REZONED QZ : ',qz1,
764 . ' ERR: ',min(abs((qz1-qz0)/qz0*hundred),99.99d0),'%'
765 WRITE(istdo,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
766 . ' INITIAL ENER: ',ener0,' REZONED ENER: ',ener1,
767 . ' ERR: ',min(abs((ener1-ener0)/ener0*hundred),99.99d0),'%'
768 WRITE(istdo,*)
769 WRITE(iout,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
770 . ' INITIAL MASS: ',mass0,' REZONED MASS: ',mass1,
771 . ' ERR: ',min(abs((mass1-mass0)/mass0*hundred),99.99d0),'%'
772 WRITE(iout,'(A18,G11.4,A15,G11.4,A6,F5.2,A1)')
773 . ' INITIAL QX : ',qx0,' REZONED QX : ',qx1,
774 . ' ERR: ',min(abs((qx1-qx0)/qx0*hundred),99.99d0),'%'
775 WRITE(IOUT,'(a18,g11.4,a15,g11.4,a6,f5.2,a1)')
776 . ' initial qy : ',QY0,' rezoned qy : ',QY1,
777 . ' err: ',MIN(ABS((QY1-QY0)/QY0*HUNDRED),99.99D0),'%'
778 WRITE(IOUT,'(a18,g11.4,a15,g11.4,a6,f5.2,a1)')
779 . ' initial qz : ',QZ0,' rezoned qz : ',QZ1,
780 . ' err: ',MIN(ABS((QZ1-QZ0)/QZ0*HUNDRED),99.99D0),'%'
781 WRITE(IOUT,'(a18,g11.4,a15,g11.4,a6,f5.2,a1)')
782 . ' initial ener: ',ENER0,' rezoned ener: ',ENER1,
783 . ' err: ',MIN(ABS((ENER1-ENER0)/ENER0*HUNDRED),99.99D0),'%'
784 WRITE(IOUT,*)
785C
786 DO I=1,NPOLH1
787 GAMA=GPOLH1(I)
788 RPOLH1(I)=MPOLH1(I)/VOLU1(I)
789 PPOLH1(I)=(GAMA-ONE)*EPOLH1(I)/VOLU1(I)
790 ENDDO
791C
792 RETURN
793 END
794C
795!||====================================================================
796!|| pinpolh ../engine/source/airbag/fvrezone.F
797!||--- called by ------------------------------------------------------
798!|| fvrezone1 ../engine/source/airbag/fvrezone.F
799!||====================================================================
800 SUBROUTINE PINPOLH(NEL , ELEM, XB , NNB , X ,
801 . IN , TOLE, XMIN, XMAX, YMIN,
802 . YMAX, ZMIN, ZMAX)
803C-----------------------------------------------
804C I m p l i c i t T y p e s
805C-----------------------------------------------
806#include "implicit_f.inc"
807C-----------------------------------------------
808C D u m m y A r g u m e n t s
809C-----------------------------------------------
810 INTEGER NEL, ELEM(3,*), NNB, IN(*)
811 my_real
812 . XB(3,*), X(3,*), TOLE, XMIN, XMAX, YMIN, YMAX,
813 . ZMIN, ZMAX
814C-----------------------------------------------
815C L o c a l V a r i a b l e s
816C-----------------------------------------------
817 INTEGER I, IEL, N1, N2, N3, II
818 my_real
819 . XX, YY, ZZ, ST, X1, Y1, Z1, X2, Y2, Z2, X3, Y3, Z3,
820 . VX1, VY1, VZ1, VX2, VY2, VZ2, VX3, VY3, VZ3, NX1, NY1,
821 . NZ1, NX2, NY2, NZ2, NX3, NY3, NZ3, RR1, RR2, RR3, SS1,
822 . SS2, SS3, SAREA, NX, NY, NZ, XC, YC, ZC, SS
823C
824 DO I=1,NNB
825 XX=XB(1,I)
826 YY=XB(2,I)
827 ZZ=XB(3,I)
828 IN(I)=0
829C
830.OR..OR. IF (XX<XMINXX>XMAX
831.OR..OR. . YY<YMINYY>YMAX
832.OR. . ZZ<ZMINZZ>ZMAX) CYCLE
833C
834 ST=ZERO
835 DO IEL=1,NEL
836 N1=ELEM(1,IEL)
837 N2=ELEM(2,IEL)
838 N3=ELEM(3,IEL)
839 X1=X(1,N1)
840 Y1=X(2,N1)
841 Z1=X(3,N1)
842 X2=X(1,N2)
843 Y2=X(2,N2)
844 Z2=X(3,N2)
845 X3=X(1,N3)
846 Y3=X(2,N3)
847 Z3=X(3,N3)
848 VX1=X1-XX
849 VY1=Y1-YY
850 VZ1=Z1-ZZ
851 VX2=X2-XX
852 VY2=Y2-YY
853 VZ2=Z2-ZZ
854 VX3=X3-XX
855 VY3=Y3-YY
856 VZ3=Z3-ZZ
857C
858 NX1=VY1*VZ2-VZ1*VY2
859 NY1=VZ1*VX2-VX1*VZ2
860 NZ1=VX1*VY2-VY1*VX2
861 NX2=VY2*VZ3-VZ2*VY3
862 NY2=VZ2*VX3-VX2*VZ3
863 NZ2=VX2*VY3-VY2*VX3
864 NX3=VY3*VZ1-VZ3*VY1
865 NY3=VZ3*VX1-VX3*VZ1
866 NZ3=VX3*VY1-VY3*VX1
867 RR1=SQRT(NX1**2+NY1**2+NZ1**2)
868 RR2=SQRT(NX2**2+NY2**2+NZ2**2)
869 RR3=SQRT(NX3**2+NY3**2+NZ3**2)
870 SS1=-(NX1*NX2+NY1*NY2+NZ1*NZ2)/RR1/RR2
871 SS2=-(NX2*NX3+NY2*NY3+NZ2*NZ3)/RR2/RR3
872 SS3=-(NX3*NX1+NY3*NY1+NZ3*NZ1)/RR3/RR1
873 IF (SS1>ONE) SS1=ONE
874 IF (SS1<-ONE) SS1=-ONE
875 IF (SS2>ONE) SS2=ONE
876 IF (SS2<-ONE) SS2=-ONE
877 IF (SS3>ONE) SS3=ONE
878 IF (SS3<-ONE) SS3=-ONE
879 SAREA=ACOS(SS1)+ACOS(SS2)+ACOS(SS3)-PI
880C
881 VX1=X2-X1
882 VY1=Y2-Y1
883 VZ1=Z2-Z1
884 VX2=X3-X1
885 VY2=Y3-Y1
886 VZ2=Z3-Z1
887 NX=VY1*VZ2-VZ1*VY2
888 NY=VZ1*VX2-VX1*VZ2
889 NZ=VX1*VY2-VY1*VX2
890 XC=THIRD*(X1+X2+X3)
891 YC=THIRD*(Y1+Y2+Y3)
892 ZC=THIRD*(Z1+Z2+Z3)
893 SS=NX*(XC-XX)+NY*(YC-YY)+NZ*(ZC-ZZ)
894 IF (SS<ZERO) SAREA=-SAREA
895 ST=ST+SAREA
896 ENDDO
897C
898 IF (ABS(ST)>TOLE) IN(I)=1
899 ENDDO
900C
901 RETURN
902 END
my_real
#define my_real
Definition cppsort.cpp:32
fvrezone1
subroutine fvrezone1(ibuf, elem, x, nel, ifvnod1, rfvnod1, ifvtri1, ifvpoly1, ifvtadr1, ifvpolh1, ifvpadr1, mpolh1, qpolh1, epolh1, ppolh1, rpolh1, gpolh1, ifvnod0, rfvnod0, ifvtri0, ifvpoly0, ifvtadr0, ifvpolh0, ifvpadr0, mpolh0, qpolh0, epolh0, ppolh0, rpolh0, gpolh0, nns1, nntr1, npolh1, nns0, nntr0, npolh0, npoly1, npoly0, nstep, id, cpapolh1, cpbpolh1, cpcpolh1, rmwpolh1, cpapolh0, cpbpolh0, cpcpolh0, rmwpolh0, ilvout, nnf, nna, ifv)
Definition fvrezone.F:150
pinpolh
subroutine pinpolh(nel, elem, xb, nnb, x, in, tole, xmin, xmax, ymin, ymax, zmin, zmax)
Definition fvrezone.F:803
fvrezone0
subroutine fvrezone0(monvol, x)
Definition fvrezone.F:33
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
dlamch
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
ymax
subroutine ymax(idn, fac, npc, pld, stiffmin, stiffmax, stiffini, stiffavg)
Definition law100_upd.F:272
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
id
initmumps id
Definition multiplerhs_example.m:4
fvbag_mod
Definition fvbag_mod.F:44
fvbag_mod::fvdata_old
type(fvbag_data), dimension(:), allocatable fvdata_old
Definition fvbag_mod.F:193
fvbag_mod::fvspmd
type(fvbag_spmd), dimension(:), allocatable fvspmd
Definition fvbag_mod.F:129
fvbag_mod::fvdata
type(fvbag_data), dimension(:), allocatable fvdata
Definition fvbag_mod.F:128
progbar_c
void progbar_c(int *icur, int *imax)
Definition progbar_c.c:30
spmd_fvb_gath
subroutine spmd_fvb_gath(ifv, x, xxx, xxxa, xxxsa, ido)
Definition spmd_fvb_gath.F:40