OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
crklayer3n_ini.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!|| crklayer3n_ini ../engine/source/elements/xfem/crklayer3n_ini.F
25!||--- called by ------------------------------------------------------
26!|| c3forc3 ../engine/source/elements/sh3n/coque3n/c3forc3.F
27!||--- calls -----------------------------------------------------
28!|| arret ../engine/source/system/arret.F
29!|| lsint4 ../engine/source/elements/xfem/crklayer4n_adv.F
30!||--- uses -----------------------------------------------------
31!|| crackxfem_mod ../engine/share/modules/crackxfem_mod.F
32!|| element_mod ../common_source/modules/elements/element_mod.F90
33!||====================================================================
34 SUBROUTINE crklayer3n_ini(
35 . NEL ,NFT ,ILAY ,NLAY ,IXTG ,
36 . ELCUTC ,ELCRKINI ,IEL_CRKTG ,INOD_CRK ,IAD_CRKTG ,
37 . NODENR ,DIR1 ,DIR2 ,NODEDGE ,CRKNODIAD ,
38 . KNOD2ELC ,CRKEDGE ,XEDGE3N ,NGL ,AREA ,
39 . XL2 ,XL3 ,YL2 ,YL3 )
40C-----------------------------------------------
41C crack initialisation, shells 3N
42C-----------------------------------------------
43C M o d u l e s
44C-----------------------------------------------
45 USE crackxfem_mod
46 use element_mod , only : nixtg
47C-----------------------------------------------
48C I m p l i c i t T y p e s
49C-----------------------------------------------
50#include "implicit_f.inc"
51C-----------------------------------------------
52C C o m m o n B l o c k s
53C-----------------------------------------------
54#include "com04_c.inc"
55#include "com_xfem1.inc"
56#include "comlock.inc"
57#include "units_c.inc"
58C-----------------------------------------------
59C D u m m y A r g u m e n t s
60C-----------------------------------------------
61 INTEGER NEL,NFT,ILAY,NLAY
62 INTEGER IXTG(NIXTG,*),IEL_CRKTG(*),NGL(NEL),
63 . INOD_CRK(*),NODENR(*),IAD_CRKTG(3,*),ELCRKINI(NLAY,*),
64 . ELCUTC(2,*),NODEDGE(2,*),CRKNODIAD(*),KNOD2ELC(*),XEDGE3N(3,*)
65 my_real, DIMENSION(NEL) :: AREA,XL2,YL2,XL3,YL3
66 my_real, DIMENSION(NLAY,NEL) :: dir1,dir2
67 TYPE (XFEM_EDGE_) , DIMENSION(*) :: CRKEDGE
68C-----------------------------------------------
69C L o c a l V a r i a b l e s
70C-----------------------------------------------
71 INTEGER I,K,ILEV,IR,P1,P2,IED,IED1,IED2,PP1,PP2,PP3,FAC,
72 . icrk,newcrk,elcrk,elcrktg,iedge,icut,sigbeta,itri,nm,np,
73 . nx1,nx2,nx3,nod1,nod2,ie1,ie2,ifi1,ifi2,ie10,ie20,iad1,iad2,iad3
74 INTEGER JCT(NEL),EDGEL(3,NEL),IADC(3),ISIGN(3),
75 . DD(3),D(6),DX(6),N(3),NN(3),IENR0(3),IENR(3)
76C----------
77 my_real, DIMENSION(NEL) :: xl1,yl1
78 my_real, DIMENSION(2,NEL) :: xin,yin
79 my_real, DIMENSION(3,NEL) :: xxl,yyl,len,fit
80 my_real beta0(3,nel),xn(3),yn(3),xmi(2),ymi(2)
81 my_real beta,xint,yint,fi,bmin,bmax,m12,mm,cross1,cross12,
82 . dir11,dir22,x1,y1,x2,y2,x3,y3,area1,area2,area3
83c------------------------
84 DATA d/1,2,2,3,1,3/
85 DATA dd/2,3,1/
86 DATA dx/1,2,3,1,2,3/
87 parameter(bmin = 0.01, bmax = 0.99)
88C=======================================================================
89 newcrk = 0
90 DO i=1,nel
91 jct(i) = 0
92 IF (elcrkini(ilay,i) == -1) THEN
93 newcrk = newcrk + 1
94 jct(newcrk) = i
95 ENDIF
96 ENDDO
97 IF (newcrk == 0) RETURN
98C-----------------------
99 pp1 = nxel*(ilay-1) + 1
100 pp2 = pp1 + 1
101 pp3 = pp1 + 2
102c
103 DO i=1,nel
104 beta0(1:3,i) = zero
105C
106 edgel(1,i)=0
107 edgel(2,i)=0
108 edgel(3,i)=0
109C
110 xin(1,i) = zero ! first inters point in local skew
111 yin(1,i) = zero
112 xin(2,i) = zero ! second inters point in local skew
113 yin(2,i) = zero
114 ENDDO
115C----------------------------------------------------------------------
116c local coords
117C----------------------------------------------------------------------
118 DO i=1,nel
119 xl1(i) = zero
120 yl1(i) = zero
121 xxl(1,i) = xl1(i)
122 yyl(1,i) = yl1(i)
123 xxl(2,i) = xl2(i)
124 yyl(2,i) = yl2(i)
125 xxl(3,i) = xl3(i)
126 yyl(3,i) = yl3(i)
127 ENDDO
128c
129 DO i=1,nel
130 len(1,i) = (xl2(i)-xl1(i))*(xl2(i)-xl1(i))
131 . + (yl2(i)-yl1(i))*(yl2(i)-yl1(i))
132 len(2,i) = (xl3(i)-xl2(i))*(xl3(i)-xl2(i))
133 . + (yl3(i)-yl2(i))*(yl3(i)-yl2(i))
134 len(3,i) = (xl1(i)-xl3(i))*(xl1(i)-xl3(i))
135 . + (yl1(i)-yl3(i))*(yl1(i)-yl3(i))
136 ENDDO
137C------------------------------------------------
138C - intersections -
139C------------------------------------------------
140 DO ir=1,newcrk
141 i=jct(ir)
142 fac = 0
143C---
144 dir11 = -dir2(ilay,i)
145 dir22 = dir1(ilay,i)
146C---
147 IF (dir11 == zero) THEN
148 DO 140 k=1,3
149 xint = zero
150 yint = zero
151 elcrktg = iel_crktg(i+nft)
152 elcrk = elcrktg + ecrkxfec
153 iedge = xedge3n(k,elcrktg)
154 nod1 = nodedge(1,iedge)
155 nod2 = nodedge(2,iedge)
156 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
157 p1 = k
158 p2 = dd(k)
159 ELSEIF (nod2 == ixtg(k+1,i).and.nod1 == ixtg(dd(k)+1,i))THEN
160 p1 = dd(k)
161 p2 = k
162 ENDIF
163C
164 IF (xxl(p1,i) == xxl(p2,i)) GOTO 140 ! no inters (parallel to dir 2)
165 m12 = xxl(p2,i)-xxl(p1,i)
166 m12 = (yyl(p2,i)-yyl(p1,i))/m12
167 xint = zero
168 yint = yyl(p1,i) - m12*xxl(p1,i)
169 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
170 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
171C
172 IF(cross12 > zero)GOTO 140
173 fac = fac + 1
174 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
175 beta = sqrt(cross1 / len(k,i))
176 beta = max(beta, bmin)
177 beta = min(beta, bmax)
178 beta0(k,i) = beta
179C---
180 xin(fac,i) = xint
181 yin(fac,i) = yint
182 edgel(k,i) = fac
183 IF(fac == 2)EXIT
184C---
185 140 CONTINUE
186C---
187 ELSEIF(dir22 == zero)THEN
188 DO 150 k=1,3
189 xint = zero
190 yint = zero
191 elcrktg = iel_crktg(i+nft)
192 elcrk = elcrktg + ecrkxfec
193 iedge = xedge3n(k,elcrktg)
194 nod1 = nodedge(1,iedge)
195 nod2 = nodedge(2,iedge)
196 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
197 p1 = k
198 p2 = dd(k)
199 ELSEIF (nod2 == ixtg(k+1,i).and.nod1 == ixtg(dd(k)+1,i))THEN
200 p1 = dd(k)
201 p2 = k
202 ENDIF
203C
204 IF (yyl(p1,i) == yyl(p2,i)) GOTO 150 ! no inters (parallel to dir 2)
205 m12 = yyl(p2,i)-yyl(p1,i)
206 m12 = (xxl(p2,i)-xxl(p1,i))/m12
207 yint = zero
208 xint = xxl(p1,i) - m12*yyl(p1,i)
209 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
210 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
211C
212 IF (cross12 > zero) GOTO 150
213 fac = fac + 1
214C
215 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
216 beta = sqrt(cross1 / len(k,i))
217 beta = max(beta, bmin)
218 beta = min(beta, bmax)
219 beta0(k,i) = beta
220C---
221 xin(fac,i) = xint
222 yin(fac,i) = yint
223 edgel(k,i) = fac
224 IF(fac == 2)EXIT
225C---
226 150 CONTINUE
227C---
228 ELSEIF (dir11 /= zero .AND. dir22 /= zero)THEN
229 DO 160 k=1,3
230 xint = zero
231 yint = zero
232 elcrktg = iel_crktg(i+nft)
233 elcrk = elcrktg + ecrkxfec
234 iedge = xedge3n(k,elcrktg)
235 nod1 = nodedge(1,iedge)
236 nod2 = nodedge(2,iedge)
237 IF (nod1 == ixtg(k+1,i) .and. nod2 == ixtg(dd(k)+1,i))THEN
238 p1 = k
239 p2 = dd(k)
240 ELSEIF (nod2 == ixtg(k+1,i).and.nod1 == ixtg(dd(k)+1,i))THEN
241 p1 = dd(k)
242 p2 = k
243 ENDIF
244C
245 IF (xxl(p1,i) == xxl(p2,i)) THEN
246 mm = dir22/dir11 ! slope of dir (2)
247 xint = xxl(p1,i) ! or = XXL(p2,I)
248 yint = mm*xint
249 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
250 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
251C
252 IF (cross12 > zero) GOTO 160
253 fac = fac + 1
254C
255 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
256 beta = sqrt(cross1 / len(k,i))
257 beta = max(beta, bmin)
258 beta = min(beta, bmax)
259 beta0(k,i) = beta
260C---
261 xin(fac,i) = xint
262 yin(fac,i) = yint
263 edgel(k,i) = fac
264 IF(fac == 2)EXIT
265 ELSE
266 mm = dir22/dir11 ! slope of dir (2)
267 m12 = xxl(p2,i)-xxl(p1,i)
268 m12 = (yyl(p2,i)-yyl(p1,i))/m12
269 IF (mm == m12) GOTO 160 ! no inters (parallel to dir 2)
270 xint = (yyl(p1,i) - m12*xxl(p1,i))/(mm - m12)
271 yint = mm*xint
272 cross12 = (xint-xxl(p1,i))*(xint-xxl(p2,i))+
273 . (yint-yyl(p1,i))*(yint-yyl(p2,i))
274C
275 IF (cross12 > zero) GOTO 160
276 fac = fac + 1
277C
278 cross1 = (xxl(p1,i) - xint)**2 + (yyl(p1,i) - yint)**2
279 beta = sqrt(cross1 / len(k,i))
280 beta = max(beta, bmin)
281 beta = min(beta, bmax)
282 beta0(k,i) = beta
283C---
284 xin(fac,i) = xint
285 yin(fac,i) = yint
286 edgel(k,i) = fac
287 IF(fac == 2)EXIT
288 ENDIF
289C---
290 160 CONTINUE
291C---
292 ENDIF
293 ENDDO
294C
295C check for getting both intersections
296C
297 DO ir=1,newcrk
298 i=jct(ir)
299 fac = 0
300 DO k=1,3
301 IF(edgel(k,i)==1.or.edgel(k,i)==2) fac=fac+1
302 ENDDO
303 IF(fac /= 2)THEN
304 WRITE(iout,*) 'ERROR IN INITIATION CRACK.NO CUT EDGES'
305 CALL arret(2)
306 ENDIF
307 ENDDO
308C---
309 DO ir=1,newcrk
310 i=jct(ir)
311 elcrktg = iel_crktg(i+nft)
312 elcrk = elcrktg + ecrkxfec
313C---
314cc ICRK = ICRK + 1 ! increase nb of cracks
315C---
316 DO k=1,3
317 ied = edgel(k,i)
318 IF (ied == 1 .or. ied == 2) THEN
319 crkedge(ilay)%IEDGETG(k,elcrktg) = ied
320 ENDIF
321 ENDDO
322 ENDDO
323C----------------------------------------------------------------------
324C SIGN DISTANCE FOR NEW CRACKED LAYER
325C----------------------------------------------------------------------
326 DO ir=1,newcrk
327 i=jct(ir)
328 fit(1,i)= zero
329 fit(2,i)= zero
330 fit(3,i)= zero
331 xn(1) = xl1(i)
332 yn(1) = yl1(i)
333 xn(2) = xl2(i)
334 yn(2) = yl2(i)
335 xn(3) = xl3(i)
336 yn(3) = yl3(i)
337 DO k=1,3
338 p1 = k
339 p2 = dd(k)
340 ied = edgel(k,i)
341 IF (ied > 0) THEN
342 xmi(ied) = half*(xn(p1)+xn(p2))
343 ymi(ied) = half*(yn(p1)+yn(p2))
344 ENDIF
345 ENDDO
346C
347 DO k=1,3
348 fi=zero
349 CALL lsint4(xmi(1),ymi(1),xmi(2),ymi(2),xn(k),yn(k),fi)
350 IF (fit(k,i)==zero) fit(k,i) = fi
351 ENDDO
352 ENDDO
353C-----------------------
354 DO ir=1,newcrk
355 i=jct(ir)
356 elcrktg = iel_crktg(i+nft)
357 elcrk = elcrktg + ecrkxfec
358C
359 DO k=1,3
360 ied = edgel(k,i)
361 IF (ied == 1 .or. ied == 2) THEN
362 iedge = xedge3n(k,elcrktg)
363
364 icut = crkedge(ilay)%ICUTEDGE(iedge)
365 IF (icut > 0) THEN ! already cut before => edge connecting two cracks
366 crkedge(ilay)%ICUTEDGE(iedge) = 3 ! only for spmd si 2 cracks sur le meme edge
367 ELSE
368 crkedge(ilay)%ICUTEDGE(iedge) = 2 ! edge cut
369 crkedge(ilay)%RATIO(iedge) = beta0(k,i)
370 ENDIF
371 ENDIF
372 ENDDO
373 ENDDO
374C-----------------------
375c Fill New Cut Phantom Element
376C-----------------------
377 DO ir=1,newcrk
378 i = jct(ir)
379 elcrktg = iel_crktg(i+nft)
380 elcrk = elcrktg + ecrkxfec
381C
382 iadc(1) = iad_crktg(1,elcrktg)
383 iadc(2) = iad_crktg(2,elcrktg)
384 iadc(3) = iad_crktg(3,elcrktg)
385C
386 n(1) = ixtg(2,i)
387 n(2) = ixtg(3,i)
388 n(3) = ixtg(4,i)
389C
390 nn(1) = inod_crk(n(1))
391 nn(2) = inod_crk(n(2))
392 nn(3) = inod_crk(n(3))
393C
394 elcutc(1,i) = 2
395 numelcrk = numelcrk + 1
396C
397 isign(1) = int(sign(one,fit(1,i)))
398 isign(2) = int(sign(one,fit(2,i)))
399 isign(3) = int(sign(one,fit(3,i)))
400C
401 IF (fit(1,i) == zero) isign(1) = 0
402 IF (fit(2,i) == zero) isign(2) = 0
403 IF (fit(3,i) == zero) isign(3) = 0
404c--------------------------------------------
405c
406 icrk = crkshell(pp1)%CRKSHELLID(elcrk)
407 itri = 0
408 nm = 0
409 np = 0
410 DO k=1,3
411 IF (isign(k) > 0) THEN
412 itri = itri + 1
413 np = k
414 ELSEIF (isign(k) < 0) THEN
415 nm = k
416 ENDIF
417 ENDDO
418 IF (itri == 1) THEN
419 itri = -1
420 nx1 = np
421 ELSEIF (itri == 2) THEN
422 itri = 1
423 nx1 = nm
424 ENDIF
425 nx2 = dx(nx1+1)
426 nx3 = dx(nx2+1)
427 xfem_phantom(ilay)%ITRI(1,elcrk) = itri
428 xfem_phantom(ilay)%ITRI(2,elcrk) = nx1
429c--------------------------------------------
430c activate enriched nodes
431c--------------------------------------------
432c
433 ienr0(1) = crknodiad(iadc(1))
434 ienr0(2) = crknodiad(iadc(2))
435 ienr0(3) = crknodiad(iadc(3))
436C
437 ienr(1) = ienr0(1) + knod2elc(nn(1))*(ilay-1)
438 ienr(2) = ienr0(2) + knod2elc(nn(2))*(ilay-1)
439 ienr(3) = ienr0(3) + knod2elc(nn(3))*(ilay-1)
440C
441 sigbeta = 0
442 DO k=1,3
443 ied = edgel(k,i)
444 IF (ied > 0) THEN
445 iedge = xedge3n(k,elcrktg)
446 nod1 = nodedge(1,iedge)
447 nod2 = nodedge(2,iedge)
448 ie10 = 0
449 ie20 = 0
450 ie10 = crkedge(ilay)%EDGEENR(1,iedge)
451 ie20 = crkedge(ilay)%EDGEENR(2,iedge)
452 IF (nod1 == n(k) .and. nod2 == n(dd(k))) THEN
453 ie1 = ienr(k)
454 ie2 = ienr(dd(k))
455 ifi1 = isign(k)
456 ifi2 = isign(dd(k))
457 sigbeta = iedge
458 ELSE IF (nod2 == n(k) .and. nod1 == n(dd(k))) THEN
459 ie1 = ienr(dd(k))
460 ie2 = ienr(k)
461 ifi1 = isign(dd(k))
462 ifi2 = isign(k)
463 sigbeta = -iedge
464 END IF
465 crkedge(ilay)%EDGEENR(1,iedge) = max(ie1,ie10)
466 crkedge(ilay)%EDGEENR(2,iedge) = max(ie2,ie20)
467 IF (crkedge(ilay)%EDGEICRK(iedge) == 0)
468 . crkedge(ilay)%EDGEICRK(iedge) = icrk
469 ENDIF
470 ENDDO
471C------------------
472 crkedge(ilay)%LAYCUT(elcrk) = -1 ! layer cut
473 xfem_phantom(ilay)%ELCUT(elcrk) = icrk
474C------------------
475 DO k=1,3
476 ied = edgel(k,i)
477 iedge = xedge3n(k,elcrktg)
478 IF (ied > 0) THEN
479 crkedge(ilay)%EDGETIP(1,iedge) = ied
480 crkedge(ilay)%EDGETIP(2,iedge) =
481 . crkedge(ilay)%EDGETIP(2,iedge) + 1
482 ENDIF
483 ENDDO
484C
485 xfem_phantom(ilay)%IFI(iadc(1)) = isign(1)
486 xfem_phantom(ilay)%IFI(iadc(2)) = isign(2)
487 xfem_phantom(ilay)%IFI(iadc(3)) = isign(3)
488C------------------
489c IXEL = 1 => positif element within ILAY (ILEV = PP1)
490C------------------
491 ilev = pp1
492C
493 crklvset(ilev)%ENR0(1,iadc(1)) = -ienr(1)
494 crklvset(ilev)%ENR0(1,iadc(2)) = -ienr(2)
495 crklvset(ilev)%ENR0(1,iadc(3)) = -ienr(3)
496C
497 IF (isign(1) > 0) crklvset(ilev)%ENR0(1,iadc(1)) = 0
498 IF (isign(2) > 0) crklvset(ilev)%ENR0(1,iadc(2)) = 0
499 IF (isign(3) > 0) crklvset(ilev)%ENR0(1,iadc(3)) = 0
500C------------------
501c IXEL = 2 => negatif element within ILAY (ILEV = PP2)
502C------------------
503 ilev = pp2
504C------------------
505 crklvset(ilev)%ENR0(1,iadc(1)) = -ienr(1)
506 crklvset(ilev)%ENR0(1,iadc(2)) = -ienr(2)
507 crklvset(ilev)%ENR0(1,iadc(3)) = -ienr(3)
508C
509 IF (isign(1) < 0) crklvset(ilev)%ENR0(1,iadc(1)) = 0
510 IF (isign(2) < 0) crklvset(ilev)%ENR0(1,iadc(2)) = 0
511 IF (isign(3) < 0) crklvset(ilev)%ENR0(1,iadc(3)) = 0
512C------------------
513c IXEL = 3 => Third phantom element (ILEV = PP3)
514C------------------
515 IF (itri < 0) THEN ! sign ILEV3 = ILEV2 < 0
516 ie2 = xedge3n(nx3,elcrktg) ! tip edge
517 IF (crkedge(ilay)%ICUTEDGE(ie2) > 1) THEN
518 sigbeta = -sigbeta
519 iad1 = iadc(nx1)
520 iad2 = iadc(nx2)
521 iad3 = iadc(nx3)
522 crklvset(pp3)%ENR0(1,iad1) = abs(crklvset(pp2)%ENR0(1,iad1))
523 crklvset(pp3)%ENR0(1,iad2) = crklvset(pp2)%ENR0(1,iad2)
524 crklvset(pp3)%ENR0(1,iad3) = crklvset(pp2)%ENR0(1,iad3)
525 crklvset(pp2)%ENR0(1,iad1) = -crknodiad(iad1) - knod2elc(nn(nx1))*(ilay-1)
526c
527c-- AREA factors for each phantom
528 x1 = xxl(nx1,i)
529 y1 = yyl(nx1,i)
530 ied = crkedge(ilay)%IEDGETG(nx1,elcrktg)
531 x2 = xin(ied,i)
532 y2 = yin(ied,i)
533 ied = crkedge(ilay)%IEDGETG(nx3,elcrktg)
534 x3 = xin(ied,i)
535 y3 = yin(ied,i)
536 area1 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
537 area1 = area1 / area(i)
538 x1 = xxl(nx2,i)
539 y1 = yyl(nx2,i)
540 x2 = xxl(nx3,i)
541 y2 = yyl(nx3,i)
542 area2 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
543 area2 = area2 / area(i)
544 area3 = one - area1 - area2
545c
546 ELSE
547 ! print*,' IMPOSSIBLE CASE'
548 ENDIF
549c
550 ELSEIF (itri > 0) THEN ! sign ILEV3 = ILEV1 > 0
551c
552 ie1 = xedge3n(nx1,elcrktg) ! tip edge
553 IF (crkedge(ilay)%ICUTEDGE(ie1) > 1) THEN
554 iad1 = iadc(nx1)
555 iad2 = iadc(nx2)
556 iad3 = iadc(nx3)
557 crklvset(pp3)%ENR0(1,iad1) = abs(crklvset(pp1)%ENR0(1,iad1))
558 crklvset(pp3)%ENR0(1,iad2) = crklvset(pp1)%ENR0(1,iad2)
559 crklvset(pp3)%ENR0(1,iad3) = crklvset(pp1)%ENR0(1,iad3)
560 crklvset(pp1)%ENR0(1,iad1) = -crknodiad(iad1) - knod2elc(nn(nx1))*(ilay-1)
561c
562c-- AREA factors for each phantom
563 ied1= crkedge(ilay)%IEDGETG(nx1,elcrktg)
564 ied2= crkedge(ilay)%IEDGETG(nx3,elcrktg)
565 x1 = xin(ied1,i)
566 y1 = yin(ied1,i)
567 x2 = xxl(nx2,i)
568 y2 = yyl(nx2,i)
569 x3 = xxl(nx3,i)
570 y3 = yyl(nx3,i)
571 area1 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
572 area1 = area1 / area(i)
573 x1 = xxl(nx1,i)
574 y1 = yyl(nx1,i)
575 x2 = xin(ied1,i)
576 y2 = yin(ied1,i)
577 x3 = xin(ied2,i)
578 y3 = yin(ied2,i)
579 area2 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
580 area2 = area2 / area(i)
581 area3 = one - area1 - area2
582 ELSE
583 ! print*,' IMPOSSIBLE CASE'
584 ENDIF
585 ENDIF ! ITRI
586c
587 crklvset(pp1)%AREA(elcrk) = area1
588 crklvset(pp2)%AREA(elcrk) = area2
589 crklvset(pp3)%AREA(elcrk) = area3
590c
591 ENDDO ! IR=1,NEWCRK
592C-------------------
593 nlevset = nlevset + 1 ! update nb of cracks
594C-------------------
595 RETURN
596 END
crklayer3n_ini
subroutine crklayer3n_ini(nel, nft, ilay, nlay, ixtg, elcutc, elcrkini, iel_crktg, inod_crk, iad_crktg, nodenr, dir1, dir2, nodedge, crknodiad, knod2elc, crkedge, xedge3n, ngl, area, xl2, xl3, yl2, yl3)
Definition crklayer3n_ini.F:40
lsint4
subroutine lsint4(y1, z1, y2, z2, y, z, fi)
Definition crklayer4n_adv.F:712
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
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::crkshell
type(xfem_shell_), dimension(:), allocatable crkshell
Definition crackxfem_mod.F:190
crackxfem_mod::crklvset
type(xfem_lvset_), dimension(:), allocatable crklvset
Definition crackxfem_mod.F:146
arret
subroutine arret(nn)
Definition arret.F:86