OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
crklayer3n_ini.F File Reference
#include "implicit_f.inc"
#include "com04_c.inc"
#include "com_xfem1.inc"
#include "comlock.inc"
#include "units_c.inc"

Go to the source code of this file.

Functions/Subroutines

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)

Function/Subroutine Documentation

◆ crklayer3n_ini()

subroutine crklayer3n_ini ( integer nel,
integer nft,
integer ilay,
integer nlay,
integer, dimension(nixtg,*) ixtg,
integer, dimension(2,*) elcutc,
integer, dimension(nlay,*) elcrkini,
integer, dimension(*) iel_crktg,
integer, dimension(*) inod_crk,
integer, dimension(3,*) iad_crktg,
integer, dimension(*) nodenr,
dimension(nlay,nel) dir1,
dimension(nlay,nel) dir2,
integer, dimension(2,*) nodedge,
integer, dimension(*) crknodiad,
integer, dimension(*) knod2elc,
type (xfem_edge_), dimension(*) crkedge,
integer, dimension(3,*) xedge3n,
integer, dimension(nel) ngl,
dimension(nel) area,
dimension(nel) xl2,
dimension(nel) xl3,
dimension(nel) yl2,
dimension(nel) yl3 )

Definition at line 34 of file crklayer3n_ini.F.

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
my_real
#define my_real
Definition cppsort.cpp:32
lsint4
subroutine lsint4(y1, z1, y2, z2, y, z, fi)
Definition crklayer4n_adv.F:712
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
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