OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
xfeconnec4n.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!|| xfeconnec4n ../engine/source/output/anim/generate/xfeconnec4n.F
25!||--- called by ------------------------------------------------------
26!|| xfecut ../engine/source/output/anim/generate/xfecut.F
27!||--- uses -----------------------------------------------------
28!|| crackxfem_mod ../engine/share/modules/crackxfem_mod.F
29!||====================================================================
30 SUBROUTINE xfeconnec4n(JFT ,JLT ,NFT ,IXC ,ELCUTC ,
31 . IEL_CRK ,IADC_CRK,ILEV ,NODEDGE,CRKEDGE,
32 . XEDGE4N )
33C-----------------------------------------------
34 USE crackxfem_mod
35C-----------------------------------------------
36C I m p l i c i t T y p e s
37C-----------------------------------------------
38#include "implicit_f.inc"
39C-----------------------------------------------
40C G l o b a l P a r a m e t e r s
41C-----------------------------------------------
42#include "mvsiz_p.inc"
43#include "comlock.inc"
44C-----------------------------------------------
45C C o m m o n B l o c k s
46C-----------------------------------------------
47#include "com_xfem1.inc"
48C-----------------------------------------------
49C D u m m y A r g u m e n t s
50C-----------------------------------------------
51 INTEGER JFT,JLT,NFT,ILEV,IXC(NIXC,*),ELCUTC(2,*),IEL_CRK(*),
52 . IADC_CRK(4,*),XEDGE4N(4,*),NODEDGE(2,*)
53 TYPE (XFEM_EDGE_) , DIMENSION(*) :: CRKEDGE
54C-----------------------------------------------
55C L o c a l V a r i a b l e s
56C-----------------------------------------------
57 INTEGER I,K,K1,K2,K3,K4,KK,p1,p2,NOD1,NOD2,IED0,IED1,IED2,
58 . IEDGE,IEDGE1,IEDGE2,EDGE,EDGE1,EDGE2,IXEL,ILAY,ITRI,ILV,
59 . fac,ok,sn,elcut,elcrk,iadc1,iadc2,iadc3,iadc4,nx1,nx2,nx3,nx4
60 INTEGER IFI0(4,MVSIZ),POS(2),
61 . IED(4),IXI(4),D(4),DX(8),NX(MVSIZ)
62 my_real
63 . xint,yint,zint,xm,ym,zm,x10,y10,z10,x20,y20,z20,beta
64 my_real
65 . xin(4,mvsiz),yin(4,mvsiz),zin(4,mvsiz),
66 . xx(4,mvsiz),yy(4,mvsiz),zz(4,mvsiz)
67c---------------------
68 DATA d /2,3,4,1/
69 DATA dx/1,2,3,4,1,2,3,4/
70C=======================================================================
71c Re-build phantom element connectivities
72C-----------------------------------------------
73 ixel = mod(ilev-1, nxel) + 1
74 ilay = (ilev-ixel)/nxel + 1
75 p1 = 0
76 p2 = 0
77c
78 DO i=jft,jlt
79 xin(1,i) = zero
80 yin(1,i) = zero
81 zin(1,i) = zero
82 xin(2,i) = zero
83 yin(2,i) = zero
84 zin(2,i) = zero
85 xin(3,i) = zero
86 yin(3,i) = zero
87 zin(3,i) = zero
88 xin(4,i) = zero
89 yin(4,i) = zero
90 zin(4,i) = zero
91 END DO
92c-----------------
93 DO i=jft,jlt
94 elcrk = iel_crk(i+nft)
95 iadc1 = iadc_crk(1,elcrk)
96 iadc2 = iadc_crk(2,elcrk)
97 iadc3 = iadc_crk(3,elcrk)
98 iadc4 = iadc_crk(4,elcrk)
99C
100 ifi0(1,i) = xfem_phantom(ilay)%IFI(iadc1)
101 ifi0(2,i) = xfem_phantom(ilay)%IFI(iadc2)
102 ifi0(3,i) = xfem_phantom(ilay)%IFI(iadc3)
103 ifi0(4,i) = xfem_phantom(ilay)%IFI(iadc4)
104C
105 ifi0(1,i) = isign(1,ifi0(1,i))
106 ifi0(2,i) = isign(1,ifi0(2,i))
107 ifi0(3,i) = isign(1,ifi0(3,i))
108 ifi0(4,i) = isign(1,ifi0(4,i))
109C--------------
110c Copy local phantom node coordinates (per ILEV)
111C--------------
112c node 1:
113 xx(1,i) = crkavx(ilev)%X(1,iadc1)
114 yy(1,i) = crkavx(ilev)%X(2,iadc1)
115 zz(1,i) = crkavx(ilev)%X(3,iadc1)
116c node 2:
117 xx(2,i) = crkavx(ilev)%X(1,iadc2)
118 yy(2,i) = crkavx(ilev)%X(2,iadc2)
119 zz(2,i) = crkavx(ilev)%X(3,iadc2)
120c node 3:
121 xx(3,i) = crkavx(ilev)%X(1,iadc3)
122 yy(3,i) = crkavx(ilev)%X(2,iadc3)
123 zz(3,i) = crkavx(ilev)%X(3,iadc3)
124c node 4:
125 xx(4,i) = crkavx(ilev)%X(1,iadc4)
126 yy(4,i) = crkavx(ilev)%X(2,iadc4)
127 zz(4,i) = crkavx(ilev)%X(3,iadc4)
128 END DO
129c-----------------------------------------------
130c calculate intersection coordinates of cut edges : XIN, YIN, ZIN
131c-----------------------------------------------
132 DO i=jft,jlt
133 elcrk = iel_crk(i+nft) ! xfem sys number
134 elcut = xfem_phantom(ilay)%ELCUT(elcrk)
135 IF (elcut /= 0) THEN
136 DO k=1,4
137 ied0 = crkedge(ilay)%IEDGEC(k,elcrk) ! = 0,1,2
138 IF (ied0 > 0) THEN
139 edge = xedge4n(k,elcrk) ! global xfem edge number
140 beta = crkedge(ilay)%RATIO(edge)
141 nod1 = nodedge(1,edge)
142 nod2 = nodedge(2,edge)
143 IF (nod1 == ixc(k+1,i+nft) .and.
144 . nod2 == ixc(d(k)+1,i+nft)) THEN
145 p1 = k
146 p2 = d(k)
147 ELSEIF (nod2 == ixc(k+1,i+nft) .and.
148 . nod1 == ixc(d(k)+1,i+nft)) THEN
149 p1 = d(k)
150 p2 = k
151 ENDIF
152 x10 = xx(p1,i)
153 y10 = yy(p1,i)
154 z10 = zz(p1,i)
155 x20 = xx(p2,i)
156 y20 = yy(p2,i)
157 z20 = zz(p2,i)
158 xin(ied0,i) = x10 + beta*(x20-x10)
159 yin(ied0,i) = y10 + beta*(y20-y10)
160 zin(ied0,i) = z10 + beta*(z20-z10)
161 END IF
162 END DO
163 END IF
164 END DO
165c
166c-----------------------------------------------
167c main loop over elements
168C SIMPLE CRACKED ELEMENT
169C only one crack inside element
170c-----------------------------------------------
171 DO i=jft,jlt
172 IF (elcutc(1,i+nft) == 0) cycle ! element standard pas coupe
173c
174 elcrk = iel_crk(i+nft)
175 elcut = xfem_phantom(ilay)%ELCUT(elcrk)
176 IF (elcut == 0) cycle ! element fantome pas coupe
177 itri = xfem_phantom(ilay)%ITRI(1,elcrk)
178 nx1 = xfem_phantom(ilay)%ITRI(2,elcrk)
179 nx2 = dx(nx1+1)
180 nx3 = dx(nx1+2)
181 nx4 = dx(nx1+3)
182c-------------------------------
183 IF (itri == 0) THEN ! element cut by two phantoms
184c-------------------------------
185 pos(1) = 0
186 pos(2) = 0
187 fac = 0
188c
189 IF (ixel == 1) THEN ! first phantom (positif)
190C---
191 DO k=1,4
192 ied0 = crkedge(ilay)%IEDGEC(k,elcrk)
193 IF (ied0 > 0) THEN
194 fac = fac + 1
195 p1 = k
196 p2 = d(p1)
197 IF (ifi0(p1,i) < 0) pos(fac) = p1 ! save negative nodes
198 IF (ifi0(p2,i) < 0) pos(fac) = p2
199 ied(fac) = p1
200 IF (fac == 2) EXIT
201 ENDIF
202 END DO
203C---
204 IF (pos(1) /= 0 .and. pos(2) /= 0) THEN
205 DO k=1,2
206 iedge = crkedge(ilay)%IEDGEC(ied(k),elcrk)
207 IF (iedge > 0) THEN
208c move negative nodes to intersection position
209 kk = crkshell(ilev)%XNODEL(pos(k),elcrk)
210 kk = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
211 crkavx(ilev)%XX(1,kk) = xin(iedge,i)
212 crkavx(ilev)%XX(2,kk) = yin(iedge,i)
213 crkavx(ilev)%XX(3,kk) = zin(iedge,i)
214 ENDIF
215 END DO
216 ENDIF
217C---
218 ELSE IF (ixel == 2) THEN ! second phantom (negatif)
219C---
220 DO k=1,4
221 ied0 = crkedge(ilay)%IEDGEC(k,elcrk)
222 IF (ied0 > 0) THEN
223 fac = fac + 1
224 p1 = k
225 p2 = d(p1)
226 IF (ifi0(p1,i) > 0) pos(fac) = p1 ! save positive nodes
227 IF (ifi0(p2,i) > 0) pos(fac) = p2
228 ied(fac) = p1
229 IF (fac == 2) EXIT
230 ENDIF
231 END DO
232C---
233 IF (pos(1) /= 0 .and. pos(2) /= 0) THEN
234 DO k=1,2
235 iedge = crkedge(ilay)%IEDGEC(ied(k),elcrk)
236 IF (iedge > 0) THEN
237c move positive nodes to intersection position
238 kk = crkshell(ilev)%XNODEL(pos(k),elcrk)
239 kk = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
240 crkavx(ilev)%XX(1,kk) = xin(iedge,i)
241 crkavx(ilev)%XX(2,kk) = yin(iedge,i)
242 crkavx(ilev)%XX(3,kk) = zin(iedge,i)
243 ENDIF
244 END DO
245 ENDIF
246C---
247 ELSE IF (ixel == 3) THEN ! third phantom not actif
248 ENDIF
249
250c---------------------------------
251 ELSEIF (itri < 0) THEN
252c---------------------------------
253 ied1 = nx1
254 ied2 = nx4
255 iadc1 = iadc_crk(nx1,elcrk)
256 iadc2 = iadc_crk(nx2,elcrk)
257 iadc3 = iadc_crk(nx3,elcrk)
258 iadc4 = iadc_crk(nx4,elcrk)
259c
260 iedge1 = crkedge(ilay)%IEDGEC(ied1,elcrk)
261 iedge2 = crkedge(ilay)%IEDGEC(ied2,elcrk)
262 edge1 = xedge4n(ied1,elcrk) ! global xfem edge number
263 edge2 = xedge4n(ied2,elcrk) ! global xfem edge number
264c
265 kk = crkshell(ilev)%XNODEL(nx1,elcrk)
266 k1 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
267 kk = crkshell(ilev)%XNODEL(nx2,elcrk)
268 k2 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
269 kk = crkshell(ilev)%XNODEL(nx3,elcrk)
270 k3 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
271 kk = crkshell(ilev)%XNODEL(nx4,elcrk)
272 k4 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
273c print*,' ELCRK,ELCUT,ITRI=',ELCRK,ELCUT,ITRI
274c print*,' ILEV,NX1=',ILEV,NX1
275c print*,' ied1,EDGE1=',ied1,EDGE1
276c print*,' ied2,EDGE2=',ied2,EDGE2
277
278c--------
279 IF (ixel == 1) THEN
280c NX1 ! unchanged
281c NX2 -> intersec( Nx1, Nx2) , edge1
282c NX3 -> intersec( Nx4, Nx1)
283c NX4 -> NX3
284 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
285 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
286 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
287 crkavx(ilev)%XX(1,k3) = xin(iedge2,i)
288 crkavx(ilev)%XX(2,k3) = yin(iedge2,i)
289 crkavx(ilev)%XX(3,k3) = zin(iedge2,i)
290 crkavx(ilev)%XX(1,k4) = crkavx(ilev)%XX(1,k3)
291 crkavx(ilev)%XX(2,k4) = crkavx(ilev)%XX(2,k3)
292 crkavx(ilev)%XX(3,k4) = crkavx(ilev)%XX(3,k3)
293c
294 ELSE IF (ixel == 2) THEN
295c NX1 -> intersec( Nx1, Nx4)
296c NX2 -> N3
297c NX3 ! unchanged
298c NX4 ! unchanged
299 crkavx(ilev)%XX(1,k1) = xin(iedge2,i)
300 crkavx(ilev)%XX(2,k1) = yin(iedge2,i)
301 crkavx(ilev)%XX(3,k1) = zin(iedge2,i)
302 crkavx(ilev)%XX(1,k2) = crkavx(ilev)%XX(1,k3)
303 crkavx(ilev)%XX(2,k2) = crkavx(ilev)%XX(2,k3)
304 crkavx(ilev)%XX(3,k2) = crkavx(ilev)%XX(3,k3)
305c
306 ELSE IF (ixel == 3) THEN
307c NX1 -> intersec( Nx1, Nx2), edge1
308c NX2 ! unchanged
309c NX3 ! unchanged
310c NX4 -> intersec( Nx4, Nx1), edge2, unchanged if IED2 = TIP
311c
312 crkavx(ilev)%XX(1,k1) = xin(iedge1,i)
313 crkavx(ilev)%XX(2,k1) = yin(iedge1,i)
314 crkavx(ilev)%XX(3,k1) = zin(iedge1,i)
315
316 crkavx(ilev)%XX(1,k4) = xin(iedge2,i)
317 crkavx(ilev)%XX(2,k4) = yin(iedge2,i)
318 crkavx(ilev)%XX(3,k4) = zin(iedge2,i)
319
320c print*,'NX1,K1=',NX1,k1
321c print*,' x,y=',CRKAVX(ILEV)%XX(1,K1),CRKAVX(ILEV)%XX(2,K1)
322c print*,'NX2,K2=',NX2,K2
323c print*,' x,y=',CRKAVX(ILEV)%XX(1,K2),CRKAVX(ILEV)%XX(2,K2)
324c print*,'NX3,K3=',NX3,K3
325c print*,' x,y=',CRKAVX(ILEV)%XX(1,K3),CRKAVX(ILEV)%XX(2,K3)
326c print*,'NX4,K4=',NX4,K4
327c print*,' x,y=',CRKAVX(ILEV)%XX(1,K4),CRKAVX(ILEV)%XX(2,K4)
328c
329 END IF ! IXEL
330
331c---------------------------------
332 ELSEIF (itri > 0) THEN
333c---------------------------------
334 ied1 = nx1
335 ied2 = nx4
336 iadc1 = iadc_crk(nx1,elcrk)
337 iadc2 = iadc_crk(nx2,elcrk)
338 iadc3 = iadc_crk(nx3,elcrk)
339 iadc4 = iadc_crk(nx4,elcrk)
340c
341 iedge1 = crkedge(ilay)%IEDGEC(ied1,elcrk)
342 iedge2 = crkedge(ilay)%IEDGEC(ied2,elcrk)
343 edge1 = xedge4n(ied1,elcrk) ! global xfem edge number
344 edge2 = xedge4n(ied2,elcrk) ! global xfem edge number
345c
346 kk = crkshell(ilev)%XNODEL(nx1,elcrk)
347 k1 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
348 kk = crkshell(ilev)%XNODEL(nx2,elcrk)
349 k2 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
350 kk = crkshell(ilev)%XNODEL(nx3,elcrk)
351 k3 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
352 kk = crkshell(ilev)%XNODEL(nx4,elcrk)
353 k4 = kk - crknod(ilev)%CRKNUMNODS * (ilev-1)
354
355c print*,' '
356c print*,' ELCRK,ELCUT,ITRI=',ELCRK,ELCUT,ITRI
357c print*,' ILEV,NX1=',ILEV,NX1
358c print*,' ied1,EDGE1=',ied1,EDGE1
359c print*,' ied2,EDGE2=',ied2,EDGE2
360c--------
361 IF (ixel == 1) THEN
362c NX1 -> intersec( Nx1, Nx2), edge1 = tip
363c NX2 = const
364c NX3 = const
365c NX4 -> N3
366 crkavx(ilev)%XX(1,k1) = xin(iedge1,i)
367 crkavx(ilev)%XX(2,k1) = yin(iedge1,i)
368 crkavx(ilev)%XX(3,k1) = zin(iedge1,i)
369 crkavx(ilev)%XX(1,k4) = crkavx(ilev)%XX(1,k3)
370 crkavx(ilev)%XX(2,k4) = crkavx(ilev)%XX(2,k3)
371 crkavx(ilev)%XX(3,k4) = crkavx(ilev)%XX(3,k3)
372c
373c--------
374 ELSE IF (ixel == 2) THEN
375c NX1 = const
376c NX2 -> intersec( Nx1, Nx2) , edge1 = tip
377c NX3=NX2
378c NX4 -> intersec( Nx4, Nx1) , edge2
379 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
380 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
381 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
382 crkavx(ilev)%XX(1,k3) = crkavx(ilev)%XX(1,k2)
383 crkavx(ilev)%XX(2,k3) = crkavx(ilev)%XX(2,k2)
384 crkavx(ilev)%XX(3,k3) = crkavx(ilev)%XX(3,k2)
385 crkavx(ilev)%XX(1,k4) = xin(iedge2,i)
386 crkavx(ilev)%XX(2,k4) = yin(iedge2,i)
387 crkavx(ilev)%XX(3,k4) = zin(iedge2,i)
388c
389c--------
390 ELSE IF (ixel == 3) THEN
391c NX1 -> intersec( Nx4, Nx1) , edge2
392c NX2 -> intersec( Nx1, Nx2) , edge1 (unchanged if Nx2 move)
393c NX3 -> unchanged
394c NX4 -> unchanged
395
396 crkavx(ilev)%XX(1,k1) = xin(iedge2,i)
397 crkavx(ilev)%XX(2,k1) = yin(iedge2,i)
398 crkavx(ilev)%XX(3,k1) = zin(iedge2,i)
399c CRKAVX(ILEV)%XX(1,K2) = CRKAVX(ILEV-2)%XX(1,K1)
400c CRKAVX(ILEV)%XX(2,K2) = CRKAVX(ILEV-2)%XX(2,K1)
401c CRKAVX(ILEV)%XX(3,K2) = CRKAVX(ILEV-2)%XX(3,K1)
402
403 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
404 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
405 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
406c
407! print*,'NX1,K1=',NX1,k1
408! print*,' x,y=',CRKAVX(ILEV)%XX(1,K1),CRKAVX(ILEV)%XX(2,K1)
409! print*,'NX2,K2=',NX2,K2
410! print*,' x,y=',CRKAVX(ILEV)%XX(1,K2),CRKAVX(ILEV)%XX(2,K2)
411! print*,'NX3,K3=',NX3,K3
412! print*,' x,y=',CRKAVX(ILEV)%XX(1,K3),CRKAVX(ILEV)%XX(2,K3)
413! print*,'NX4,K4=',NX4,K4
414! print*,' x,y=',CRKAVX(ILEV)%XX(1,K4),CRKAVX(ILEV)%XX(2,K4)
415
416 END IF ! IXEL
417C---
418 END IF ! ITRI
419C-----------------
420 ENDDO ! I=JFT,JLT
421C-----------------
422 RETURN
423 END
my_real
#define my_real
Definition cppsort.cpp:32
crackxfem_mod
Definition crackxfem_mod.F:121
crackxfem_mod::xfem_phantom
type(xfem_phantom_), dimension(:), allocatable xfem_phantom
Definition crackxfem_mod.F:139
crackxfem_mod::crknod
type(xfem_nodes_), dimension(:), allocatable crknod
Definition crackxfem_mod.F:218
crackxfem_mod::crkavx
type(xfem_avx_), dimension(:), allocatable crkavx
Definition crackxfem_mod.F:163
crackxfem_mod::crkshell
type(xfem_shell_), dimension(:), allocatable crkshell
Definition crackxfem_mod.F:190
xfeconnec4n
subroutine xfeconnec4n(jft, jlt, nft, ixc, elcutc, iel_crk, iadc_crk, ilev, nodedge, crkedge, xedge4n)
Definition xfeconnec4n.F:33