OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i11trivox.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!|| i11trivox ../engine/source/interfaces/intsort/i11trivox.F
25!||--- called by ------------------------------------------------------
26!|| i11buce_vox ../engine/source/interfaces/intsort/i11buce.F
27!||--- calls -----------------------------------------------------
28!|| i11sto_vox ../engine/source/interfaces/intsort/i11sto.F
29!|| ireallocate ../engine/share/modules/realloc_mod.F
30!|| my_barrier ../engine/source/system/machine.F
31!||--- uses -----------------------------------------------------
32!|| realloc_mod ../engine/share/modules/realloc_mod.F
33!|| tri11 ../engine/share/modules/tri11_mod.F
34!|| tri7box ../engine/share/modules/tri7box.F
35!||====================================================================
36 SUBROUTINE i11trivox(
37 1 IRECTS, IRECTM , X , NRTM ,NRTSR ,
38 2 XYZM , II_STOK, CAND_S , CAND_M ,NSN4 ,
39 3 NOINT , TZINF , I_MEM , ESHIFT ,ADDCM ,
40 5 CHAINE, NRTS , ITAB , STFS ,STFM ,
41 6 IAUTO , VOXEL , NBX , NBY ,NBZ ,
42 7 ITASK , IFPEN , IFORM , GAPMIN ,DRAD ,
43 8 MARGE ,GAP_S , GAP_M , GAP_S_L ,GAP_M_L,
44 9 BGAPSMX, IGAP ,GAP ,FLAGREMNODE,KREMNODE,
45 1 REMNODE,DGAPLOAD)
46C============================================================================
47C M o d u l e s
48C-----------------------------------------------
49 USE realloc_mod
50 USE tri7box
51 USE tri11
52C-----------------------------------------------
53C I m p l i c i t T y p e s
54C-----------------------------------------------
55#include "implicit_f.inc"
56#include "comlock.inc"
57C-----------------------------------------------
58C G l o b a l P a r a m e t e r s
59C-----------------------------------------------
60#include "mvsiz_p.inc"
61C-----------------------------------------------
62C C o m m o n B l o c k s
63C-----------------------------------------------
64#include "param_c.inc"
65C-----------------------------------------------
66C M e s s a g e P a s s i n g
67C-----------------------------------------------
68#ifdef MPI
69#endif
70!-----------------------------------------------
71! SUBROUTINE AIM
72! ==============
73! VOXEL SEARCH to find couple (edge,edge) with penetration among all possible couples defined by secnd and main side.
74! Temporary found candidate are written in Temporary array PROV_S and PROV_M in order to optimise OpenMP performances.
75! There is no order.
76! PROV_S(i),PROV_M(i) : is a potential candidate couple because edges are geometrically near each other.
77! I11STO subroutine will compute if penetration is positive and if couple was not already stoked, in this case, candidate is stoked in CAND(S(i), CAND_M(i)
78!
79C-----------------------------------------------
80C D u m m y A r g u m e n t s
81C
82C NOM DESCRIPTION E/S
83C
84C ADD(2,*) TABLEAU DES ADRESSES E/S
85C 1.........ADRESSES NOEUDS
86C 2.........ADRESSES ELEMENTS
87C ZYZM(6,*) TABLEAU DES XYZMIN E/S
88C 1.........XMIN BOITE
89C 2.........YMIN BOITE
90C 3.........ZMIN BOITE
91C 4.........XMAX BOITE
92C 5.........YMAX BOITE
93C 6.........ZMAX BOITE
94C IRECTM(2,*) TABLEAU DES CONEC E
95C 1.........NODE 1 main EDGE
96C 2.........NODE 2 main EDGE
97C IRECTS(2,*) TABLEAU DES CONEC E
98C 1.........NODE 1 SECND EDGE
99C 2.........NODE 2 SECND EDGE
100C X(3,*) COORDONNEES NODALES E
101C II_STOK niveau de stockage des couples
102C candidats impact E/S
103C CAND_S boites resultats noeuds
104C CAND_M adresses des boites resultat elements
105C NOINT NUMERO USER DE L'INTERFACE
106C TZINF TAILLE ZONE INFLUENCE
107C VOXEL(*,*,*) VOXEL PARTIONNEMENT DE l'ESPACE (NBX+2,NBY+2,NBZ+2)
108C Stocke dans chaque voxel les edge secnd.
109C En pratique designe la premier arrete d'une liste chainee
110C MAX_ADD adresse maximum pour les tableaux chaines
111C NSN4 4*NSN TAILLE MAX ADMISE MAINTENANT POUR LES
112C COUPLES NOEUDS,ELT CANDIDATS
113C-----------------------------------------------
114C D u m m y A r g u m e n t s
115C-----------------------------------------------
116 INTEGER ::
117 . NRTM,NRTSR,ESHIFT,NRTS,IGAP,
118 . NSN4,NOINT,ITAB(*),NBX,NBY,NBZ,IAUTO,
119 . IRECTS(2,NRTS),IRECTM(2,NRTM)
120 INTEGER ITASK,IFORM
121 INTEGER, INTENT(INOUT) ::
122 . CAND_S(*),CAND_M(*),ADDCM(*),CHAINE(2,*),
123 . VOXEL(1:NBX+2,1:NBY+2,1:NBZ+2), I_MEM,IFPEN(*),II_STOK,
124 . FLAGREMNODE,KREMNODE(*),REMNODE(*)
125 my_real
126 . ,INTENT(IN) ::
127 . x(3,*),xyzm(6,*),
128 . stfs(nrts),stfm(nrtm), tzinf, gap
129 my_real , INTENT(IN) :: dgapload,drad
130 my_real
131 . gapmin,marge,bgapsmx,
132 . gap_s(*),gap_m(*), gap_s_l(*), gap_m_l(*)
133C-----------------------------------------------
134C L o c a l V a r i a b l e s
135C-----------------------------------------------
136 INTEGER
137 . I,J,SS1,SS2,IBUG,
138 . n1,n2,mm1,mm2, in1, in2, im1, im2, k,l,
139 . prov_s(2*mvsiz),prov_m(2*mvsiz), !tableau provisoire de candidats envoye a I11STOK
140 . ix1,iy1,iz1,ix2,iy2,iz2,
141 . ix,iy,iz, first_add,
142 . i_stok, i_stok_bak, iedg,
143 . prev_add, chain_add, current_add, !pour le balayage des tableaux chaines
144 . nedg, deja , max_add ,ii_stok0, m,remove_remote
145 INTEGER, DIMENSION(3) :: TMIN,TMAX
146 my_real
147 . XX1, XX2,
148 . XMIN, XMAX,YMIN, YMAX,ZMIN, ZMAX,
149 . YY1,YY2,ZZ1,ZZ2,
150 . AAA, DD,
151 . XMAX_EDGS(NRTS+NRTSR), XMIN_EDGS(NRTS+NRTSR), !cotes min/max des arrete seconds et mains
152 . YMAX_EDGS(NRTS+NRTSR), YMIN_EDGS(NRTS+NRTSR),
153 . ZMAX_EDGS(NRTS+NRTSR), ZMIN_EDGS(NRTS+NRTSR),
154 . xmax_edgm(nrtm), xmin_edgm(nrtm),
155 . ymax_edgm(nrtm), ymin_edgm(nrtm),
156 . zmax_edgm(nrtm), zmin_edgm(nrtm),
157 . xminb,yminb,zminb,xmaxb,ymaxb,zmaxb
158C-----------------------------------------------
159 INTEGER, DIMENSION(:), ALLOCATABLE :: TAGREMLINE
160C-----------------------------------------------
161C
162 IF(FLAGREMNODE==2) then
163 ALLOCATE(tagremline(nrts))
164 tagremline(1:nrts) = 0
165 ENDIF
166C
167 aaa = zero
168 !ATTENTION A POPTIONMISER POUR NE PAS FAIRE ONE RAZ COMPLET SI NRTS ==0
169 min_ix=nbx+2
170 min_iy=nby+2
171 min_iz=nbz+2
172 max_ix=1
173 max_iy=1
174 max_iz=1
175
176 !---------------------------------------------------------!
177 ! Allocation des tableaux chaines !
178 !---------------------------------------------------------!
179 IF(itask == 0)THEN
180 max_add = max(1,4*(nrts+nrtsr))
181 ALLOCATE(lchain_elem(1:max_add))
182 ALLOCATE(lchain_next(1:max_add))
183 ALLOCATE(lchain_last(1:max_add))
184 END IF
185
186 CALL my_barrier !all threads wait for allocation
187
188 IF(nrtm==0.OR.nrts==0)THEN
189 !ne pas reinitiliser tout le voxel s'il n'y a pas de candidat
190 min_ix=1
191 min_iy=1
192 min_iz=1
193 END IF
194
195 !---------------------------------------------------------!
196 ! recuperation des bornes du domaine !
197 !---------------------------------------------------------!
198 xmin = xyzm(1,1)
199 ymin = xyzm(2,1)
200 zmin = xyzm(3,1)
201 xmax = xyzm(4,1)
202 ymax = xyzm(5,1)
203 zmax = xyzm(6,1)
204c dev future: xminb plus grand que xmin...
205 xminb = xmin
206 yminb = ymin
207 zminb = zmin
208 xmaxb = xmax
209 ymaxb = ymax
210 zmaxb = zmax
211C=======================================================================
212C 1 Pour chaque edge, on marque les voxels occupes.
213C Le nombre d edge dans un voxels etant variable, on
214C utilise un tableau chaine.
215C Ces voxels representent le voisinage de l edge.
216C On cherchera ensuite toutes les entites interfacees
217C dans ce voisinage.
218C=======================================================================
219 IF(itask == 0)THEN
220
221 current_add=1 ! premiere adresse
222
223 DO i = 1,nrts !si besoin on peut inverser Main/Secnd
224
225 IF(stfs(i)==zero)cycle !on ne retient pas les facettes detruites
226
227 !-------------------------------------------!
228 ! Nodes ID for edge (N1,N2) !
229 !-------------------------------------------!
230 n1=irects(1,i)
231 n2=irects(2,i)
232 !-------------------------------------------!
233 ! Coordinates of the two nodes !
234 ! +Optimisation // recherche les noeuds !
235 !compris dans xmin xmax des elements du proc!
236 !-------------------------------------------!
237 xx1=x(1,n1)
238 xx2=x(1,n2)
239 xmax_edgs(i)=max(xx1,xx2); IF(xmax_edgs(i) < xmin) cycle
240 xmin_edgs(i)=min(xx1,xx2); IF(xmin_edgs(i) > xmax) cycle
241 yy1=x(2,n1)
242 yy2=x(2,n2)
243 ymax_edgs(i)=max(yy1,yy2); IF(ymax_edgs(i) < ymin) cycle
244 ymin_edgs(i)=min(yy1,yy2); IF(ymin_edgs(i) > ymax) cycle
245 zz1=x(3,n1)
246 zz2=x(3,n2)
247 zmax_edgs(i)=max(zz1,zz2); IF(zmax_edgs(i) < zmin) cycle
248 zmin_edgs(i)=min(zz1,zz2); IF(zmin_edgs(i) > zmax) cycle
249
250 !-------------------------------------------!
251 ! VOXEL OCCUPIED BY THE EDGE !
252 !-------------------------------------------!
253 !Voxel_lower_left_bound for this edge
254 ix1=int(nbx*(xmin_edgs(i)-xminb)/(xmaxb-xminb))
255 iy1=int(nby*(ymin_edgs(i)-yminb)/(ymaxb-yminb))
256 iz1=int(nbz*(zmin_edgs(i)-zminb)/(zmaxb-zminb))
257 ix1=max(1,2+min(nbx,ix1))
258 iy1=max(1,2+min(nby,iy1))
259 iz1=max(1,2+min(nbz,iz1))
260 !Voxel_upper_right_bound for this edge
261 ix2=int(nbx*(xmax_edgs(i)-xminb)/(xmaxb-xminb))
262 iy2=int(nby*(ymax_edgs(i)-yminb)/(ymaxb-yminb))
263 iz2=int(nbz*(zmax_edgs(i)-zminb)/(zmaxb-zminb))
264 ix2=max(1,2+min(nbx,ix2))
265 iy2=max(1,2+min(nby,iy2))
266 iz2=max(1,2+min(nbz,iz2))
267
268 !pour reset des voxel
269 min_ix = min(min_ix,ix1)
270 min_iy = min(min_iy,iy1)
271 min_iz = min(min_iz,iz1)
272 max_ix = max(max_ix,ix2)
273 max_iy = max(max_iy,iy2)
274 max_iz = max(max_iz,iz2)
275
276 !----------------------------------------------!
277 ! EDGE STORAGE FOR EACH VOXEL (CHAINED ARRAY) !
278 !----------------------------------------------!
279C
280C VOXEL(i,j,k) LCHAIN_LAST(FIRST)
281C +-----------+------------+
282C | =>FIRST | =>LAST |
283C +--+--------+--+---------+
284C | |
285C | |
286C | |
287C | | LCHAIN_ELEM(*) LCHAIN_NEXT(*)
288C | | +------------+-----------+
289C +-------------->| edge_id | iadd 3 | 1:FIRST --+
290C | +------------+-----------+ |
291C | | | | 2 |
292C | +------------+-----------+ |
293C | | edge_id | iadd 4 | 3 <-------+
294C | +------------+-----------+ |
295C | | edge_id | iadd 6 | 4 <-------+
296C | +------------+-----------+ |
297C | | | | 5 |
298C | +------------+-----------+ |
299C +-->| edge_id | 0 | 6:LAST <--+
300C +------------+-----------+
301C | | | MAX_ADD
302C +------------+-----------+
303C
304 !Pour tous les voxels qu'occupe la brique
305 DO iz = iz1,iz2
306 DO iy = iy1,iy2
307 DO ix = ix1,ix2
308
309 first_add = voxel(ix,iy,iz)
310
311 IF(first_add == 0)THEN
312 !voxel encore vide
313 voxel(ix,iy,iz) = current_add ! adresse dans le tableau chaine de la premiere eddge trouvee occupant le voxel
314 lchain_last(current_add) = current_add ! dernier=adresse pour l edge courante
315 lchain_elem(current_add) = i ! edge ID
316 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste !
317 ELSE
318 !voxel contenant deja une edge
319 prev_add = lchain_last(first_add) ! devient l'avant-dernier
320 lchain_last(first_add) = current_add ! maj du dernier
321 lchain_elem(current_add) = i ! edge ID
322 lchain_next(prev_add) = current_add ! maj du suivant 0 -> CURRENT_ADD
323 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste
324 ENDIF
325
326 current_add = current_add+1
327
328 IF( current_add>=max_add)THEN
329 !OPTIMISATION : suprresion du deallocate/GOTO debut.
330 !REALLOCATE SI PAS ASSEZ DE PLACE : inutile de recommencer de 1 a MAX_ADD-1, on poursuit de MAX_ADD a 2*MAX_ADD
331 max_add = 2 * max_add
332 !print *, "reallocate"
336 ENDIF
337
338 ENDDO !IX
339 ENDDO !IY
340 ENDDO !IZ
341
342 ENDDO !DO I=1,NRTS
343
344C=======================================================================
345C 2 Traiter les edges remote. Recuperer les 2 noeuds des
346C edge remotes qui sont dans les memes voxels
347C
348C a faire
349C
350C=======================================================================
351 DO i = nrts+1,nrts+nrtsr !si besoin on peut inverser Main/Secnd
352c IF(STFS(I)==ZERO)CYCLE !on ne retient pas les facettes detruites, deja fait dans SPMD_MACH::spmd_tri11vox
353 j=i-nrts
354 !-------------------------------------------!
355 ! Coordinates of the two nodes !
356 ! +Optimisation // recherche les noeuds !
357 !compris dans xmin xmax des elements du proc!
358 !-------------------------------------------!
359 xx1=xrem(1,j)
360 xx2=xrem(8,j)
361 xmax_edgs(i)=max(xx1,xx2) ; IF(xmax_edgs(i) < xmin) cycle
362 xmin_edgs(i)=min(xx1,xx2) ; IF(xmin_edgs(i) > xmax) cycle
363 yy1=xrem(2,j)
364 yy2=xrem(9,j)
365 ymax_edgs(i)=max(yy1,yy2) ; IF(ymax_edgs(i) < ymin) cycle
366 ymin_edgs(i)=min(yy1,yy2) ; IF(ymin_edgs(i) > ymax) cycle
367 zz1=xrem(3,j)
368 zz2=xrem(10,j)
369 zmax_edgs(i)=max(zz1,zz2) ; IF(zmax_edgs(i) < zmin) cycle
370 zmin_edgs(i)=min(zz1,zz2) ; IF(zmin_edgs(i) > zmax) cycle
371
372 !-------------------------------------------!
373 ! VOXEL OCCUPIED BY THE EDGE !
374 !-------------------------------------------!
375 !Voxel_lower_left_bound for this edge
376 ix1=int(nbx*(xmin_edgs(i)-xminb)/(xmaxb-xminb))
377 iy1=int(nby*(ymin_edgs(i)-yminb)/(ymaxb-yminb))
378 iz1=int(nbz*(zmin_edgs(i)-zminb)/(zmaxb-zminb))
379 ix1=max(1,2+min(nbx,ix1))
380 iy1=max(1,2+min(nby,iy1))
381 iz1=max(1,2+min(nbz,iz1))
382 !Voxel_upper_right_bound for this edge
383 ix2=int(nbx*(xmax_edgs(i)-xminb)/(xmaxb-xminb))
384 iy2=int(nby*(ymax_edgs(i)-yminb)/(ymaxb-yminb))
385 iz2=int(nbz*(zmax_edgs(i)-zminb)/(zmaxb-zminb))
386 ix2=max(1,2+min(nbx,ix2))
387 iy2=max(1,2+min(nby,iy2))
388 iz2=max(1,2+min(nbz,iz2))
389
390 !pour reset des voxel
391 min_ix = min(min_ix,ix1)
392 min_iy = min(min_iy,iy1)
393 min_iz = min(min_iz,iz1)
394 max_ix = max(max_ix,ix2)
395 max_iy = max(max_iy,iy2)
396 max_iz = max(max_iz,iz2)
397
398 !----------------------------------------------!
399 ! EDGE STORAGE FOR EACH VOXEL (CHAINED ARRAY) !
400 !----------------------------------------------!
401C
402C VOXEL(i,j,k) LCHAIN_LAST(FIRST)
403C +-----------+------------+
404C | =>FIRST | =>LAST |
405C +--+--------+--+---------+
406C | |
407C | |
408C | |
409C | | LCHAIN_ELEM(*) LCHAIN_NEXT(*)
410C | | +------------+-----------+
411C +-------------->| edge_id | iadd 3 | 1:FIRST --+
412C | +------------+-----------+ |
413C | | | | 2 |
414C | +------------+-----------+ |
415C | | edge_id | iadd 4 | 3 <-------+
416C | +------------+-----------+ |
417C | | edge_id | iadd 6 | 4 <-------+
418C | +------------+-----------+ |
419C | | | | 5 |
420C | +------------+-----------+ |
421C +-->| edge_id | 0 | 6:LAST <--+
422C +------------+-----------+
423C | | | MAX_ADD
424C +------------+-----------+
425C
426 !Pour tous les voxels qu'occupe la brique
427 DO iz = iz1,iz2
428 DO iy = iy1,iy2
429 DO ix = ix1,ix2
430 first_add = voxel(ix,iy,iz)
431 IF(first_add == 0)THEN
432 !voxel encore vide
433 voxel(ix,iy,iz) = current_add ! adresse dans le tableau chaine de la premiere eddge trouvee occupant le voxel
434 lchain_last(current_add) = current_add ! dernier=adresse pour l edge courante
435 lchain_elem(current_add) = i ! edge ID
436 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste !
437 ELSE
438 !voxel contenant deja une edge
439 prev_add = lchain_last(first_add) ! devient l'avant-dernier
440 lchain_last(first_add) = current_add ! maj du dernier
441 lchain_elem(current_add) = i ! edge ID
442 lchain_next(prev_add) = current_add ! maj du suivant 0 -> CURRENT_ADD
443 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste
444 ENDIF
445 current_add = current_add+1
446 IF( current_add>=max_add)THEN
447 !OPTIMISATION : suprresion du deallocate/GOTO debut.
448 !REALLOCATE SI PAS ASSEZ DE PLACE : inutile de recommencer de 1 a MAX_ADD-1, on poursuit de MAX_ADD a 2*MAX_ADD
449 max_add = 2 * max_add
450 !print *, "reallocate remote"
454 ENDIF
455 ENDDO !IX
456 ENDDO !IY
457 ENDDO !IZ
458
459 ENDDO !DO NRTS+1,NRTS+NRTSR
460
461
462 END IF !(ITASK==0)
463
464 CALL my_barrier !le tableau voxel doit etre rempli avant de continuer
465 !max_add doit etre le meme pour tout le monde.
466C=======================================================================
467C 3 A partir des voxels occupes par une edge main, on est en
468C mesure de connaitre toutes les edges escalves dans ce voisinage.
469C Ce qui permet de creer des couples cancidats pour le contact
470C Si la penetration est positive.
471C=======================================================================
472 nedg = 0
473 i_stok = 0
474 marge = tzinf - max(gap+dgapload,drad)
475
476 DO iedg=1,nrtm
477
478 IF(stfm(iedg) == zero)cycle ! on ne retient pas les facettes detruites
479
480c AAA = ZERO !MARGE
481 aaa = tzinf
482 IF(igap == 0)THEN
483 aaa = tzinf
484 ELSE
485 aaa = marge+
486 . max(max(gapmin,bgapsmx+gap_m(iedg))+dgapload,drad)
487 ENDIF
488
489
490 !-------------------------------------------!
491 ! (N1,N2) is the current main edge !
492 !-------------------------------------------!
493 n1 = irectm(1,iedg)
494 n2 = irectm(2,iedg)
495 mm1 = itab(n1)
496 mm2 = itab(n2)
497
498 !-------------------------------------------!
499 ! X-coordinates of the four nodes !
500 !-------------------------------------------!
501 xx1=x(1,n1)
502 xx2=x(1,n2)
503 yy1=x(2,n1)
504 yy2=x(2,n2)
505 zz1=x(3,n1)
506 zz2=x(3,n2)
507
508 xmax_edgm(iedg)=max(xx1,xx2) ! +TZINF
509 xmin_edgm(iedg)=min(xx1,xx2) ! -TZINF
510 ymax_edgm(iedg)=max(yy1,yy2) ! +TZINF
511 ymin_edgm(iedg)=min(yy1,yy2) ! -TZINF
512 zmax_edgm(iedg)=max(zz1,zz2) ! +TZINF
513 zmin_edgm(iedg)=min(zz1,zz2) ! -TZINF
514
515 !-------------------------------------------!
516 ! VOXEL OCCUPIED BY THE BRICK !
517 !-------------------------------------------!
518 !Voxel_lower_left_bound for this element---+
519 ix1=int(nbx*(xmin_edgm(iedg)-aaa-xminb)/(xmaxb-xminb))
520 iy1=int(nby*(ymin_edgm(iedg)-aaa-yminb)/(ymaxb-yminb))
521 iz1=int(nbz*(zmin_edgm(iedg)-aaa-zminb)/(zmaxb-zminb))
522 ix1=max(1,2+min(nbx,ix1))
523 iy1=max(1,2+min(nby,iy1))
524 iz1=max(1,2+min(nbz,iz1))
525 !Voxel_upper_right_bound for this element---+
526 ix2=int(nbx*(xmax_edgm(iedg)+aaa-xminb)/(xmaxb-xminb))
527 iy2=int(nby*(ymax_edgm(iedg)+aaa-yminb)/(ymaxb-yminb))
528 iz2=int(nbz*(zmax_edgm(iedg)+aaa-zminb)/(zmaxb-zminb))
529 ix2=max(1,2+min(nbx,ix2))
530 iy2=max(1,2+min(nby,iy2))
531 iz2=max(1,2+min(nbz,iz2))
532
533 deja = 0 ! l edge n est pas encore candidate.
534 i_stok_bak = i_stok
535C
536C--- IREMGAP - tag of deactivated lines
537 IF(flagremnode==2)THEN
538 k = kremnode(2*(iedg-1)+1)
539 l = kremnode(2*(iedg-1)+2)-1
540 DO m=k,l
541 tagremline(remnode(m)) = 1
542 ENDDO
543 ENDIF
544C
545 !ON PARCOURS A NOUVEAU LES EDGES SECND DANS LE VOISINAGE DE L EDGE main IEDG
546 !ON CONSTITUE ICI ONE COUPLE
547
548 DO iz = iz1,iz2
549 DO iy = iy1,iy2
550 DO ix = ix1,ix2
551
552 chain_add = voxel(ix,iy,iz) ! adresse dans le tableau chaine de la premiere edge stoquee dans le voxel
553 DO WHILE(chain_add /= 0) ! BOUCLE SUR LES EDGES DU VOXEL COURANT
554 i = lchain_elem(chain_add) ! numeros des edge_id balayes dans le voxel courant
555
556 !secnd edge nodes, exclure couples avec noeud commun
557 IF (i<=nrts)THEN
558 ss1=itab(irects(1,i))
559 ss2=itab(irects(2,i))
560 ELSE
561 ss1=irem(2,i-nrts)
562 ss2=irem(3,i-nrts)
563 END IF
564
565 IF( (ss1==mm1).OR.(ss1==mm2).OR.
566 . (ss2==mm1).OR.(ss2==mm2) )THEN
567 chain_add = lchain_next(chain_add)
568 cycle
569 END IF
570
571 !unicite des couples
572 IF(iauto==1 .AND. mm1<ss1 )THEN
573 chain_add = lchain_next(chain_add)
574 cycle
575 END IF
576
577C IREMPGAP
578 IF (flagremnode == 2) THEN
579 IF (i <= nrts) THEN
580C- Local Taged lines are removed
581 IF(tagremline(i)==1) THEN
582 chain_add = lchain_next(chain_add)
583 cycle
584 ENDIF
585 ELSE
586C- Remote lines are identified by nodes
587 k = kremnode(2*(iedg-1)+2)
588 l = kremnode(2*(iedg-1)+3)-1
589 remove_remote = 0
590 DO m=k,l,2
591 IF ((ss1==remnode(m)).AND.(ss2==remnode(m+1))) remove_remote = 1
592 ENDDO
593 IF (remove_remote==1) THEN
594 chain_add = lchain_next(chain_add)
595 cycle
596 ENDIF
597 ENDIF
598 ENDIF
599
600 i_stok = i_stok + 1 !on dispose d'un candidat
601 prov_s(i_stok) = i !edge secnd
602 prov_m(i_stok) = iedg !edge main
603
604 !print *, "candidat:", IEDG, I
605 IF(deja==0) nedg = nedg + 1 !nombre d edges candidate au calcul de contact (debug)
606 deja=1 !l edge main IEDG fait l'objet d'une ecriture de candidat. On compte les edges main faisant l'objet dun couple candidate : on ne doit plus incrementer NEDG pour les autres edge secnd testees.
607 chain_add = lchain_next(chain_add)
608C-----------------------------------------------------
609 IF(i_stok>=nvsiz)THEN
610 CALL i11sto_vox(
611 1 nvsiz ,irects,irectm,x ,ii_stok,
612 2 cand_s,cand_m,nsn4 ,noint ,marge,
613 3 i_mem ,prov_s,prov_m,eshift,addcm ,
614 4 chaine,nrts ,itab ,ifpen ,iform,
615 5 gapmin,drad ,igap, gap_s, gap_m,
616 7 gap_s_l, gap_m_l ,dgapload)
617
618 IF(i_mem==2) THEN
619 !print *, "too much candidates"
620c IF (ITASK==0)II_STOK=ZERO
621 GOTO 1000
622 END if!(I_MEM==2)
623 i_stok = i_stok-nvsiz
624C !DIR$ ASSUME (I_STOK < NVSIZ)
625 DO j=1,i_stok
626 prov_s(j) = prov_s(j+nvsiz)
627 prov_m(j) = prov_m(j+nvsiz)
628 ENDDO
629 ENDIF
630C-----------------------------------------------------
631
632 ENDDO !NEXT WHILE(CHAIN_ADD /= 0)
633 ENDDO !NEXT IZ
634 ENDDO !NEXT IY
635 ENDDO !NEXT IZ
636
637C--- IREMGAP - clean of tagremline
638 IF(flagremnode==2)THEN
639 k = kremnode(2*(iedg-1)+1)
640 l = kremnode(2*(iedg-1)+2)-1
641 DO m=k,l
642 tagremline(remnode(m)) = 0
643 ENDDO
644 ENDIF
645
646 ENDDO !NEXT IEDG
647
648C-------------------------------------------------------------------------
649C FIN DU TRI
650C-------------------------------------------------------------------------
651
652 IF(i_stok/=0)CALL i11sto_vox(
653 1 i_stok,irects,irectm,x ,ii_stok,
654 2 cand_s,cand_m,nsn4 ,noint ,marge ,
655 3 i_mem ,prov_s,prov_m,eshift,addcm ,
656 4 chaine,nrts ,itab ,ifpen ,iform ,
657 5 gapmin,drad ,igap, gap_s ,gap_m ,
658 7 gap_s_l, gap_m_l ,dgapload)
659
660
661C=======================================================================
662C 4 remise a zero des noeuds dans les boites et desallocation
663C=======================================================================
664
665 1000 CONTINUE
666
667 CALL my_barrier !ne pas desalloue tant que les autres threads n'ont pas fini de travailler
668
669
670 ! peut etre otpimisee : ne pas reinitialiser tout le tableau (plusieurs solutions possibles)
671 !! VOXEL( MIN_IX:MAX_IX, MIN_IY:MAX_IY, MIN_IZ:MAX_IZ ) = 0
672 tmin(1) = min_ix
673 tmin(2) = min_iy
674 tmin(3) = min_iz
675
676 tmax(1) = max_ix
677 tmax(2) = max_iy
678 tmax(3) = max_iz
679
680 IF (itask==0)THEN
681 !RESET VOXEL WITHIN USED RANGE ONLY
682 DO k= tmin(3),tmax(3)
683 DO j= tmin(2),tmax(2)
684 DO i= tmin(1),tmax(1)
685 voxel(i,j,k) = 0
686 END DO
687 END DO
688 END DO
689 !CHAINED LIST DEALLOCATION
690 DEALLOCATE(lchain_next)
691 DEALLOCATE(lchain_elem)
692 DEALLOCATE(lchain_last)
693 IF(flagremnode==2) DEALLOCATE(tagremline)
694 END IF
695
696C___________________________________________________________________________________________________________
697
698 RETURN
699 END
700
701
702
#define my_real
Definition cppsort.cpp:32
if(complex_arithmetic) id
subroutine i11sto_vox(j_stok, irects, irectm, x, ii_stok, cand_s, cand_m, nsn4, noint, marge, i_mem, prov_s, prov_m, eshift, addcm, chaine, nrts, itab, ifpen, iform, gapmin, drad, igap, gap_s, gap_m, gap_s_l, gap_m_l, dgapload)
Definition i11sto.F:39
subroutine i11trivox(irects, irectm, x, nrtm, nrtsr, xyzm, ii_stok, cand_s, cand_m, nsn4, noint, tzinf, i_mem, eshift, addcm, chaine, nrts, itab, stfs, stfm, iauto, voxel, nbx, nby, nbz, itask, ifpen, iform, gapmin, drad, marge, gap_s, gap_m, gap_s_l, gap_m_l, bgapsmx, igap, gap, flagremnode, kremnode, remnode, dgapload)
Definition i11trivox.F:46
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
integer function, dimension(:), pointer ireallocate(ptr, new_size)
Definition realloc_mod.F:39
integer, dimension(:), pointer lchain_elem
Definition tri11_mod.F:37
integer max_iz
Definition tri11_mod.F:33
integer min_ix
Definition tri11_mod.F:33
integer, dimension(:), pointer lchain_last
Definition tri11_mod.F:39
integer min_iz
Definition tri11_mod.F:33
integer min_iy
Definition tri11_mod.F:33
integer, dimension(:), pointer lchain_next
Definition tri11_mod.F:38
integer max_iy
Definition tri11_mod.F:33
integer max_ix
Definition tri11_mod.F:33
integer, dimension(:,:), allocatable irem
Definition tri7box.F:339
subroutine my_barrier
Definition machine.F:31