OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
comput_coinknot.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!|| comput_coinknot ../starter/source/elements/ige3d/comput_coinknot.F
25!||--- called by ------------------------------------------------------
26!|| prerafig3d ../starter/source/elements/ige3d/prerafig3d.F
27!||--- calls -----------------------------------------------------
28!||====================================================================
29 SUBROUTINE comput_coinknot(IEL ,IXIG3D,KXIG3D ,MESHIGE,PTANG1 ,PTANG2 ,IDDIR,
30 . IDTANG1,IDTANG2,NELDIR,NELTANG1,NELTANG2,DIR,
31 . TAB_COINKNOT,L_TAB_COINKNOT,TAB_ELCUT,L_TAB_ELCUT,TAB_NEWEL,L_TAB_NEWEL,
32 . KNOT,IAD_KNOT,NKNOT1,NKNOT2,NKNOT3,IDFILS,KNOTLOCEL,NEWKNOT,IPARTIG3D,
33 . TAB_OLDIDCUT,IDCUT,FLAG)
34C----------------------------------------------------------------------
35C routine that defines the shapes of the meshsurf according to the elements
36C neighbors to refine in the cutting plane
37C this routine also fills all the elementary data of the
38C elements and the new elements (kxig3d + ixig3d)
39C----------------------------------------------------------------------
40C-----------------------------------------------
41C I m p l i c i t T y p e s
42C-----------------------------------------------
43#include "implicit_f.inc"
44C-----------------------------------------------
45C C o m m o n B l o c k s
46C-----------------------------------------------
47#include "param_c.inc"
48#include "ige3d_c.inc"
49C-----------------------------------------------
50C D u m m y A r g u m e n t s
51C-----------------------------------------------
52 INTEGER IXIG3D(*),KXIG3D(NIXIG3D,*),MESHIGE(NELTANG1,NELTANG2,NELDIR),
53 . IDDIR(*),IDTANG1(*),IDTANG2(*),TAB_ELCUT(L_TAB_ELCUT),
54 . TAB_NEWEL(L_TAB_NEWEL),IDFILS(NBFILSMAX,*),IPARTIG3D(*),TAB_OLDIDCUT(3,*)
55 INTEGER IDCUT,IDCUT_VOISIN,IEL,DIR,PTANG1,PTANG2,NELDIR,NELTANG1,NELTANG2,
56 . l_tab_elcut,l_tab_newel,l_tab_coinknot,flag,iad_knot,
57 . nknot1,nknot2,nknot3
58 my_real tab_coinknot(2,*),knot(*)
59 my_real knotlocel(2,3,*),newknot
60C-----------------------------------------------
61C L o c a l V a r i a b l e s
62C-----------------------------------------------
63 INTEGER I,J,K,L,M,IDMESHIGE,NBSEGMENTS,DIRTANG1,DIRTANG2
64 INTEGER IDNBCUT,ID1KNOT1,ID2KNOT1,ID1KNOT2,ID2KNOT2, IDNEXTG,
65 . IDNEXTD, IDNEXTEL,IAD_IXIG3D
66 INTEGER COINDEP(2),COIN(2),DIRECTION,WORK(70000),
67 . l_tabworkel, coinenglob(2,2),nb_newfils
68 INTEGER L_TABWORK, NEXT, REFNBCUT, IOUT, INTERSEC
69 INTEGER, DIMENSION(:), ALLOCATABLE :: INDEX, TABWORKEL, TABWORKELTRI
70 INTEGER OFFSET_KNOT,OFFSET_TANG1,OFFSET_TANG2,BORDHAUT,BORDDROIT
71 my_real COIN_TMP(2,20), DET, T1, T2, XA(5),YA(5),
72 . xb, yb, xc, yc, xd, yd, tol
73C=======================================================================
74C----------------------------------------------------------
75C ELEMENT ISO GEOMETRIQUE
76C-----------------------------------------------
77C KXIG3D(6,*) : index of 1st knot in the Xknot vector corresponding to the element
78C KXIG3D(7,*) : index of 1st knot in the Yknot vector corresponding to the element
79C KXIG3D(8,*) : index of 1st knot in the Zknot vector corresponding to the element
80C KXIG3D(9,*) : index of 2nd knot in the Xknot vector corresponding to the element
81C KXIG3D(10,*) : index of 2nd knot in the Yknot vector corresponding to the element
82C KXIG3D(11,*) : index of 2nd knot in the Zknot vector corresponding to the element
83C KXIG3D(12,*) : number of element's cuts needed in the X direction
84C KXIG3D(13,*) : number of element's cuts needed in the Y direction
85C KXIG3D(14,*) : number of element's cuts needed in the Z direction
86C=======================================================================
87C
88C in the parametric space
89C
90C RAFFINEMENT X : RAFFINEMENT Y : RAFFINEMENT Z :
91C
92C Z X Y
93C | | |
94C |_ _ Y |_ _ Z |_ _ X
95C / / /
96C X Y Z
97C
98C=======================================================================
99C
100C TANG2 _______
101C | | |
102C |_ _ TANG1 | |
103C / |_______|
104C DIR COIN(1,2) /
105C
106C=======================================================================
107C----------------------------------------------------------------------
108C initialization of the variables according to the cutting direction
109C----------------------------------------------------------------------
110C
111 IF(dir==1) THEN
112 dirtang1 = 2
113 dirtang2 = 3
114 idnbcut=12
115 id1knot1=7 ! index of 1st knot in the Yknot vector corresponding to the element
116 id2knot1=10 ! index of 2nd knot in the Yknot vector corresponding to the element
117 id1knot2=8 ! index of 1st knot in the Zknot vector corresponding to the element
118 id2knot2=11 ! index of 2nd knot in the Zknot vector corresponding to the element
119 idmeshige=iddir(kxig3d(6,iel))
120 offset_knot = iad_knot
121 offset_tang1 = iad_knot + nknot1
122 offset_tang2 = iad_knot + nknot1 + nknot2
123 ELSEIF(dir==2) THEN
124 dirtang1 = 3
125 dirtang2 = 1
126 idnbcut=13
127 id1knot1=8 ! index of 1st knot in the Zknot vector corresponding to the element
128 id2knot1=11 ! index of 2nd knot in the Zknot vector corresponding to the element
129 id1knot2=6 ! index of 1st knot in the Xknot vector corresponding to the element
130 id2knot2=9 ! index of 2nd knot in the Xknot vector corresponding to the element
131 idmeshige=iddir(kxig3d(7,iel))
132 offset_knot = iad_knot + nknot1
133 offset_tang1 = iad_knot + nknot1 + nknot2
134 offset_tang2 = iad_knot
135 ELSEIF(dir==3) THEN
136 dirtang1 = 1
137 dirtang2 = 2
138 idnbcut=14
139 id1knot1=6 ! index of 1st knot in the Xknot vector corresponding to the element
140 id2knot1=9 ! index of 2nd knot in the xknot vector corresponding to the element
141 id1knot2=7 ! index of 1st knot in the Yknot vector corresponding to the element
142 id2knot2=10 ! index of 2nd knot in the Yknot vector corresponding to the element
143 idmeshige=iddir(kxig3d(8,iel))
144 offset_knot = iad_knot + nknot1 + nknot2
145 offset_tang1 = iad_knot
146 offset_tang2 = iad_knot + nknot1
147 ENDIF
148
149 l_tabwork = 50000
150 tol = em06
151
152C----------------------------------------------------------------------
153C FIRST STEP : ON PART DE L'ELEMENT INITIAL A RAFFINER
154C----------------------------------------------------------------------
155
156 coindep(1)=idtang1(kxig3d(id1knot1,iel))
157 coindep(2)=idtang2(kxig3d(id1knot2,iel))
158
159 coin(:)=coindep(:)
160
161 coinenglob(1,1)= min(coin(1),10000)
162 coinenglob(2,1)= min(coin(2),10000)
163 coinenglob(1,2)= max(coin(1),0)
164 coinenglob(2,2)= max(coin(2),0)
165
166 nbsegments = 0
167 l_tab_coinknot = 0
168 l_tab_elcut = 0
169 l_tab_newel = 0
170 l_tabworkel = 0
171
172 ALLOCATE(tabworkel(l_tabwork))
173 ALLOCATE(tabworkeltri(l_tabwork))
174 tabworkeltri(:) = 0
175 tabworkel(:) = ep06
176 coin_tmp(:,:) = 0
177
178 refnbcut = kxig3d(idnbcut,iel)
179
180c
181c PREMIER COIN
182c
183 l_tab_coinknot = l_tab_coinknot + 1
184
185 coinenglob(1,1)= min(coin(1),coinenglob(1,1))
186 coinenglob(2,1)= min(coin(2),coinenglob(2,1))
187 coinenglob(1,2)= max(coin(1),coinenglob(1,2))
188 coinenglob(2,2)= max(coin(2),coinenglob(2,2))
189 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
190 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
191 IF(flag==1) THEN
192 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
193 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
194 ENDIF
195c
196 borddroit = 0
197 bordhaut = 0
198
199C----------------------------------------------------------------------
200C always start from the bottom left cell, so we search upwards
201C----------------------------------------------------------------------
202c
203 direction = 1
204
205 DO WHILE (direction==1)
206
207c extension of the segment upwards
208
209 idnextd=meshige(coin(1),coin(2),idmeshige)
210 IF(kxig3d(idnbcut,idnextd)==1) THEN
211 idcut_voisin = 0
212 ELSE
213 idcut_voisin = tab_oldidcut(dir,idnextd)-kxig3d(idnbcut,idnextd)+1
214 ENDIF
215 IF(idcut_voisin/=idcut) THEN
216 direction=2 ! We have to turn right
217 cycle
218 ELSE
219 l_tabworkel = l_tabworkel + 1
220 tabworkel(l_tabworkel) = idnextd
221 ENDIF
222
223 IF(coin(2)<neltang2) THEN
224 IF(coin(1)>1) THEN
225 idnextg=meshige(coin(1)-1,coin(2),idmeshige)
226 IF(kxig3d(idnbcut,idnextg)==1) THEN
227 idcut_voisin = 0
228 ELSE
229 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
230 ENDIF
231 IF(idcut_voisin==idcut) THEN
232 direction=4 ! We have to turn left
233 cycle
234 ELSE
235 coin(2) = coin(2)+1 ! on monte
236 ENDIF
237 ELSE
238 coin(2) = coin(2)+1 ! on monte
239 ENDIF
240 ELSE
241 IF(coin(1)>1) THEN
242 idnextg=meshige(coin(1)-1,coin(2),idmeshige)
243 IF(kxig3d(idnbcut,idnextg)==1) THEN
244 idcut_voisin = 0
245 ELSE
246 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
247 ENDIF
248 IF(idcut_voisin==idcut) THEN
249 direction=4 ! We have to turn left
250 cycle
251 ELSE
252 bordhaut=1
253 direction=2 ! We have to turn right
254 cycle
255 ENDIF
256 ELSE
257 bordhaut=1
258 direction=2 ! We have to turn right
259 cycle
260 ENDIF
261 ENDIF
262
263 ENDDO
264c
265c NOUVEAU COIN
266c
267 nbsegments = nbsegments + 1
268 l_tab_coinknot = l_tab_coinknot + 1
269
270 coinenglob(1,1)= min(coin(1),coinenglob(1,1))
271 coinenglob(2,1)= min(coin(2),coinenglob(2,1))
272 coinenglob(1,2)= max(coin(1),coinenglob(1,2))
273 coinenglob(2,2)= max(coin(2),coinenglob(2,2))
274 IF(bordhaut==0) THEN
275 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
276 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
277 IF(flag==1) THEN
278 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
279 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
280 ENDIF
281 ELSE
282 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
283 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
284 IF(flag==1) THEN
285 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
286 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
287 ENDIF
288 ENDIF
289C
290 DO WHILE (.NOT.(coin(1)==coindep(1).AND.coin(2)==coindep(2)).OR.l_tab_coinknot<5)
291
292 SELECT CASE (direction)
293
294 CASE(1)
295
296 DO WHILE (direction==1)
297
298c extension of the segment upwards
299
300 idnextd=meshige(coin(1),coin(2),idmeshige)
301 IF(kxig3d(idnbcut,idnextd)==1) THEN
302 idcut_voisin = 0
303 ELSE
304 idcut_voisin = tab_oldidcut(dir,idnextd)-kxig3d(idnbcut,idnextd)+1
305 ENDIF
306 IF(idcut_voisin==idcut) THEN
307 l_tabworkel = l_tabworkel + 1
308 tabworkel(l_tabworkel) = idnextd
309 ELSE
310 direction=2 ! We have to turn right
311 cycle
312 ENDIF
313
314 IF(coin(2)<neltang2) THEN
315
316 IF(coin(1)>1) THEN
317 idnextg=meshige(coin(1)-1,coin(2),idmeshige)
318 IF(kxig3d(idnbcut,idnextg)==1) THEN
319 idcut_voisin = 0
320 ELSE
321 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
322 ENDIF
323 IF(idcut_voisin==idcut) THEN
324 direction=4 ! We have to turn left
325 cycle
326 ELSE
327 coin(2) = coin(2)+1
328 ENDIF
329 ELSE
330 coin(2) = coin(2)+1
331 ENDIF
332 ELSE
333 IF(coin(1)>1) THEN
334 idnextg=meshige(coin(1)-1,coin(2),idmeshige)
335 IF(kxig3d(idnbcut,idnextg)==1) THEN
336 idcut_voisin = 0
337 ELSE
338 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
339 ENDIF
340 IF(idcut_voisin==idcut) THEN
341 direction=4 ! We have to turn left
342 cycle
343 ELSE
344 bordhaut=1
345 direction=2 ! We have to turn right
346 cycle
347 ENDIF
348 ELSE
349 bordhaut=1
350 direction=2 ! We have to turn right
351 cycle
352 ENDIF
353 ENDIF
354
355 ENDDO
356c
357c NOUVEAU COIN
358c
359 nbsegments = nbsegments + 1
360 l_tab_coinknot = l_tab_coinknot + 1
361
362 IF(bordhaut==0) THEN
363 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
364 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
365 IF(flag==1) THEN
366 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
367 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
368 ENDIF
369 ELSE
370 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
371 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
372 IF(flag==1) THEN
373 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
374 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
375 ENDIF
376 ENDIF
377c
378 CASE(2)
379
380 DO WHILE (direction==2)
381
382c extension of the segment to the right
383
384 IF(bordhaut==1) THEN ! If we are at the top edge
385 idnextg=meshige(coin(1),coin(2),idmeshige)
386 IF(kxig3d(idnbcut,idnextg)==1) THEN
387 idcut_voisin = 0
388 ELSE
389 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
390 ENDIF
391 IF(idcut_voisin==idcut) THEN
392 l_tabworkel = l_tabworkel + 1
393 tabworkel(l_tabworkel) = idnextg
394 ELSE
395 direction=3 ! We have to turn downstairs
396 cycle
397 ENDIF
398 IF(coin(1)<neltang1) THEN ! If we are not on the right edge
399 coin(1)=coin(1)+1
400 cycle
401 ELSE
402 direction=3 ! We have to turn downstairs
403 borddroit=1
404 cycle
405 ENDIF
406 ENDIF
407
408 IF(bordhaut==0) THEN
409
410 IF(coin(2)>1) THEN
411 idnextd=meshige(coin(1),coin(2)-1,idmeshige)
412 IF(kxig3d(idnbcut,idnextd)==1) THEN
413 idcut_voisin = 0
414 ELSE
415 idcut_voisin = tab_oldidcut(dir,idnextd)-kxig3d(idnbcut,idnextd)+1
416 ENDIF
417 IF(idcut_voisin==idcut) THEN
418 l_tabworkel = l_tabworkel + 1
419 tabworkel(l_tabworkel) = idnextd
420 ELSE
421 direction=3 ! We have to turn downstairs
422 cycle
423 ENDIF
424
425 IF(coin(1)<neltang1) THEN ! If we are not on the right edge
426 idnextg=meshige(coin(1),coin(2),idmeshige)
427 IF(kxig3d(idnbcut,idnextg)==1) THEN
428 idcut_voisin = 0
429 ELSE
430 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
431 ENDIF
432 IF(idcut_voisin==idcut) THEN
433 l_tabworkel = l_tabworkel + 1
434 tabworkel(l_tabworkel) = idnextg
435 direction=1 ! We have to turn at the top
436 cycle
437 ELSE
438 coin(1)=coin(1)+1 ! On avance
439 cycle
440 ENDIF
441 ELSE ! If we are on the right edge
442 idnextg=meshige(coin(1),coin(2),idmeshige)
443 IF(kxig3d(idnbcut,idnextg)==1) THEN
444 idcut_voisin = 0
445 ELSE
446 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
447 ENDIF
448 IF(idcut_voisin==idcut) THEN
449 l_tabworkel = l_tabworkel + 1
450 tabworkel(l_tabworkel) = idnextg
451 direction=1 ! We have to turn downstairs
452 cycle
453 ELSE
454 direction=3 ! We have to turn downstairs
455 borddroit=1
456 cycle
457 ENDIF
458 ENDIF
459
460 ENDIF
461
462 ENDIF
463
464 ENDDO
465c
466c NOUVEAU COIN
467c
468 nbsegments = nbsegments + 1
469 l_tab_coinknot = l_tab_coinknot + 1
470
471 IF(borddroit==0) THEN
472 IF(bordhaut==0) THEN
473 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
474 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
475 IF(flag==1) THEN
476 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
477 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
478 ENDIF
479 ELSE
480 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
481 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
482 IF(flag==1) THEN
483 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
484 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
485 ENDIF
486 ENDIF
487 ELSE
488 IF(bordhaut==0) THEN
489 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
490 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
491 IF(flag==1) THEN
492 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
493 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
494 ENDIF
495 ELSE
496 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
497 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
498 IF(flag==1) THEN
499 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
500 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id2knot2,meshige(coin(1),coin(2),idmeshige)))
501 ENDIF
502 ENDIF
503 ENDIF
504c
505 CASE(3)
506
507 DO WHILE (direction==3)
508
509c extension of the segment downwards
510
511 IF(bordhaut==1) THEN
512 IF(borddroit==1) THEN
513 bordhaut=0
514 cycle
515 ELSE
516 IF(coin(2)>1) THEN
517 bordhaut=0
518 cycle
519 ELSE
520 direction=4 ! We have to turn left
521 bordhaut=0
522 cycle
523 ENDIF
524 ENDIF
525 ELSE
526
527 IF(borddroit==1) THEN
528 IF(coin(2)>1) THEN
529 idnextg=meshige(coin(1),coin(2)-1,idmeshige)
530 IF(kxig3d(idnbcut,idnextg)==1) THEN
531 idcut_voisin = 0
532 ELSE
533 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
534 ENDIF
535 IF(idcut_voisin==idcut) THEN
536 l_tabworkel = l_tabworkel + 1
537 tabworkel(l_tabworkel) = idnextg
538 coin(2)=coin(2)-1
539 cycle
540 ELSE
541 direction=4 ! We have to turn to the left
542 cycle
543 ENDIF
544 ELSE
545 direction=4
546 cycle
547 ENDIF
548 ENDIF
549
550 IF(coin(2)>1) THEN
551 IF(coin(1)>1) THEN
552 idnextd=meshige(coin(1)-1,coin(2)-1,idmeshige)
553 IF(kxig3d(idnbcut,idnextd)==1) THEN
554 idcut_voisin = 0
555 ELSE
556 idcut_voisin = tab_oldidcut(dir,idnextd)-kxig3d(idnbcut,idnextd)+1
557 ENDIF
558 IF(idcut_voisin==idcut) THEN
559 l_tabworkel = l_tabworkel + 1
560 tabworkel(l_tabworkel) = idnextd
561 ELSE
562 direction=4 ! We have to turn to the left
563 cycle
564 ENDIF
565 idnextg=meshige(coin(1),coin(2)-1,idmeshige)
566 IF(kxig3d(idnbcut,idnextg)==1) THEN
567 idcut_voisin = 0
568 ELSE
569 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
570 ENDIF
571 IF(idcut_voisin==idcut) THEN
572 IF(borddroit==0) THEN
573 direction=2 ! We have to turn to the right
574 cycle
575 ELSE
576 l_tabworkel = l_tabworkel + 1
577 tabworkel(l_tabworkel) = idnextg
578 coin(2) = coin(2)-1
579 bordhaut=0
580 ENDIF
581 ELSE
582 coin(2) = coin(2)-1
583 bordhaut=0
584 ENDIF
585 ELSE ! We are on the left
586 idnextg=meshige(coin(1),coin(2)-1,idmeshige)
587 IF(kxig3d(idnbcut,idnextg)==1) THEN
588 idcut_voisin = 0
589 ELSE
590 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
591 ENDIF
592 IF(idcut_voisin/=idcut) THEN
593 direction=2 ! We have to turn to the right
594 cycle
595 ELSE
596 coin(2) = coin(2)-1
597 bordhaut=0
598 ENDIF
599 ENDIF
600 ELSE ! We are at the bottom edge
601 direction=4 ! We have to turn to the left
602 ENDIF
603 ENDIF
604 ENDDO
605c
606c NOUVEAU COIN
607c
608 nbsegments = nbsegments + 1
609 l_tab_coinknot = l_tab_coinknot + 1
610
611 IF(borddroit==0) THEN
612 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
613 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
614 IF(flag==1) THEN
615 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
616 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
617 ENDIF
618 ELSE
619 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
620 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
621 IF(flag==1) THEN
622 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id2knot1,meshige(coin(1),coin(2),idmeshige)))
623 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
624 ENDIF
625 ENDIF
626c
627 CASE(4)
628
629 DO WHILE (direction==4)
630
631c extension of the segment to the left
632
633 IF(borddroit==1) THEN
634 IF(coin(1)>1) THEN
635 borddroit=0
636 cycle
637 ELSE
638 borddroit=0
639 direction=1
640 cycle
641 ENDIF
642 ENDIF
643
644 IF(coin(1)>1) THEN
645 idnextd=meshige(coin(1)-1,coin(2),idmeshige)
646 IF(kxig3d(idnbcut,idnextd)==1) THEN
647 idcut_voisin = 0
648 ELSE
649 idcut_voisin = tab_oldidcut(dir,idnextd)-kxig3d(idnbcut,idnextd)+1
650 ENDIF
651 IF(idcut_voisin==idcut) THEN
652 l_tabworkel = l_tabworkel + 1
653 tabworkel(l_tabworkel) = idnextd
654 ELSE
655 direction=1 ! We have to turn up
656 cycle
657 ENDIF
658
659 IF(coin(2)>1) THEN
660 idnextg=meshige(coin(1)-1,coin(2)-1,idmeshige)
661 IF(kxig3d(idnbcut,idnextg)==1) THEN
662 idcut_voisin = 0
663 ELSE
664 idcut_voisin = tab_oldidcut(dir,idnextg)-kxig3d(idnbcut,idnextg)+1
665 ENDIF
666 IF(idcut_voisin==idcut) THEN
667 direction=3 ! We have to turn down
668 cycle
669 ELSE
670 coin(1) = coin(1)-1
671 borddroit=0
672 ENDIF
673 ELSE
674 coin(1) = coin(1)-1
675 borddroit=0
676 ENDIF
677 ELSE
678 direction=1 ! We have to turn up
679 ENDIF
680
681 ENDDO
682c
683c NOUVEAU COIN
684c
685 nbsegments = nbsegments + 1
686 l_tab_coinknot = l_tab_coinknot + 1
687
688 coin_tmp(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
689 coin_tmp(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
690 IF(flag==1) THEN
691 tab_coinknot(1,l_tab_coinknot) = knot(offset_tang1+kxig3d(id1knot1,meshige(coin(1),coin(2),idmeshige)))
692 tab_coinknot(2,l_tab_coinknot) = knot(offset_tang2+kxig3d(id1knot2,meshige(coin(1),coin(2),idmeshige)))
693 ENDIF
694c
695 END SELECT
696
697 coinenglob(1,1)= min(coin(1),coinenglob(1,1))
698 coinenglob(2,1)= min(coin(2),coinenglob(2,1))
699 coinenglob(1,2)= max(coin(1),coinenglob(1,2))
700 coinenglob(2,2)= max(coin(2),coinenglob(2,2))
701
702 ENDDO
703
704C----------------------------------------------------------------------
705C processing of the elements of the defined contour and those inside
706C----------------------------------------------------------------------
707
708 ALLOCATE(index(2*l_tabwork))
709 CALL my_orders(0, work, tabworkel, index, l_tabwork , 1)
710
711 DO i=1,l_tabwork
712 tabworkeltri(i)=tabworkel(index(i))
713 ENDDO
714
715 DEALLOCATE(index)
716
717 IF(flag==0) THEN
718 DO i=1,l_tabwork
719 IF(tabworkeltri(i)==ep06) EXIT
720 IF(i/=1) THEN
721 IF(tabworkeltri(i-1)==tabworkeltri(i)) cycle
722 ENDIF
723c
724 nb_newfils=0
725c
726 IF(newknot>=knotlocel(1,dir,tabworkeltri(i)).AND.
727 . newknot<=knotlocel(2,dir,tabworkeltri(i))) THEN
728 l_tab_elcut = l_tab_elcut + 1
729 l_tab_newel = l_tab_newel + 1
730 ENDIF
731 DO j=1,idfils(1,tabworkeltri(i))-nb_newfils
732 IF(newknot>=knotlocel(1,dir,idfils(j+1,tabworkeltri(i))).AND.
733 . newknot<=knotlocel(2,dir,idfils(j+1,tabworkeltri(i)))) THEN
734 l_tab_elcut = l_tab_elcut + 1
735 l_tab_newel = l_tab_newel + 1
736 ENDIF
737 ENDDO
738c
739 ENDDO
740 ELSE
741 DO i=1,l_tabwork
742 IF(tabworkeltri(i)==ep06) EXIT
743 IF(i/=1) THEN
744 IF(tabworkeltri(i-1)==tabworkeltri(i)) cycle
745 ENDIF
746c
747 nb_newfils=0
748c
749 IF(newknot>=knotlocel(1,dir,tabworkeltri(i)).AND.
750 . newknot<=knotlocel(2,dir,tabworkeltri(i))) THEN
751 l_tab_elcut = l_tab_elcut + 1
752 tab_elcut(l_tab_elcut) = tabworkeltri(i)
753
754 addelig3d=addelig3d+1
755
756 l_tab_newel = l_tab_newel + 1
757 tab_newel(l_tab_newel) = numelig3d0+addelig3d
758
759 idfils(1,tabworkeltri(i))=idfils(1,tabworkeltri(i))+1
760 idfils(idfils(1,tabworkeltri(i))+1,tabworkeltri(i))=numelig3d0+addelig3d
761C counting the new children
762 nb_newfils=nb_newfils+1
763
764 kxig3d(:,numelig3d0+addelig3d) = kxig3d(:,tabworkeltri(i))
765 iad_ixig3d = sixig3d + addsixig3d + 1
766 kxig3d(4,numelig3d0+addelig3d) = iad_ixig3d
767 DO k=1,kxig3d(3,tabworkeltri(i))
768 ixig3d(iad_ixig3d+k-1) = ixig3d(kxig3d(4,tabworkeltri(i))+k-1)
769 ENDDO
770 addsixig3d = addsixig3d + kxig3d(3,tabworkeltri(i))
771
772 kxig3d(5,numelig3d0+addelig3d) = numelig3d0+addelig3d ! Do not prevent the double id
773 ipartig3d(numelig3d0+addelig3d) = ipartig3d(tabworkeltri(i))
774
775 kxig3d(15,numelig3d0+addelig3d) = inod_ige
776 inod_ige = inod_ige + 64
777
778 kxig3d(idnbcut,numelig3d0+addelig3d) = 1
779
780 knotlocel(:,:,numelig3d0+addelig3d) = knotlocel(:,:,tabworkeltri(i))
781 knotlocel(1,dir,numelig3d0+addelig3d) = newknot
782
783 knotlocel(2,dir,tabworkeltri(i)) = newknot
784
785 ENDIF
786
787C loop over the old children of this element
788 DO j=1,idfils(1,tabworkeltri(i))-nb_newfils
789 IF(newknot>=knotlocel(1,dir,idfils(j+1,tabworkeltri(i))).AND.
790 . newknot<=knotlocel(2,dir,idfils(j+1,tabworkeltri(i)))) THEN
791 l_tab_elcut = l_tab_elcut + 1
792 tab_elcut(l_tab_elcut) = idfils(j+1,tabworkeltri(i))
793
794 addelig3d=addelig3d+1
795
796 l_tab_newel = l_tab_newel + 1
797 tab_newel(l_tab_newel) = numelig3d0+addelig3d
798
799 idfils(1,tabworkeltri(i))=idfils(1,tabworkeltri(i))+1
800 idfils(idfils(1,tabworkeltri(i))+1,tabworkeltri(i))=numelig3d0+addelig3d
801 nb_newfils=nb_newfils+1
802
803
804 kxig3d(:,numelig3d0+addelig3d) = kxig3d(:,idfils(j+1,tabworkeltri(i)))
805 iad_ixig3d = sixig3d + addsixig3d + 1
806 kxig3d(4,numelig3d0+addelig3d) = iad_ixig3d
807 DO k=1,kxig3d(3,idfils(j+1,tabworkeltri(i)))
808 ixig3d(iad_ixig3d+k-1) = ixig3d(kxig3d(4,idfils(j+1,tabworkeltri(i)))+k-1)
809 ENDDO
810 addsixig3d = addsixig3d + kxig3d(3,idfils(j+1,tabworkeltri(i)))
811
812 kxig3d(5,numelig3d0+addelig3d) = numelig3d0+addelig3d ! Do not prevent the double id
813 ipartig3d(numelig3d0+addelig3d) = ipartig3d(idfils(j+1,tabworkeltri(i)))
814
815 kxig3d(15,numelig3d0+addelig3d) = inod_ige
816 inod_ige = inod_ige + 64
817
818 kxig3d(idnbcut,numelig3d0+addelig3d) = 1
819
820 knotlocel(:,:,numelig3d0+addelig3d) = knotlocel(:,:,idfils(j+1,tabworkeltri(i)))
821 knotlocel(1,dir,numelig3d0+addelig3d) = newknot
822
823 knotlocel(2,dir,idfils(j+1,tabworkeltri(i))) = newknot
824
825 ENDIF
826 ENDDO
827c
828 kxig3d(idnbcut,tabworkeltri(i)) = kxig3d(idnbcut,tabworkeltri(i)) - 1 ! it is the parent that carries the number of refinements
829
830 ENDDO
831 ENDIF
832
833 DEALLOCATE(tabworkel)
834
835C----------------------------------------------------------------------
836C additional processing for the middle elements, not bordering the segments
837C----------------------------------------------------------------------
838
839 DO i=coinenglob(1,1),coinenglob(1,2)
840 DO j=coinenglob(2,1),coinenglob(2,2)
841 idnextel=meshige(i,j,idmeshige)
842 next = 0
843 DO k=1,l_tabwork ! check that the element is not already listed
844 IF(idnextel==tabworkeltri(k)) THEN
845 next = 1
846 EXIT
847 ENDIF
848 ENDDO
849 IF(next==1) cycle
850
851 iout = 0
852
853 xa(1) = knotlocel(1,dirtang1,idnextel) + tol
854 xa(2) = knotlocel(2,dirtang1,idnextel) - tol
855 xa(3) = knotlocel(2,dirtang1,idnextel) - tol
856 xa(4) = knotlocel(1,dirtang1,idnextel) + tol
857 xa(5) = xa(1)
858
859 ya(1) = knotlocel(1,dirtang2,idnextel) + tol
860 ya(2) = knotlocel(1,dirtang2,idnextel) + tol
861 ya(3) = knotlocel(2,dirtang2,idnextel) - tol
862 ya(4) = knotlocel(2,dirtang2,idnextel) - tol
863 ya(5) = ya(1)
864
865
866 xb=xa(1)-1000
867 yb=ya(1)-2000
868
869 DO k=1,4 ! loop over the 4 corners of the extent of the element
870 intersec=0
871 DO l=1,l_tab_coinknot-1
872
873 xc=coin_tmp(1,l)
874 yc=coin_tmp(2,l)
875 xd=coin_tmp(1,l+1)
876 yd=coin_tmp(2,l+1)
877
878 det = (xb-xa(k))*(yc-yd) - (xc-xd)*(yb-ya(k))
879 IF(det==0) THEN
880 ELSE
881 t1 = ((xc-xa(k))*(yc-yd)-(xc-xd)*(yc-ya(k)))/det
882 t2 = ((xb-xa(k))*(yc-ya(k))-(xc-xa(k))*(yb-ya(k)))/det
883 IF(t1>1.OR.t1<0.OR.t2>1.OR.t2<=0) THEN
884 ELSE
885 intersec = intersec + 1
886 ENDIF
887 ENDIF
888 ENDDO
889 IF(mod(intersec,2)==0) iout=1
890 ENDDO
891
892 IF(iout==1) cycle
893
894 IF(kxig3d(idnbcut,idnextel)==1) THEN
895 idcut_voisin = 0
896 ELSE
897 idcut_voisin = tab_oldidcut(dir,idnextel)-kxig3d(idnbcut,idnextel)+1
898 ENDIF
899 IF(idcut_voisin==idcut) THEN
900
901 IF(flag==0) THEN
902
903 nb_newfils=0
904
905 IF(newknot>=knotlocel(1,dir,idnextel).AND.
906 . newknot<=knotlocel(2,dir,idnextel)) THEN
907 l_tab_elcut = l_tab_elcut + 1
908 l_tab_newel = l_tab_newel + 1
909 ENDIF
910 DO k=1,idfils(1,idnextel)-nb_newfils
911 IF(newknot>=knotlocel(1,dir,idfils(k+1,idnextel)).AND.
912 . newknot<=knotlocel(2,dir,idfils(k+1,idnextel))) THEN
913 l_tab_elcut = l_tab_elcut + 1
914 l_tab_newel = l_tab_newel + 1
915 ENDIF
916 ENDDO
917
918 ELSE
919
920 nb_newfils=0
921
922 IF(newknot>=knotlocel(1,dir,idnextel).AND.
923 . newknot<=knotlocel(2,dir,idnextel)) THEN
924 l_tab_elcut = l_tab_elcut + 1
925 tab_elcut(l_tab_elcut) = idnextel
926
927 addelig3d=addelig3d+1
928
929 l_tab_newel = l_tab_newel + 1
930 tab_newel(l_tab_newel) = numelig3d0+addelig3d
931
932 idfils(1,idnextel)=idfils(1,idnextel)+1
933 idfils(idfils(1,idnextel)+1,idnextel)=numelig3d0+addelig3d
934 nb_newfils=nb_newfils+1
935
936 kxig3d(:,numelig3d0+addelig3d) = kxig3d(:,idnextel)
937 iad_ixig3d = sixig3d + addsixig3d + 1
938 kxig3d(4,numelig3d0+addelig3d) = iad_ixig3d
939 DO m=1,kxig3d(3,idnextel)
940 ixig3d(iad_ixig3d+m-1) = ixig3d(kxig3d(4,idnextel)+m-1)
941 ENDDO
942 addsixig3d = addsixig3d + kxig3d(3,idnextel)
943
944 kxig3d(5,numelig3d0+addelig3d) = numelig3d0+addelig3d ! Do not prevent the double id
945 ipartig3d(numelig3d0+addelig3d) = ipartig3d(idnextel)
946
947 kxig3d(15,numelig3d0+addelig3d) = inod_ige
948 inod_ige = inod_ige + 64
949
950
951 kxig3d(idnbcut,numelig3d0+addelig3d) = 1
952
953 knotlocel(:,:,numelig3d0+addelig3d) = knotlocel(:,:,idnextel)
954 knotlocel(1,dir,numelig3d0+addelig3d) = newknot
955
956 knotlocel(2,dir,idnextel) = newknot
957
958 ENDIF
959
960 DO k=1,idfils(1,idnextel)-nb_newfils
961 IF(newknot>=knotlocel(1,dir,idfils(k+1,idnextel)).AND.
962 . newknot<=knotlocel(2,dir,idfils(k+1,idnextel))) THEN
963 l_tab_elcut = l_tab_elcut + 1
964 tab_elcut(l_tab_elcut) = idfils(k+1,idnextel)
965
966 addelig3d=addelig3d+1
967
968 l_tab_newel = l_tab_newel + 1
969 tab_newel(l_tab_newel) = numelig3d0+addelig3d
970
971 idfils(1,idnextel)=idfils(1,idnextel)+1
972 idfils(idfils(1,idnextel)+1,idnextel)=numelig3d0+addelig3d
973 nb_newfils=nb_newfils+1
974
975 kxig3d(:,numelig3d0+addelig3d) = kxig3d(:,idfils(k+1,idnextel))
976 iad_ixig3d = sixig3d + addsixig3d + 1
977 kxig3d(4,numelig3d0+addelig3d) = iad_ixig3d
978 DO m=1,kxig3d(3,idfils(k+1,idnextel))
979 ixig3d(iad_ixig3d+m-1) = ixig3d(kxig3d(4,idfils(k+1,idnextel))+m-1)
980 ENDDO
981 addsixig3d = addsixig3d + kxig3d(3,idfils(k+1,idnextel))
982
983 kxig3d(5,numelig3d0+addelig3d) = numelig3d0+addelig3d ! n'empeche pas les id doubles
984 ipartig3d(numelig3d0+addelig3d) = ipartig3d(idfils(k+1,idnextel))
985
986 kxig3d(15,numelig3d0+addelig3d) = inod_ige
987 inod_ige = inod_ige + 64
988
989 kxig3d(idnbcut,numelig3d0+addelig3d) = 1
990
991 knotlocel(:,:,numelig3d0+addelig3d) = knotlocel(:,:,idfils(k+1,idnextel))
992 knotlocel(1,dir,numelig3d0+addelig3d) = newknot
993
994 knotlocel(2,dir,idfils(k+1,idnextel)) = newknot
995
996 ENDIF
997 ENDDO
998
999 kxig3d(idnbcut,idnextel) = kxig3d(idnbcut,idnextel) - 1 ! it is the parent that carries the number of refinements
1000
1001 ENDIF
1002 ENDIF
1003 ENDDO
1004 ENDDO
1005
1006 DEALLOCATE(tabworkeltri)
1007
1008 RETURN
1009 END
subroutine comput_coinknot(iel, ixig3d, kxig3d, meshige, ptang1, ptang2, iddir, idtang1, idtang2, neldir, neltang1, neltang2, dir, tab_coinknot, l_tab_coinknot, tab_elcut, l_tab_elcut, tab_newel, l_tab_newel, knot, iad_knot, nknot1, nknot2, nknot3, idfils, knotlocel, newknot, ipartig3d, tab_oldidcut, idcut, flag)
#define my_real
Definition cppsort.cpp:32
end diagonal values have been computed in the(sparse) matrix id.SOL
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
void my_orders(int *mode, int *iwork, int *data, int *index, int *n, int *irecl)
Definition my_orders.c:82