OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i25trivox_edg.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!|| i25trivox_edg ../starter/source/interfaces/inter3d1/i25trivox_edg.F
25!||--- called by ------------------------------------------------------
26!|| i25buce_edg ../starter/source/interfaces/inter3d1/i25buce_edg.F
27!||--- calls -----------------------------------------------------
28!|| bitget ../starter/source/interfaces/inter3d1/bitget.F
29!|| i25sto_e2s ../starter/source/interfaces/inter3d1/i25sto_e2s.F
30!|| i25sto_edg ../starter/source/interfaces/inter3d1/i25sto_edg.f
31!|| ireallocate ../starter/share/modules1/realloc_mod.F
32!||--- uses -----------------------------------------------------
33!|| realloc_mod ../starter/share/modules1/realloc_mod.F
34!|| tri11 ../starter/share/modules1/tri11_mod.F
35!|| tri7box ../starter/share/modules1/tri7box.F
36!||====================================================================
37 SUBROUTINE i25trivox_edg(
38 1 I_MEM ,VMAXDT ,INACTI ,
39 2 IRECT ,X ,STF ,XYZM ,
40 3 II_STOK ,CANDS_E2E ,ESHIFT ,NEDGE_T ,CANDM_E2E ,
41 4 MULNSNE,NOINT ,BGAPEMX ,SSHIFT ,NRTM_T ,
42 5 VOXEL ,NBX ,NBY ,NBZ ,
43 6 IGAP ,GAP_M ,GAP_M_L ,DRAD ,MARGE ,
44 7 ITASK ,ITAB ,LL_STOK ,MULNSNS ,
45 8 MBINFLG ,EBINFLG,ILEV ,CAND_A ,CAND_P ,
46 9 FLAGREMNODE,KREMNODE,REMNODE ,S_REMNODE_EDG,
47 A IEDGE ,NEDGE ,LEDGE ,MSEGTYP ,IGAP0 ,
48 B ADMSR,EDG_BISECTOR,VTX_BISECTOR,
49 C CANDM_E2S,CANDS_E2S,CAND_B,CAND_PS,GAPE ,
50 D GAP_E_L ,DGAPLOAD,FLAG_REMOVED_NODE,
51 E S_KREMNODE_E2S,S_REMNODE_E2S,KREMNODE_E2S,REMNODE_E2S,
52 F S_KREMNODE_EDG)
53C============================================================================
54C M o d u l e s
55C-----------------------------------------------
56 USE realloc_mod
57 USE tri7box
58 USE tri11
59C-----------------------------------------------
60C I m p l i c i t T y p e s
61C-----------------------------------------------
62#include "implicit_f.inc"
63C-----------------------------------------------
64C G l o b a l P a r a m e t e r s
65C-----------------------------------------------
66#include "mvsiz_p.inc"
67c parameter setting the size for the vector (orig version is 128)
68 INTEGER NVECSZ
69 PARAMETER (NVECSZ = mvsiz)
70C-----------------------------------------------
71C C o m m o n B l o c k s
72C-----------------------------------------------
73#include "param_c.inc"
74C-----------------------------------------------
75C ROLE DE LA ROUTINE:
76C ===================
77C CLASSE LES EDGES DANS DES BOITES
78C RECHERCHE POUR CHAQUE FACETTE DES BOITES CONCERNES
79C RECHERCHE DES CANDIDATS
80C-----------------------------------------------
81C D u m m y A r g u m e n t s
82C-----------------------------------------------
83 INTEGER, INTENT(in) :: S_REMNODE_EDG
84 LOGICAL, INTENT(in) :: FLAG_REMOVED_NODE !< flag to remove some S node from the list of candidates
85 INTEGER, INTENT(in) :: S_KREMNODE_E2S !< size of KREMNODE_E2S array
86 INTEGER, INTENT(in) :: S_REMNODE_E2S !< size of REMNODE_E2S array
87 INTEGER, INTENT(in) :: S_KREMNODE_EDG
88
89 INTEGER I_MEM(2),INACTI,ITASK,IGAP,IEDGE,NEDGE,ESHIFT,NEDGE_T,SSHIFT,NRTM_T,IGAP0,
90 . MULNSNE,MULNSNS,NOINT,NBX,NBY,NBZ,
91 . CANDS_E2E(*),CANDM_E2E(*),
92 . IRECT(4,*), VOXEL(NBX+2,NBY+2,NBZ+2),II_STOK,LL_STOK,ITAB(*),
93 . MBINFLG(*),EBINFLG(*),ILEV,CAND_A(*),LEDGE(NLEDGE,*),ADMSR(4,*),MSEGTYP(*),
94 . CANDM_E2S(*),CANDS_E2S(*),CAND_B(*),
95 . FLAGREMNODE,KREMNODE(S_KREMNODE_EDG),REMNODE(*)
96 INTEGER, DIMENSION(S_KREMNODE_E2S), INTENT(in) :: KREMNODE_E2S !< adress of forbidden S edge
97 INTEGER, DIMENSION(S_REMNODE_E2S), INTENT(in) :: REMNODE_E2S !< list of forbidden S edge
98C REAL
99 my_real
100 . x(3,*),xyzm(6),stf(*),gap_m(*),gap_m_l(*),gape(*),gap_e_l(*),cand_p(*), cand_ps(*),
101 . marge,bgapemx,vmaxdt,drad
102 my_real , INTENT(IN) :: dgapload
103 real*4 edg_bisector(3,4,*), vtx_bisector(3,2,*)
104C-----------------------------------------------
105C L o c a l V a r i a b l e s
106C-----------------------------------------------
107 INTEGER I,J,I_STOK, SOL_EDGE, SH_EDGE,
108 . N1,N2,NN,NE,K,L,J_STOK,II,JJ,NA,NB,
109 . PROV_S(MVSIZ),PROV_M(MVSIZ),
110 . M,NS1,NS2,NSE,NS,SIZE,Z_FIRST,Z_LAST
111C REAL
112 my_real
113 . dx,dy,dz,xs,ys,zs,xx,sx,sy,sz,s2,
114 . xmin, xmax,ymin, ymax,zmin, zmax, tz, gapsmx, gapl,
115 . xx1,xx2,xx3,xx4,yy1,yy2,yy3,yy4,zz1,zz2,zz3,zz4,
116 . d1x,d1y,d1z,d2x,d2y,d2z,dd1,dd2,d2,a2,gs
117c provisoire
118 INTEGER TAGEDG(NEDGE)
119 INTEGER IX,IY,IZ,IEDG,
120 . M1, M2, M3, M4, MM1,MM2,MM3,MM4,SS1,SS2,
121 . IMS1,IMS2,ISS1,ISS2,
122 . am1,am2,as1,as2,
123 . ix1,iy1,iz1,ix2,iy2,iz2
124 INTEGER, DIMENSION(3) :: TMIN,TMAX
125 my_real
126 . XMINB,YMINB,ZMINB,XMAXB,YMAXB,ZMAXB,AAA,DRAD2,
127 . xmax_edgs, xmin_edgs, !cotes min/max des aretes seconds et mains
128 . ymax_edgs, ymin_edgs,
129 . zmax_edgs, zmin_edgs,
130 . xmax_edgm, xmin_edgm,
131 . ymax_edgm, ymin_edgm,
132 . zmax_edgm, zmin_edgm
133 INTEGER FIRST_ADD, PREV_ADD, CHAIN_ADD, CURRENT_ADD, MAX_ADD
134 INTEGER BITGET
135 INTEGER WORK(70000)
136 INTEGER, DIMENSION(:,:),ALLOCATABLE :: IIXYZ,LEDGE_TMP
137 INTEGER, DIMENSION(:),ALLOCATABLE :: INDEX
138 EXTERNAL BITGET
139C-----------------------------------------------
140C
141 INTEGER, DIMENSION(:), ALLOCATABLE :: TAGREMLINE
142C-----------------------------------------------
143C
144 DRAD2 = zero
145C--- IREMGAP - array for tag of deactivated lines
146 IF(flag_removed_node)THEN
147 ALLOCATE(tagremline(nedge))
148 tagremline(1:nedge) = 0
149 ENDIF
150
151 sol_edge =iedge/10 ! solids
152 sh_edge =iedge-10*sol_edge ! shells
153
154 min_ix=nbx+2
155 min_iy=nby+2
156 min_iz=nbz+2
157 max_ix=1
158 max_iy=1
159 max_iz=1
160
161 !---------------------------------------------------------!
162 ! Allocation des tableaux chaines !
163 !---------------------------------------------------------!
164 max_add = max(1,4*nedge)
165 ALLOCATE(lchain_elem(1:max_add))
166 ALLOCATE(lchain_next(1:max_add))
167 ALLOCATE(lchain_last(1:max_add))
168
169 ALLOCATE(iixyz(4,nedge),ledge_tmp(nledge,nedge),index(2*nedge))
170C Barrier to wait init voxel and allocation
171C CALL MY_BARRIER
172C
173 xmin = xyzm(1)
174 ymin = xyzm(2)
175 zmin = xyzm(3)
176 xmax = xyzm(4)
177 ymax = xyzm(5)
178 zmax = xyzm(6)
179
180c dev future: xminb plus grand que xmin...
181 xminb = xmin
182 yminb = ymin
183 zminb = zmin
184 xmaxb = xmax
185 ymaxb = ymax
186 zmaxb = zmax
187
188C=======================================================================
189C 1 mise des edges dans les boites
190C=======================================================================
191 IF(itask == 0)THEN
192
193 iixyz(1:4,1:nedge)=0
194 DO i=1,nedge
195 index(i) =i
196 END DO
197
198 DO i=1,nedge
199
200 ne =ledge(1,i)
201 IF(stf(ne)==0) cycle ! on ne retient pas les facettes detruites
202
203C IF(LEDGE(3,I)/=0) CYCLE ! on ne retient que les aretes de bord
204
205 n1 = ledge(5,i)
206 n2 = ledge(6,i)
207
208 !-------------------------------------------!
209 xx1=x(1,n1)
210 xx2=x(1,n2)
211 yy1=x(2,n1)
212 yy2=x(2,n2)
213 zz1=x(3,n1)
214 zz2=x(3,n2)
215 IF(igap0 == 0)THEN
216 xmax_edgs=max(xx1,xx2); IF(xmax_edgs < xmin) cycle
217 xmin_edgs=min(xx1,xx2); IF(xmin_edgs > xmax) cycle
218 ymax_edgs=max(yy1,yy2); IF(ymax_edgs < ymin) cycle
219 ymin_edgs=min(yy1,yy2); IF(ymin_edgs > ymax) cycle
220 zmax_edgs=max(zz1,zz2); IF(zmax_edgs < zmin) cycle
221 zmin_edgs=min(zz1,zz2); IF(zmin_edgs > zmax) cycle
222 ELSE
223 xmax_edgs=max(xx1,xx2)+gap_m(ne); IF(xmax_edgs < xmin) cycle
224 xmin_edgs=min(xx1,xx2)-gap_m(ne); IF(xmin_edgs > xmax) cycle
225 ymax_edgs=max(yy1,yy2)+gap_m(ne); IF(ymax_edgs < ymin) cycle
226 ymin_edgs=min(yy1,yy2)-gap_m(ne); IF(ymin_edgs > ymax) cycle
227 zmax_edgs=max(zz1,zz2)+gap_m(ne); IF(zmax_edgs < zmin) cycle
228 zmin_edgs=min(zz1,zz2)-gap_m(ne); IF(zmin_edgs > zmax) cycle
229 END IF
230
231 !-------------------------------------------!
232 ! VOXEL OCCUPIED BY THE EDGE !
233 !-------------------------------------------!
234 !Voxel_lower_left_bound for this edge
235 ix1=int(nbx*(xmin_edgs-xminb)/(xmaxb-xminb))
236 iy1=int(nby*(ymin_edgs-yminb)/(ymaxb-yminb))
237 iz1=int(nbz*(zmin_edgs-zminb)/(zmaxb-zminb))
238 ix1=max(1,2+min(nbx,ix1))
239 iy1=max(1,2+min(nby,iy1))
240 iz1=max(1,2+min(nbz,iz1))
241 !Voxel_upper_right_bound for this edge
242 ix2=int(nbx*(xmax_edgs-xminb)/(xmaxb-xminb))
243 iy2=int(nby*(ymax_edgs-yminb)/(ymaxb-yminb))
244 iz2=int(nbz*(zmax_edgs-zminb)/(zmaxb-zminb))
245 ix2=max(1,2+min(nbx,ix2))
246 iy2=max(1,2+min(nby,iy2))
247 iz2=max(1,2+min(nbz,iz2))
248
249 iixyz(1,i)=min(n1,n2)
250 iixyz(2,i)=0
251 iixyz(3,i)=0
252 iixyz(4,i)=0
253
254 index(i) =i
255
256 END DO
257
258 CALL my_orders(0,work,iixyz,index,nedge,4)
259 ledge_tmp(1:nledge,1:nedge)=ledge(1:nledge,1:nedge)
260 DO i=1,nedge
261 k=index(i)
262 ledge(1:nledge,i)=ledge_tmp(1:nledge,k)
263c write(500,'(3I10)') I,LEDGE(5,I),LEDGE(6,I)
264 END DO
265 DEALLOCATE(iixyz,index,ledge_tmp)
266
267 current_add=1 ! premiere adresse
268
269 DO i=1,nedge
270
271 ne =ledge(1,i)
272 IF(stf(ne)==zero) cycle ! on ne retient pas les facettes detruites
273
274 IF(ledge(7,i) < 0) cycle ! larete nest pas une arete second
275
276 n1 = ledge(5,i)
277 n2 = ledge(6,i)
278
279 !-------------------------------------------!
280 xx1=x(1,n1)
281 xx2=x(1,n2)
282 yy1=x(2,n1)
283 yy2=x(2,n2)
284 zz1=x(3,n1)
285 zz2=x(3,n2)
286 IF(igap0 == 0)THEN
287 xmax_edgs=max(xx1,xx2); IF(xmax_edgs < xmin) cycle
288 xmin_edgs=min(xx1,xx2); IF(xmin_edgs > xmax) cycle
289 ymax_edgs=max(yy1,yy2); IF(ymax_edgs < ymin) cycle
290 ymin_edgs=min(yy1,yy2); IF(ymin_edgs > ymax) cycle
291 zmax_edgs=max(zz1,zz2); IF(zmax_edgs < zmin) cycle
292 zmin_edgs=min(zz1,zz2); IF(zmin_edgs > zmax) cycle
293 ELSE
294 xmax_edgs=max(xx1,xx2)+gape(i); IF(xmax_edgs < xmin) cycle
295 xmin_edgs=min(xx1,xx2)-gape(i); IF(xmin_edgs > xmax) cycle
296 ymax_edgs=max(yy1,yy2)+gape(i); IF(ymax_edgs < ymin) cycle
297 ymin_edgs=min(yy1,yy2)-gape(i); IF(ymin_edgs > ymax) cycle
298 zmax_edgs=max(zz1,zz2)+gape(i); IF(zmax_edgs < zmin) cycle
299 zmin_edgs=min(zz1,zz2)-gape(i); IF(zmin_edgs > zmax) cycle
300 END IF
301
302 !-------------------------------------------!
303 ! VOXEL OCCUPIED BY THE EDGE !
304 !-------------------------------------------!
305 !Voxel_lower_left_bound for this edge
306 ix1=int(nbx*(xmin_edgs-xminb)/(xmaxb-xminb))
307 iy1=int(nby*(ymin_edgs-yminb)/(ymaxb-yminb))
308 iz1=int(nbz*(zmin_edgs-zminb)/(zmaxb-zminb))
309 ix1=max(1,2+min(nbx,ix1))
310 iy1=max(1,2+min(nby,iy1))
311 iz1=max(1,2+min(nbz,iz1))
312 !Voxel_upper_right_bound for this edge
313 ix2=int(nbx*(xmax_edgs-xminb)/(xmaxb-xminb))
314 iy2=int(nby*(ymax_edgs-yminb)/(ymaxb-yminb))
315 iz2=int(nbz*(zmax_edgs-zminb)/(zmaxb-zminb))
316 ix2=max(1,2+min(nbx,ix2))
317 iy2=max(1,2+min(nby,iy2))
318 iz2=max(1,2+min(nbz,iz2))
319
320 !pour reset des voxel
321 min_ix = min(min_ix,ix1)
322 min_iy = min(min_iy,iy1)
323 min_iz = min(min_iz,iz1)
324 max_ix = max(max_ix,ix2)
325 max_iy = max(max_iy,iy2)
326 max_iz = max(max_iz,iz2)
327
328 !----------------------------------------------!
329 ! EDGE STORAGE FOR EACH VOXEL (CHAINED ARRAY) !
330 !----------------------------------------------!
331C
332C VOXEL(i,j,k) LCHAIN_LAST(FIRST)
333C +-----------+------------+
334C | =>FIRST | =>LAST |
335C +--+--------+--+---------+
336C | |
337C | |
338C | |
339C | | LCHAIN_ELEM(*) LCHAIN_NEXT(*)
340C | | +------------+-----------+
341C +-------------->| edge_id | iadd 3 | 1:FIRST --+
342C | +------------+-----------+ |
343C | | | | 2 |
344C | +------------+-----------+ |
345C | | edge_id | iadd 4 | 3 <-------+
346C | +------------+-----------+ |
347C | | edge_id | iadd 6 | 4 <-------+
348C | +------------+-----------+ |
349C | | | | 5 |
350C | +------------+-----------+ |
351C +-->| edge_id | 0 | 6:LAST <--+
352C +------------+-----------+
353C | | | MAX_ADD
354C +------------+-----------+
355 DO iz = iz1,iz2
356 DO iy = iy1,iy2
357 DO ix = ix1,ix2
358
359 first_add = voxel(ix,iy,iz)
360
361 IF(first_add == 0)THEN
362 !voxel encore vide
363 voxel(ix,iy,iz) = current_add ! adresse dans le tableau chaine de la premiere eddge trouvee occupant le voxel
364 lchain_last(current_add) = current_add ! dernier=adresse pour l edge courante
365 lchain_elem(current_add) = i ! edge ID
366 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste !
367 ELSE
368 !voxel contenant deja une edge
369 prev_add = lchain_last(first_add) ! devient l'avant-dernier
370 lchain_last(first_add) = current_add ! maj du dernier
371 lchain_elem(current_add) = i ! edge ID
372 lchain_next(prev_add) = current_add ! maj du suivant 0 -> CURRENT_ADD
373 lchain_next(current_add) = 0 ! pas de suivant car dernier de la liste
374 ENDIF
375
376 current_add = current_add+1
377
378 IF( current_add>=max_add)THEN
379 !OPTIMISATION : suprresion du deallocate/GOTO debut.
380 !REALLOCATE SI PAS ASSEZ DE PLACE : inutile de recommencer de 1 a MAX_ADD-1, on poursuit de MAX_ADD a 2*MAX_ADD
381 max_add = 2 * max_add
382 !print *, "reallocate"
383 lchain_next => ireallocate(lchain_next, max_add)
384 lchain_elem => ireallocate(lchain_elem, max_add)
385 lchain_last => ireallocate(lchain_last, max_add)
386 ENDIF
387
388 ENDDO !IX
389 ENDDO !IY
390 ENDDO !IZ
391
392 ENDDO
393
394 END IF
395
396 tagedg(1:nedge)=0
397C Barrier to wait task0 treatment
398C CALL MY_BARRIER
399C=======================================================================
400C Sorting vs main shell edges
401C=======================================================================
402 IF(sh_edge==0) GOTO 300
403C=======================================================================
404C 3 A partir des voxels occupes par une edge main, on est en
405C mesure de connaitre toutes les edges escalves dans ce voisinage.
406C Ce qui permet de creer des couples cancidats pour le contact
407C Si la penetration est positive.
408C=======================================================================
409
410 j_stok = 0
411
412 DO i=1,nedge_t
413
414 iedg=eshift+i
415
416 ne=ledge(1,iedg)
417 IF(stf(ne)==zero) cycle ! on ne retient pas les facettes detruites
418
419 IF(iabs(ledge(7,iedg))==1) cycle ! Main solid edge
420
421 !-------------------------------------------!
422 ! (N1,N2) is the current main edge !
423 !-------------------------------------------!
424
425 aaa = marge+bgapemx+gape(iedg)+dgapload
426
427 n1 = ledge(5,iedg)
428 n2 = ledge(6,iedg)
429 mm1 = itab(n1)
430 mm2 = itab(n2)
431 am1 = min(mm1,mm2)
432 am2 = max(mm1,mm2)
433
434 IF(ilev==2)THEN
435 ims1 = bitget(ebinflg(iedg),0)
436 ims2 = bitget(ebinflg(iedg),1)
437 END IF
438
439 !-------------------------------------------!
440 ! X-coordinates of the four nodes !
441 !-------------------------------------------!
442
443 xx1=x(1,n1)
444 xx2=x(1,n2)
445 yy1=x(2,n1)
446 yy2=x(2,n2)
447 zz1=x(3,n1)
448 zz2=x(3,n2)
449 xmax_edgm=max(xx1,xx2)+gape(iedg) ! +TZINF
450 xmin_edgm=min(xx1,xx2)-gape(iedg) ! -TZINF
451 ymax_edgm=max(yy1,yy2)+gape(iedg) ! +TZINF
452 ymin_edgm=min(yy1,yy2)-gape(iedg) ! -TZINF
453 zmax_edgm=max(zz1,zz2)+gape(iedg) ! +TZINF
454 zmin_edgm=min(zz1,zz2)-gape(iedg) ! -TZINF
455 !-------------------------------------------!
456 ! VOXEL OCCUPIED BY THE BRICK !
457 !-------------------------------------------!
458 !Voxel_lower_left_bound for this element---+
459 ix1=int(nbx*(xmin_edgm-aaa-xminb)/(xmaxb-xminb))
460 iy1=int(nby*(ymin_edgm-aaa-yminb)/(ymaxb-yminb))
461 iz1=int(nbz*(zmin_edgm-aaa-zminb)/(zmaxb-zminb))
462 ix1=max(1,2+min(nbx,ix1))
463 iy1=max(1,2+min(nby,iy1))
464 iz1=max(1,2+min(nbz,iz1))
465 !Voxel_upper_right_bound for this element---+
466 ix2=int(nbx*(xmax_edgm+aaa-xminb)/(xmaxb-xminb))
467 iy2=int(nby*(ymax_edgm+aaa-yminb)/(ymaxb-yminb))
468 iz2=int(nbz*(zmax_edgm+aaa-zminb)/(zmaxb-zminb))
469 ix2=max(1,2+min(nbx,ix2))
470 iy2=max(1,2+min(nby,iy2))
471 iz2=max(1,2+min(nbz,iz2))
472
473C--- IREMGAP - tag of deactivated lines
474 IF(flag_removed_node .AND. s_kremnode_edg > 0)THEN
475 k = kremnode(iedg)
476 l = kremnode(iedg+1)-1
477 DO m=k,l
478 tagremline(remnode(m)) = 1
479 ENDDO
480 ENDIF
481
482 DO iz = iz1,iz2
483 DO iy = iy1,iy2
484 DO ix = ix1,ix2
485
486 chain_add = voxel(ix,iy,iz) ! adresse dans le tableau chaine de la premiere edge stoquee dans le voxel
487 DO WHILE(chain_add /= 0) ! BOUCLE SUR LES EDGES DU VOXEL COURANT
488 jj = lchain_elem(chain_add) ! numeros des edge_id balayes dans le voxel courant
489
490 IF(tagedg(jj)/=0)THEN ! edge deja traitee vs cette arete main
491 chain_add = lchain_next(chain_add)
492 cycle
493 END IF
494
495 !secnd edge nodes, exclure couples avec noeud commun
496 IF (jj<=nedge)THEN
497 ss1= itab(ledge(5,jj))
498 ss2= itab(ledge(6,jj))
499 END IF
500
501 IF( (ss1==mm1).OR.(ss1==mm2).OR.
502 . (ss2==mm1).OR.(ss2==mm2) )THEN
503 chain_add = lchain_next(chain_add)
504 cycle
505 END IF
506
507 IF(ilev==2)THEN
508 IF(jj <= nedge) THEN
509 iss1=bitget(ebinflg(jj),0)
510 iss2=bitget(ebinflg(jj),1)
511 ENDIF
512 IF(.NOT.((ims1 == 1 .and. iss2==1).or.(ims2 == 1 .and. iss1==1)))THEN
513 chain_add = lchain_next(chain_add)
514 cycle
515 ENDIF
516 ENDIF
517
518
519 IF(iabs(ledge(7,iedg))/=1 .AND. ledge(7,jj)/=1)THEN
520 ! attention les traitements de i25dst3e pour les solides
521 ! ne sont pas symetriques main second.
522 as1 = min(ss1,ss2)
523 as2 = max(ss1,ss2)
524 ! unicite des couples
525 IF(am1 < as1 .OR. (am1 == as1 .AND. am2 < as2))THEN
526 chain_add = lchain_next(chain_add)
527 cycle
528 ENDIF
529 ENDIF
530
531C IREMPGAP - Tagged lines are removed
532 IF(flag_removed_node)THEN
533 IF(tagremline(jj)==1) THEN
534 chain_add = lchain_next(chain_add)
535 cycle
536 ENDIF
537 ENDIF
538
539 tagedg(jj)=1
540
541 j_stok = j_stok + 1 !on dispose d'un candidat
542 prov_s(j_stok) = jj !edge secnd
543 prov_m(j_stok) = iedg !edge main
544
545c IF(DEJA==0) NEDG = NEDG + 1 !nombre d edges candidate au calcul de contact (debug)
546c 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.
547 chain_add = lchain_next(chain_add)
548C-----------------------------------------------------
549 IF(j_stok==nvsiz)THEN
550 CALL i25sto_edg(
551 1 nvsiz ,irect ,x ,ii_stok,inacti,
552 2 cands_e2e,candm_e2e ,mulnsne,noint ,marge ,
553 3 i_mem(1) ,prov_s ,prov_m ,igap0,cand_a,
554 4 nedge ,ledge ,itab ,drad2 ,igap ,
555 5 gape ,gap_e_l,admsr ,edg_bisector,vtx_bisector ,
556 6 cand_p,dgapload)
557 IF(i_mem(1)/=0) GOTO 300
558 j_stok = 0
559 ENDIF
560C-----------------------------------------------------
561
562 ENDDO ! WHILE(CHAIN_ADD /= 0)
563
564 ENDDO !NEXT IZ
565 ENDDO !NEXT IY
566 ENDDO !NEXT IZ
567
568C Reset TAGEDG
569 DO iz = iz1,iz2
570 DO iy = iy1,iy2
571 DO ix = ix1,ix2
572
573 chain_add = voxel(ix,iy,iz)
574 DO WHILE(chain_add /= 0) ! BOUCLE SUR LES EDGES DU VOXEL COURANT
575
576 jj = lchain_elem(chain_add) ! numeros des edge_id balayes dans le voxel courant
577 tagedg(jj)=0
578
579 chain_add = lchain_next(chain_add)
580
581 END DO
582
583 ENDDO !NEXT IZ
584 ENDDO !NEXT IY
585 ENDDO !NEXT IZ
586
587C--- IREMGAP - clean of tagremline
588 IF(flag_removed_node .AND. s_kremnode_edg > 0)THEN
589 k = kremnode(iedg)
590 l = kremnode(iedg+1)-1
591 DO m=k,l
592 tagremline(remnode(m)) = 0
593 ENDDO
594 ENDIF
595
596
597 ENDDO !NEXT IEDG
598
599
600C-------------------------------------------------------------------------
601C FIN DU TRI vs Main shell edges
602C-------------------------------------------------------------------------
603 IF(j_stok/=0)CALL i25sto_edg(
604 1 j_stok ,irect ,x ,ii_stok,inacti,
605 2 cands_e2e,candm_e2e ,mulnsne,noint ,marge ,
606 3 i_mem(1) ,prov_s ,prov_m ,igap0,cand_a,
607 4 nedge ,ledge ,itab ,drad2 ,igap ,
608 5 gape ,gap_e_l,admsr ,edg_bisector,vtx_bisector ,
609 6 cand_p,dgapload)
610
611 300 CONTINUE
612C=======================================================================
613C Sorting vs main solid edges
614C=======================================================================
615 IF(sol_edge==0) GOTO 400
616C=======================================================================
617C 3bis A partir des voxels occupes par une edge main, on est en
618C mesure de connaitre toutes les edges escalves dans ce voisinage.
619C Ce qui permet de creer des couples cancidats pour le contact
620C Si la penetration est positive.
621C=======================================================================
622
623 j_stok = 0
624
625 DO i=1,nrtm_t
626
627 ne =sshift+i
628
629 IF(msegtyp(ne)/=0) cycle ! not a solid edge
630 IF(stf(ne)==zero) cycle ! on ne retient pas les facettes detruites
631
632 m1 = irect(1,ne)
633 m2 = irect(2,ne)
634 m3 = irect(3,ne)
635 m4 = irect(4,ne)
636 mm1= itab(m1)
637 mm2= itab(m2)
638 mm3= itab(m3)
639 mm4= itab(m4)
640
641 xx1=x(1,m1)
642 yy1=x(2,m1)
643 zz1=x(3,m1)
644 xx2=x(1,m2)
645 yy2=x(2,m2)
646 zz2=x(3,m2)
647 xx3=x(1,m3)
648 yy3=x(2,m3)
649 zz3=x(3,m3)
650 xx4=x(1,m4)
651 yy4=x(2,m4)
652 zz4=x(3,m4)
653
654 xmax_edgm=max(xx1,xx2,xx3,xx4) ! +TZINF
655 xmin_edgm=min(xx1,xx2,xx3,xx4) ! -TZINF
656 ymax_edgm=max(yy1,yy2,yy3,yy4) ! +TZINF
657 ymin_edgm=min(yy1,yy2,yy3,yy4) ! -TZINF
658 zmax_edgm=max(zz1,zz2,zz3,zz4) ! +TZINF
659 zmin_edgm=min(zz1,zz2,zz3,zz4) ! -TZINF
660
661 dx=em02*(xmax_edgm-xmin_edgm)
662 dy=em02*(ymax_edgm-ymin_edgm)
663 dz=em02*(zmax_edgm-zmin_edgm)
664
665 xmax_edgm=xmax_edgm+dx
666 xmin_edgm=xmin_edgm-dx
667 ymax_edgm=ymax_edgm+dy
668 ymin_edgm=ymin_edgm-dy
669 zmax_edgm=zmax_edgm+dz
670 zmin_edgm=zmin_edgm-dz
671
672 aaa = marge+bgapemx +dgapload ! filtrer vs GAPE(JJ) dans i25pen3_edg !
673
674 !-------------------------------------------!
675 ! VOXEL OCCUPIED BY THE BRICK !
676 !-------------------------------------------!
677 !Voxel_lower_left_bound for this element---+
678 ix1=int(nbx*(xmin_edgm-aaa-xminb)/(xmaxb-xminb))
679 iy1=int(nby*(ymin_edgm-aaa-yminb)/(ymaxb-yminb))
680 iz1=int(nbz*(zmin_edgm-aaa-zminb)/(zmaxb-zminb))
681 ix1=max(1,2+min(nbx,ix1))
682 iy1=max(1,2+min(nby,iy1))
683 iz1=max(1,2+min(nbz,iz1))
684 !Voxel_upper_right_bound for this element---+
685 ix2=int(nbx*(xmax_edgm+aaa-xminb)/(xmaxb-xminb))
686 iy2=int(nby*(ymax_edgm+aaa-yminb)/(ymaxb-yminb))
687 iz2=int(nbz*(zmax_edgm+aaa-zminb)/(zmaxb-zminb))
688 ix2=max(1,2+min(nbx,ix2))
689 iy2=max(1,2+min(nby,iy2))
690 iz2=max(1,2+min(nbz,iz2))
691
692 IF(ilev==2)THEN
693 ims1 = bitget(mbinflg(ne),0)
694 ims2 = bitget(mbinflg(ne),1)
695 END IF
696
697 IF(flag_removed_node .AND. s_remnode_e2s > 0)THEN
698 k = kremnode_e2s(2*(ne-1)+1)
699 l = kremnode_e2s(2*(ne-1)+2)-1
700 DO m=k,l
701 tagremline(remnode_e2s(m)) = 1
702 ENDDO
703 ENDIF
704
705 DO iz = iz1,iz2
706 DO iy = iy1,iy2
707 DO ix = ix1,ix2
708
709 chain_add = voxel(ix,iy,iz) ! adresse dans le tableau chaine de la premiere edge stoquee dans le voxel
710 DO WHILE(chain_add /= 0) ! BOUCLE SUR LES EDGES DU VOXEL COURANT
711 jj = lchain_elem(chain_add) ! numeros des edge_id balayes dans le voxel courant
712
713 IF(tagedg(jj)/=0)THEN ! edge deja traitee vs cette arete main
714 chain_add = lchain_next(chain_add)
715 cycle
716 END IF
717
718 !secnd edge nodes, exclure couples avec noeud commun
719 IF (jj<=nedge)THEN
720 ss1= itab(ledge(5,jj))
721 ss2= itab(ledge(6,jj))
722 END IF
723
724 IF((ss1==mm1).OR.(ss1==mm2).OR.(ss1==mm3).OR.(ss1==mm4).OR.
725 . (ss2==mm1).OR.(ss2==mm2).OR.(ss2==mm3).OR.(ss2==mm4))THEN
726 chain_add = lchain_next(chain_add)
727 cycle
728 END IF
729
730 IF(ilev==2)THEN
731 IF(jj <= nedge) THEN
732 iss1=bitget(ebinflg(jj),0)
733 iss2=bitget(ebinflg(jj),1)
734 ENDIF
735
736 IF(.NOT.((ims1 == 1 .and. iss2==1).or.(ims2 == 1 .and. iss1==1)))THEN
737 chain_add = lchain_next(chain_add)
738 cycle
739 ENDIF
740 ENDIF
741
742 IF (flag_removed_node) THEN
743 IF(tagremline(jj)==1) THEN
744 chain_add = lchain_next(chain_add)
745 cycle
746 ENDIF
747 ENDIF
748 tagedg(jj)=1
749
750
751C-----------------------------------------------------
752
753 j_stok = j_stok + 1 !on dispose d'un candidat
754 prov_s(j_stok) = jj !edge secnd
755 prov_m(j_stok) = ne !segment main
756
757C-----------------------------------------------------
758 IF(j_stok==nvsiz)THEN
759 CALL i25sto_e2s(
760 1 nvsiz ,irect ,x ,ll_stok,inacti,
761 2 cands_e2s,candm_e2s,mulnsns,noint ,marge ,
762 3 i_mem(2) ,prov_s ,prov_m ,igap0,cand_b,
763 4 nedge ,ledge ,itab ,drad2 ,igap ,
764 5 gap_m ,gap_m_l,gape ,gap_e_l,admsr ,
765 6 edg_bisector,vtx_bisector ,cand_ps,dgapload)
766 IF(i_mem(2)/=0) GOTO 400
767 j_stok = 0
768 ENDIF
769C-----------------------------------------------------
770
771 chain_add = lchain_next(chain_add) ! Next RTM
772
773 ENDDO ! WHILE(CHAIN_ADD /= 0)
774
775 ENDDO !NEXT IZ
776 ENDDO !NEXT IY
777 ENDDO !NEXT IZ
778
779C Reset TAGEDG
780 DO iz = iz1,iz2
781 DO iy = iy1,iy2
782 DO ix = ix1,ix2
783
784 chain_add = voxel(ix,iy,iz)
785 DO WHILE(chain_add /= 0) ! BOUCLE SUR LES EDGES DU VOXEL COURANT
786
787 jj = lchain_elem(chain_add) ! numeros des edge_id balayes dans le voxel courant
788 tagedg(jj)=0
789
790 chain_add = lchain_next(chain_add)
791
792 END DO
793
794 ENDDO !NEXT IZ
795 ENDDO !NEXT IY
796 ENDDO !NEXT IZ
797
798 IF(flag_removed_node.AND.s_remnode_e2s>0)THEN
799 k = kremnode_e2s(2*(ne-1)+1)
800 l = kremnode_e2s(2*(ne-1)+2)-1
801 DO m=k,l
802 tagremline(remnode_e2s(m)) = 0
803 ENDDO
804 ENDIF
805
806 ENDDO !NEXT IEDG
807C-------------------------------------------------------------------------
808C FIN DU TRI vs Main solid edges
809C-------------------------------------------------------------------------
810 IF(j_stok/=0)CALL i25sto_e2s(
811 1 j_stok ,irect ,x ,ll_stok,inacti,
812 2 cands_e2s,candm_e2s,mulnsns,noint ,marge ,
813 3 i_mem(2) ,prov_s ,prov_m ,igap0,cand_b,
814 4 nedge ,ledge ,itab ,drad2 ,igap ,
815 5 gap_m ,gap_m_l,gape ,gap_e_l,admsr ,
816 6 edg_bisector,vtx_bisector ,cand_ps,dgapload)
817
818C=======================================================================
819C 4 remise a zero des noeuds dans les boites
820C=======================================================================
821 400 CONTINUE
822
823C Barrier to avoid reinitialization before end of sorting
824C CALL MY_BARRIER
825
826 tmin(1) = min_ix
827 tmin(2) = min_iy
828 tmin(3) = min_iz
829
830 tmax(1) = max_ix
831 tmax(2) = max_iy
832 tmax(3) = max_iz
833
834 IF (itask==0)THEN
835 !RESET VOXEL WITHIN USED RANGE ONLY
836 DO k= tmin(3),tmax(3)
837 DO j= tmin(2),tmax(2)
838 DO i= tmin(1),tmax(1)
839 voxel(i,j,k) = 0
840 END DO
841 END DO
842 END DO
843 !CHAINED LIST DEALLOCATION
844 DEALLOCATE(lchain_next)
845 DEALLOCATE(lchain_elem)
846 DEALLOCATE(lchain_last)
847 IF(flag_removed_node) DEALLOCATE(tagremline)
848 ENDIF
849C=======================================================================
850
851 RETURN
852 END
my_real
#define my_real
Definition cppsort.cpp:32
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
my_orders
void my_orders(int *mode, int *iwork, int *data, int *index, int *n, int *irecl)
Definition my_orders.c:82
realloc_mod
Definition realloc_mod.F:29
realloc_mod::ireallocate
integer function, dimension(:), pointer ireallocate(ptr, new_size)
Definition realloc_mod.F:39
tri11::lchain_elem
integer, dimension(:), pointer lchain_elem
Definition tri11_mod.F:37
tri11::max_iz
integer max_iz
Definition tri11_mod.F:33
tri11::min_ix
integer min_ix
Definition tri11_mod.F:33
tri11::lchain_last
integer, dimension(:), pointer lchain_last
Definition tri11_mod.F:39
tri11::min_iz
integer min_iz
Definition tri11_mod.F:33
tri11::min_iy
integer min_iy
Definition tri11_mod.F:33
tri11::lchain_next
integer, dimension(:), pointer lchain_next
Definition tri11_mod.F:38
tri11::max_iy
integer max_iy
Definition tri11_mod.F:33
tri11::max_ix
integer max_ix
Definition tri11_mod.F:33
i25sto_e2s
subroutine i25sto_e2s(j_stok, irect, x, ii_stok, inacti, cand_s, cand_m, mulnsn, noint, marge, i_mem, prov_s, prov_m, igap0, cand_a, nedge, ledge, itab, drad, igap, gap_m, gap_m_l, gape, gap_e_l, admsr, edg_bisector, vtx_bisector, cand_p, dgapload)
Definition i25sto_e2s.F:39
i25sto_edg
subroutine i25sto_edg(j_stok, irect, x, ii_stok, inacti, cand_s, cand_m, mulnsn, noint, marge, i_mem, prov_s, prov_m, igap0, cand_a, nedge, ledge, itab, drad, igap, gape, gap_e_l, admsr, edg_bisector, vtx_bisector, cand_p, dgapload)
Definition i25sto_edg.F:39
i25trivox_edg
subroutine i25trivox_edg(i_mem, vmaxdt, inacti, irect, x, stf, xyzm, ii_stok, cands_e2e, eshift, nedge_t, candm_e2e, mulnsne, noint, bgapemx, sshift, nrtm_t, voxel, nbx, nby, nbz, igap, gap_m, gap_m_l, drad, marge, itask, itab, ll_stok, mulnsns, mbinflg, ebinflg, ilev, cand_a, cand_p, flagremnode, kremnode, remnode, s_remnode_edg, iedge, nedge, ledge, msegtyp, igap0, admsr, edg_bisector, vtx_bisector, candm_e2s, cands_e2s, cand_b, cand_ps, gape, gap_e_l, dgapload, flag_removed_node, s_kremnode_e2s, s_remnode_e2s, kremnode_e2s, remnode_e2s, s_kremnode_edg)
Definition i25trivox_edg.F:53
starter
program starter
Definition starter.F:39