OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
crklen3n_adv.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!|| crklen3n_adv ../engine/source/elements/xfem/crklen3n_adv.F
25!||--- called by ------------------------------------------------------
26!|| c3forc3 ../engine/source/elements/sh3n/coque3n/c3forc3.F
27!||--- calls -----------------------------------------------------
28!|| arret ../engine/source/system/arret.f
29!||--- uses -----------------------------------------------------
30!|| crackxfem_mod ../engine/share/modules/crackxfem_mod.F
31!|| element_mod ../common_source/modules/elements/element_mod.F90
32!||====================================================================
33 SUBROUTINE crklen3n_adv(
34 . NEL ,NFT ,ILAY ,NLAY ,IXTG ,
35 . CRKLEN ,ELCRKINI ,IEL_CRKTG,DIR1 ,DIR2 ,
36 . NODEDGE ,CRKEDGE ,XEDGE3N ,NGL ,XL2 ,
37 . XL3 ,YL2 ,YL3 ,ALDT )
38c-----------------------------------------------
39C crack advancement, shells 3N
40c-----------------------------------------------
41C M o d u l e s
42C-----------------------------------------------
44 use element_mod , only : nixtg
45C-----------------------------------------------
46C I m p l i c i t T y p e s
47C-----------------------------------------------
48#include "implicit_f.inc"
49C-----------------------------------------------
50C C o m m o n B l o c k s
51C-----------------------------------------------
52#include "units_c.inc"
53#include "com_xfem1.inc"
54C-----------------------------------------------
55C D u m m y A r g u m e n t s
56C-----------------------------------------------
57 INTEGER NEL,NFT,ILAY,NLAY
58 INTEGER IXTG(NIXTG,*),NGL(NEL),IEL_CRKTG(*),ELCRKINI(NLAY,*),
59 . NODEDGE(2,*),XEDGE3N(3,*)
60 my_real DIR1(NLAY,NEL),DIR2(NLAY,NEL),CRKLEN(NEL),ALDT(NEL)
61 TYPE (XFEM_EDGE_) , DIMENSION(*) :: CRKEDGE
62 my_real, DIMENSION(NEL) :: xl2,yl2,xl3,yl3
63C-----------------------------------------------
64C L o c a l V a r i a b l e s
65C-----------------------------------------------
66 INTEGER I,K,IR,p1,p2,NEWCRK,IED,IED1,IED2,FAC,OK,ICRK,
67 . NOD1,NOD2,ELCRK,ELCRKTG,IEDGE,ICUT
68 INTEGER JCT(NEL),EDGEL(3,NEL),TIP(NEL)
69 INTEGER DD(3),D(6),INV(2)
70C----------
71 my_real, DIMENSION(NEL) :: xl1,yl1
72 my_real, DIMENSION(2,NEL) :: xin,yin
73 my_real, DIMENSION(3,NEL) :: xxl,yyl,len
74 my_real beta0(3,nel)
75 my_real beta,xint,yint,bmin,bmax,x10,y10,x20,y20,
76 . m12,mm,cross1,cross12,xint0,yint0,dir11,dir22
77c
78 DATA d/1,2,2,3,1,3/
79 DATA dd/2,3,1/
80 parameter(bmin = 0.01, bmax = 0.99)
81C=======================================================================
82 newcrk = 0
83 DO i=1,nel
84 jct(i) = 0
85 IF (elcrkini(ilay,i) == 5) THEN ! crack propagation
86 newcrk = newcrk + 1
87 jct(newcrk) = i
88 elcrkini(ilay,i) = 2 ! reset for propagation
89 ELSEIF (elcrkini(ilay,i) == -5) THEN ! crack initialization
90 crklen(i) = aldt(i)
91 elcrkini(ilay,i) = -1 ! reset for initialization
92 ENDIF
93 ENDDO
94 IF (newcrk == 0) RETURN
95C---
96 DO i=1,nel
97 beta0(1:3,i) = zero
98 tip(i) = 0
99 edgel(1,i)=0
100 edgel(2,i)=0
101 edgel(3,i)=0
102 xin(1,i) = zero ! first inters point in local skew
103 yin(1,i) = zero
104 xin(2,i) = zero ! second inters point in local skew
105 yin(2,i) = zero
106 ENDDO
107C
108 inv(1) = 2
109 inv(2) = 1
110C
111C---
112C advance crack inside uncut elements in the layer
113C---
114 DO ir=1,newcrk ! loop over elements with advancing cracks
115 i = jct(ir)
116 elcrktg = iel_crktg(i+nft)
117 ok = 0
118 icut = 0
119 ied = 0
120 DO k=1,3
121 iedge = xedge3n(k,elcrktg)
122 icut = crkedge(ilay)%ICUTEDGE(iedge)
123 nod1 = nodedge(1,iedge)
124 nod2 = nodedge(2,iedge)
125 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i)) THEN
126 p1 = k
127 p2 = dd(k)
128 ELSE IF (nod2 == ixtg(k+1,i) .and. nod1 == ixtg(dd(k)+1,i)) THEN
129 p1 = dd(k)
130 p2 = k
131 ENDIF
132 IF (icut == 1) THEN
133 ok = ok + 1
134 ied = k
135c tag
136 icrk = crkedge(ilay)%EDGEICRK(iedge)
137 EXIT
138 ENDIF ! IF(ICUT == 1)THEN
139 ENDDO ! DO K=1,4
140C---
141 IF (ok /= 1) THEN
142 WRITE(iout,*) 'ERROR IN ADVANCING CRACK --- CHECK CRACK TIP'
143 CALL arret(2)
144 ENDIF
145C---
146 edgel(ied,i) = 1
147 iedge = xedge3n(ied,elcrktg)
148 tip(i) = crkedge(ilay)%EDGETIP(1,iedge)
149C---
150 END DO ! DO IR=1,NEWCRK
151C----------------------------------------------------------------------
152c local coords
153C----------------------------------------------------------------------
154 DO i=1,nel
155 xl1(i) = zero
156 yl1(i) = zero
157 xxl(1,i)=xl1(i)
158 yyl(1,i)=yl1(i)
159 xxl(2,i)=xl2(i)
160 yyl(2,i)=yl2(i)
161 xxl(3,i)=xl3(i)
162 yyl(3,i)=yl3(i)
163 ENDDO
164C---
165 DO i=1,nel
166 len(1,i) = (xl2(i)-xl1(i))*(xl2(i)-xl1(i))
167 . + (yl2(i)-yl1(i))*(yl2(i)-yl1(i))
168 len(2,i) = (xl3(i)-xl2(i))*(xl3(i)-xl2(i))
169 . + (yl3(i)-yl2(i))*(yl3(i)-yl2(i))
170 len(3,i) = (xl1(i)-xl3(i))*(xl1(i)-xl3(i))
171 . + (yl1(i)-yl3(i))*(yl1(i)-yl3(i))
172 ENDDO
173C------------------------------------------------
174C - intersections -
175C------------------------------------------------
176 DO ir=1,newcrk
177 i=jct(ir)
178 elcrktg = iel_crktg(i+nft)
179 elcrk = elcrktg + ecrkxfec
180 ied1 = 0
181 ied2 = 0
182 DO k=1,3
183 IF(edgel(k,i) > 0)THEN
184 ied1 = edgel(k,i)
185 ied2 = inv(ied1)
186 EXIT
187 END IF
188 END DO
189 DO k=1,3
190 iedge = xedge3n(k,elcrktg)
191 IF (iedge > 0 .and. edgel(k,i) == 1) THEN
192 icut = crkedge(ilay)%ICUTEDGE(iedge)
193 IF (icut == 1) THEN
194 beta = crkedge(ilay)%RATIO(iedge)
195C
196 IF (beta > one .or. beta == zero) THEN
197 WRITE(*,*) 'ERROR NEGATIV BETA, NO INTERSECTION!'
198 CALL arret(2)
199 ENDIF
200C
201 nod1 = nodedge(1,iedge)
202 nod2 = nodedge(2,iedge)
203 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i)) THEN
204 p1 = k
205 p2 = dd(k)
206 ELSEIF (nod2 == ixtg(k+1,i).and.nod1==ixtg(dd(k)+1,i)) THEN
207 p1 = dd(k)
208 p2 = k
209 ENDIF
210 x10 = xxl(p1,i)
211 y10 = yyl(p1,i)
212 x20 = xxl(p2,i)
213 y20 = yyl(p2,i)
214C
215 xint = x10+beta*(x20-x10)
216 yint = y10+beta*(y20-y10)
217 xin(ied1,i) = xint
218 yin(ied1,i) = yint
219 ENDIF
220 ENDIF
221 ENDDO
222C---
223 IF (ied1 == 0 .or. ied2 == 0) GOTO 130
224 xint0 = xin(ied1,i)
225 yint0 = yin(ied1,i)
226C---
227 dir11 = -dir2(ilay,i)
228 dir22 = dir1(ilay,i)
229C---
230 IF (dir11 == zero) THEN
231 DO 140 k=1,3
232 xint = zero
233 yint = zero
234 elcrktg = iel_crktg(i+nft)
235 elcrk = elcrktg + ecrkxfec
236 iedge = xedge3n(k,elcrktg)
237 nod1 = nodedge(1,iedge)
238 nod2 = nodedge(2,iedge)
239 IF(nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
240 p1 = k
241 p2 = dd(k)
242 ELSE IF(nod2 == ixtg(k+1,i).and.nod1==ixtg(dd(k)+1,i))THEN
243 p1 = dd(k)
244 p2 = k
245 ENDIF
246C
247 IF (edgel(k,i) == ied1) GOTO 140
248 IF (xxl(p1,i) == xxl(p2,i)) GOTO 140 ! no inters (parallel to dir 2)
249 m12 = xxl(p2,i)-xxl(p1,i)
250 m12 = (yyl(p2,i)-yyl(p1,i))/m12
251 xint = xint0
252 yint = yyl(p1,i)+m12*(xint-xxl(p1,i))
253 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
254 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
255 IF (cross12 > zero) GOTO 140
256c
257 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
258 beta = sqrt(cross1 / len(k,i))
259 beta = max(beta, bmin)
260 beta = min(beta, bmax)
261 beta0(k,i) = beta
262C
263 xin(ied2,i) = xint
264 yin(ied2,i) = yint
265 edgel(k,i) = ied2
266 EXIT
267 140 CONTINUE
268 ELSEIF(dir22 == zero)THEN
269 DO 150 k=1,3
270 xint = zero
271 yint = zero
272 elcrktg = iel_crktg(i+nft)
273 elcrk = elcrktg + ecrkxfec
274 iedge = xedge3n(k,elcrktg)
275 nod1 = nodedge(1,iedge)
276 nod2 = nodedge(2,iedge)
277 IF(nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
278 p1 = k
279 p2 = dd(k)
280 ELSE IF(nod2 == ixtg(k+1,i).and.nod1==ixtg(dd(k)+1,i))THEN
281 p1 = dd(k)
282 p2 = k
283 ENDIF
284C
285 IF (edgel(k,i) == ied1) GOTO 150
286 IF (yyl(p1,i) == yyl(p2,i)) GOTO 150 ! no inters (parallel to dir 2)
287 m12 = yyl(p2,i)-yyl(p1,i)
288 m12 = (xxl(p2,i)-xxl(p1,i))/m12
289 yint = yint0
290 xint = xxl(p1,i)+m12*(yint-yyl(p1,i))
291 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
292 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
293 IF (cross12 > zero) GOTO 150
294C
295 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
296 beta = sqrt(cross1 / len(k,i))
297 beta = max(beta, bmin)
298 beta = min(beta, bmax)
299 beta0(k,i) = beta
300C
301 xin(ied2,i) = xint
302 yin(ied2,i) = yint
303 edgel(k,i) = ied2
304 EXIT
305 150 CONTINUE
306 ELSEIF(dir11 /= zero .AND. dir22 /= zero)THEN
307 DO 160 k=1,3
308 xint = zero
309 yint = zero
310 elcrktg = iel_crktg(i+nft)
311 elcrk = elcrktg + ecrkxfec
312 iedge = xedge3n(k,elcrktg)
313 nod1 = nodedge(1,iedge)
314 nod2 = nodedge(2,iedge)
315 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
316 p1 = k
317 p2 = dd(k)
318 ELSE IF (nod2 == ixtg(k+1,i).and.nod1==ixtg(dd(k)+1,i))THEN
319 p1 = dd(k)
320 p2 = k
321 ENDIF
322C
323 IF (edgel(k,i) == ied1) GOTO 160
324 IF (xxl(p1,i) == xxl(p2,i)) THEN
325 mm = dir22/dir11
326 xint = xxl(p1,i) ! or = XXL(p2,I)
327 yint = yint0+mm*(xint-xint0)
328 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
329 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
330 IF (cross12 > zero) GOTO 160
331C
332 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
333 beta = sqrt(cross1 / len(k,i))
334 beta = max(beta, bmin)
335 beta = min(beta, bmax)
336 beta0(k,i) = beta
337C
338 xin(ied2,i) = xint
339 yin(ied2,i) = yint
340 edgel(k,i) = ied2
341 EXIT
342 ELSE
343 mm = dir22/dir11
344 m12 = xxl(p2,i)-xxl(p1,i)
345 m12 = (yyl(p2,i)-yyl(p1,i))/m12
346 IF (mm == m12) GOTO 160 ! no inters (parallel to dir 2)
347 xint = (yint0-yyl(p1,i)+m12*xxl(p1,i)-mm*xint0)/(m12-mm)
348 yint = yint0+mm*(xint-xint0)
349 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
350 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
351 IF (cross12 > zero) GOTO 160
352C
353 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
354 beta = sqrt(cross1 / len(k,i))
355 beta = max(beta, bmin)
356 beta = min(beta, bmax)
357 beta0(k,i) = beta
358C
359 xin(ied2,i) = xint
360 yin(ied2,i) = yint
361 edgel(k,i) = ied2
362 EXIT
363 ENDIF
364 160 CONTINUE
365 ENDIF
366 130 CONTINUE
367 ENDDO
368C----------------------------------------------------------------------
369C check for getting both intersections
370C
371 DO ir=1,newcrk
372 i = jct(ir)
373 fac = 0
374 DO k=1,3
375 IF (edgel(k,i)==1 .or. edgel(k,i)==2) fac=fac+1
376 ENDDO
377 IF (fac /= 2) THEN
378 WRITE(iout,*) 'ERROR IN ADVANCING CRACK.NO CUT EDGES'
379 CALL arret(2)
380 ENDIF
381 crklen(i) = sqrt((xin(2,i) - xin(1,i))**2 + (yin(2,i) - yin(1,i))**2)
382 ENDDO
383c-----------
384 RETURN
385 END
subroutine crklen3n_adv(nel, nft, ilay, nlay, ixtg, crklen, elcrkini, iel_crktg, dir1, dir2, nodedge, crkedge, xedge3n, ngl, xl2, xl3, yl2, yl3, aldt)
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
real function second()
SECOND Using ETIME
subroutine arret(nn)
Definition arret.F:86