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