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