OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
enrichtg_ini.F File Reference
#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "com_xfem1.inc"
#include "units_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine enrichtg_ini (elbuf_str, ixtg, nft, jft, jlt, nxlay, iad_crktg, iel_crktg, inod_crk, elcutc, nodedge, crknodiad, knod2elc, x, crkedge, xedge3n)

Function/Subroutine Documentation

◆ enrichtg_ini()

subroutine enrichtg_ini ( type (elbuf_struct_) elbuf_str,
integer, dimension(nixtg,*) ixtg,
integer nft,
integer jft,
integer jlt,
integer nxlay,
integer, dimension(3,*) iad_crktg,
integer, dimension(*) iel_crktg,
integer, dimension(*) inod_crk,
integer, dimension(2,*) elcutc,
integer, dimension(2,*) nodedge,
integer, dimension(*) crknodiad,
integer, dimension(*) knod2elc,
x,
type (xfem_edge_), dimension(*) crkedge,
integer, dimension(3,*) xedge3n )

Definition at line 35 of file enrichtg_ini.F.

39C-----------------------------------------------
40 USE crackxfem_mod
41 USE elbufdef_mod
42C-----------------------------------------------
43C I m p l i c i t T y p e s
44C-----------------------------------------------
45#include "implicit_f.inc"
46C-----------------------------------------------
47C G l o b a l P a r a m e t e r s
48C-----------------------------------------------
49#include "mvsiz_p.inc"
50C-----------------------------------------------
51C C o m m o n B l o c k s
52C-----------------------------------------------
53#include "com_xfem1.inc"
54#include "units_c.inc"
55C-----------------------------------------------
56C D u m m y A r g u m e n t s
57C-----------------------------------------------
58 INTEGER IXTG(NIXTG,*),NFT,JFT,JLT,INOD_CRK(*),KNOD2ELC(*),
59 . IAD_CRKTG(3,*),IEL_CRKTG(*),ELCUTC(2,*),NODEDGE(2,*),
60 . CRKNODIAD(*),XEDGE3N(3,*),NXLAY
61 my_real
62 . x(3,*)
63 TYPE (ELBUF_STRUCT_) :: ELBUF_STR
64 TYPE (XFEM_EDGE_) , DIMENSION(*) :: CRKEDGE
65C-----------------------------------------------
66C L o c a l V a r i a b l e s
67C-----------------------------------------------
68 INTEGER I,K,ELCRK,ELCRKTG,ELCUT,pp1,pp2,pp3,IADC(3),IENR0(3),
69 . IENR(3),IED,IEDGE,r,IE10,IE20,IE1,IE2,NOD1,NOD2,N(3),DX(6),
70 . NX(3),dd(3),LAYCUT,ISIGN1,ISIGN2,ISIGN3,SIGBETA,IAD1,IAD2,
71 . ISIGN0(NXEL,3),p1,p2,ICUTEDGE,IBOUNDEDGE,
72 . NTAG(3),EDGEENR(3),ENR(3),ILAY,ITRI,NX1,NX2,NX3,NM,NP
73 my_real
74 . x1g(mvsiz),x2g(mvsiz),x3g(mvsiz),y1g(mvsiz),y2g(mvsiz),
75 . y3g(mvsiz),z1g(mvsiz),z2g(mvsiz),z3g(mvsiz),area(mvsiz),
76 . lxyz0(2),rx(mvsiz),ry(mvsiz),rz(mvsiz),
77 . sx(mvsiz),sy(mvsiz),sz(mvsiz),r11(mvsiz),r12(mvsiz),
78 . r13(mvsiz),r21(mvsiz),r22(mvsiz),r23(mvsiz),r31(mvsiz),
79 . r32(mvsiz),r33(mvsiz),xl1(mvsiz),yl1(mvsiz),xl2(mvsiz),
80 . yl2(mvsiz),xl3(mvsiz),yl3(mvsiz),
81 . fit(3,mvsiz),offg(mvsiz),xin(2,mvsiz),yin(2,mvsiz),xxx,yyy,zzz,fi,
82 . xxl(3,mvsiz),yyl(3,mvsiz),xn(3),yn(3),xm(2),ym(2)
83 my_real x1,y1,x2,y2,x3,y3,x10,y10,z10,x20,y20,z20,
84 . beta,area1,area2,area3
85 DATA dd/2,3,1/
86 DATA dx/1,2,3,1,2,3/
87C=======================================================================
88C COORDONNEES (GLOBAL)
89C-----------------------
90 p1 = 0
91 p2 = 0
92 nx1 = 0
93 nx2 = 0
94 area1 = 0
95 area2 = 0
96 area3 = 0
97 DO i=jft,jlt
98 x1g(i)=x(1,ixtg(2,i+nft))
99 y1g(i)=x(2,ixtg(2,i+nft))
100 z1g(i)=x(3,ixtg(2,i+nft))
101 x2g(i)=x(1,ixtg(3,i+nft))
102 y2g(i)=x(2,ixtg(3,i+nft))
103 z2g(i)=x(3,ixtg(3,i+nft))
104 x3g(i)=x(1,ixtg(4,i+nft))
105 y3g(i)=x(2,ixtg(4,i+nft))
106 z3g(i)=x(3,ixtg(4,i+nft))
107 ENDDO
108C----------------------------
109C LOCAL SYSTEM
110C----------------------------
111 DO i=jft,jlt
112 rx(i) = x2g(i)-x1g(i)
113 ry(i) = y2g(i)-y1g(i)
114 rz(i) = z2g(i)-z1g(i)
115 sx(i) = x3g(i)-x1g(i)
116 sy(i) = y3g(i)-y1g(i)
117 sz(i) = z3g(i)-z1g(i)
118 offg(i) = elbuf_str%GBUF%OFF(i)
119 ENDDO
120 k = 0
121 CALL clskew3(jft,jlt,k,
122 . rx, ry, rz,
123 . sx, sy, sz,
124 . r11,r12,r13,r21,r22,r23,r31,r32,r33,area,offg )
125C--------------------------
126C GLOBAL-->LOCAL TRANSFORMATION
127C--------------------------
128 DO i=jft,jlt
129 xl1(i) = zero
130 yl1(i) = zero
131 xl2(i)=r11(i)*rx(i)+r21(i)*ry(i)+r31(i)*rz(i)
132 yl2(i)=r12(i)*rx(i)+r22(i)*ry(i)+r32(i)*rz(i)
133 xl3(i)=r11(i)*sx(i)+r21(i)*sy(i)+r31(i)*sz(i)
134 yl3(i)=r12(i)*sx(i)+r22(i)*sy(i)+r32(i)*sz(i)
135 area(i)=half*area(i)
136 ENDDO
137C
138C LOCAL TO CENTER OF ELEM (possible warped element)
139C
140 DO i=jft,jlt
141 lxyz0(1)=third*(xl1(i)+xl2(i)+xl3(i))
142 lxyz0(2)=third*(yl1(i)+yl2(i)+yl3(i))
143C
144 xl1(i)=xl1(i)-lxyz0(1)
145 yl1(i)=yl1(i)-lxyz0(2)
146 xl2(i)=xl2(i)-lxyz0(1)
147 yl2(i)=yl2(i)-lxyz0(2)
148 xl3(i)=xl3(i)-lxyz0(1)
149 yl3(i)=yl3(i)-lxyz0(2)
150 ENDDO
151c-----------------------------------------------------
152c Loop over layers
153c-----------------------------------------------------
154 DO ilay=1,nxlay
155 pp1 = nxel*(ilay-1)+1
156 pp2 = pp1 + 1
157 pp3 = pp2 + 1
158C
159 DO i=jft,jlt
160 fit(1,i)=zero
161 fit(2,i)=zero
162 fit(3,i)=zero
163 ENDDO
164C
165 DO i=jft,jlt
166 elcrktg = iel_crktg(i+nft)
167 elcrk = elcrktg + ecrkxfec
168 laycut = crkedge(ilay)%LAYCUT(elcrk)
169 IF (laycut /= 0) THEN
170 xn(1) = xl1(i)
171 yn(1) = yl1(i)
172 xn(2) = xl2(i)
173 yn(2) = yl2(i)
174 xn(3) = xl3(i)
175 yn(3) = yl3(i)
176 xxl(1,i)=xl1(i)
177 yyl(1,i)=yl1(i)
178 xxl(2,i)=xl2(i)
179 yyl(2,i)=yl2(i)
180 xxl(3,i)=xl3(i)
181 yyl(3,i)=yl3(i)
182 DO k=1,3
183 p1 = k
184 p2 = dd(k)
185 ied = crkedge(ilay)%IEDGETG(k,elcrktg)
186 IF (ied > 0) THEN
187 iedge = xedge3n(k,elcrktg)
188 beta = crkedge(ilay)%RATIO(iedge)
189 nod1 = nodedge(1,iedge)
190 nod2 = nodedge(2,iedge)
191 IF (nod1 == ixtg(k+1,i+nft) .and. nod2 == ixtg(dd(k)+1,i+nft)) THEN
192 p1 = k
193 p2 = dd(k)
194 ELSEIF (nod2 == ixtg(k+1,i+nft).and.nod1 == ixtg(dd(k)+1,i+nft))THEN
195 p1 = dd(k)
196 p2 = k
197 ENDIF
198 xin(ied,i) = xn(p1) + beta*(xn(p2) - xn(p1))
199 yin(ied,i) = yn(p1) + beta*(yn(p2) - yn(p1))
200 xm(ied) = half*(xn(p1)+xn(p2))
201 ym(ied) = half*(yn(p1)+yn(p2))
202 ENDIF
203 ENDDO
204C
205 DO k=1,3
206 fi = zero
207 CALL lsint4(xm(1),ym(1),xm(2),ym(2),xn(k),yn(k),fi )
208 IF (fit(k,i)==zero) fit(k,i) = fi
209 ENDDO
210 ENDIF
211 ENDDO
212C
213 IF (ilay == 1) THEN
214 DO i=jft,jlt
215 elcrktg = iel_crktg(i+nft)
216 elcrk = elcrktg + ecrkxfec
217 elcut = crkedge(ilay)%LAYCUT(elcrk)
218 IF (elcut /= 0) THEN
219 elcutc(1,i+nft) = 1
220 elcutc(2,i+nft) = 1
221 ENDIF
222 ENDDO
223 ENDIF
224c--------------------------------------------------------------------
225 DO i=jft,jlt
226 elcrktg = iel_crktg(i+nft)
227 elcrk = elcrktg + ecrkxfec
228 laycut = crkedge(ilay)%LAYCUT(elcrk)
229 IF (laycut /= 0) THEN
230c
231 iadc(1) = iad_crktg(1,elcrktg)
232 iadc(2) = iad_crktg(2,elcrktg)
233 iadc(3) = iad_crktg(3,elcrktg)
234C
235 ienr0(1) = crknodiad(iadc(1))
236 ienr0(2) = crknodiad(iadc(2))
237 ienr0(3) = crknodiad(iadc(3))
238C
239 n(1) = ixtg(2,i+nft)
240 n(2) = ixtg(3,i+nft)
241 n(3) = ixtg(4,i+nft)
242C
243 nx(1) = inod_crk(n(1))
244 nx(2) = inod_crk(n(2))
245 nx(3) = inod_crk(n(3))
246C
247 ienr(1) = ienr0(1) + knod2elc(nx(1))*(ilay-1)
248 ienr(2) = ienr0(2) + knod2elc(nx(2))*(ilay-1)
249 ienr(3) = ienr0(3) + knod2elc(nx(3))*(ilay-1)
250C
251 ntag(1:3) = 0
252 edgeenr(1:3) = 0
253 enr(1:3) = 0
254C
255 DO k=1,3
256 ied = crkedge(ilay)%IEDGETG(k,elcrktg)
257 IF (ied > 0) THEN
258 ntag(k) = ntag(k) + 1
259 ntag(dd(k)) = ntag(dd(k)) + 1
260 iedge = xedge3n(k,elcrktg)
261 nod1 = nodedge(1,iedge)
262 nod2 = nodedge(2,iedge)
263 ie10 = crkedge(ilay)%EDGEENR(1,iedge)
264 ie20 = crkedge(ilay)%EDGEENR(2,iedge)
265 IF (nod1 == n(k) .and. nod2 == n(dd(k))) THEN
266 p1 = k
267 p2 = dd(k)
268 ELSEIF (nod2 == n(k) .and. nod1 == n(dd(k))) THEN
269 p1 = dd(k)
270 p2 = k
271 ENDIF
272 edgeenr(p1) = ie10
273 edgeenr(p2) = ie20
274 ENDIF
275 ENDDO
276C
277 DO k=1,3
278 IF (ntag(k) > 0) THEN
279 enr(k) = edgeenr(k)
280 ELSE
281 enr(k) = ienr(k)
282 ENDIF
283 ENDDO
284C
285 DO k=1,3
286 IF (ienr(k) > ienrnod) THEN
287 WRITE(iout,*) 'ERROR CRACK INITIATION,ENRICHMENT NODE EXCEEDED'
288 CALL arret(2)
289 ENDIF
290 ENDDO
291C
292 isign1 = int(sign(one,fit(1,i)))
293 isign2 = int(sign(one,fit(2,i)))
294 isign3 = int(sign(one,fit(3,i)))
295C
296 IF (fit(1,i) == zero) isign1 = 0
297 IF (fit(2,i) == zero) isign2 = 0
298 IF (fit(3,i) == zero) isign3 = 0
299C
300 DO k=1,nxel
301 isign0(k,1) = isign1
302 isign0(k,2) = isign2
303 isign0(k,3) = isign3
304 ENDDO
305c------------------------------------------
306 itri = 0
307 nm = 0
308 np = 0
309 DO k=1,3
310 IF (isign0(1,k) > 0) THEN
311 itri = itri + 1
312 np = k
313 ELSEIF (isign0(1,k) < 0) THEN
314 nm = k
315 ENDIF
316 ENDDO
317 IF (itri == 1) THEN
318 itri = -1
319 nx1 = np
320 ELSEIF (itri == 2) THEN
321 itri = 1
322 nx1 = nm
323 ENDIF
324 nx2 = dx(nx1+1)
325 nx3 = dx(nx2+1)
326 xfem_phantom(ilay)%ITRI(1,elcrk) = itri
327 xfem_phantom(ilay)%ITRI(2,elcrk) = nx1
328C
329 xfem_phantom(ilay)%IFI(iadc(1)) = isign0(1,1)
330 xfem_phantom(ilay)%IFI(iadc(2)) = isign0(1,2)
331 xfem_phantom(ilay)%IFI(iadc(3)) = isign0(1,3)
332C-------------------
333c IXEL = 1 => positif element within ILAY (ILEV = PP1)
334C-------------------
335 crklvset(pp1)%ENR0(1,iadc(1)) = enr(1)
336 crklvset(pp1)%ENR0(1,iadc(2)) = enr(2)
337 crklvset(pp1)%ENR0(1,iadc(3)) = enr(3)
338C
339 IF (isign0(1,1) > 0) crklvset(pp1)%ENR0(1,iadc(1)) = 0
340 IF (isign0(1,2) > 0) crklvset(pp1)%ENR0(1,iadc(2)) = 0
341 IF (isign0(1,3) > 0) crklvset(pp1)%ENR0(1,iadc(3)) = 0
342C
343 crklvset(pp1)%ENR0(2,iadc(1)) = crklvset(pp1)%ENR0(1,iadc(1))
344 crklvset(pp1)%ENR0(2,iadc(2)) = crklvset(pp1)%ENR0(1,iadc(2))
345 crklvset(pp1)%ENR0(2,iadc(3)) = crklvset(pp1)%ENR0(1,iadc(3))
346C-------------------
347c IXEL = 2 => negatif element within ILAY (ILEV = PP2)
348C-------------------
349 crklvset(pp2)%ENR0(1,iadc(1)) = enr(1)
350 crklvset(pp2)%ENR0(1,iadc(2)) = enr(2)
351 crklvset(pp2)%ENR0(1,iadc(3)) = enr(3)
352C
353 IF (isign0(2,1) < 0) crklvset(pp2)%ENR0(1,iadc(1)) = 0
354 IF (isign0(2,2) < 0) crklvset(pp2)%ENR0(1,iadc(2)) = 0
355 IF (isign0(2,3) < 0) crklvset(pp2)%ENR0(1,iadc(3)) = 0
356C
357 crklvset(pp2)%ENR0(2,iadc(1)) = crklvset(pp2)%ENR0(1,iadc(1))
358 crklvset(pp2)%ENR0(2,iadc(2)) = crklvset(pp2)%ENR0(1,iadc(2))
359 crklvset(pp2)%ENR0(2,iadc(3)) = crklvset(pp2)%ENR0(1,iadc(3))
360C-------------------
361c IXEL = 3 => Third phantom element (ILEV = PP3)
362C-------------------
363 IF (itri < 0) THEN ! sign ILEV3 = ILEV2 < 0
364 ie1 = xedge3n(nx1,elcrktg)
365 ie2 = xedge3n(nx3,elcrktg)
366c IF (CRKEDGE(ILAY)%ICUTEDGE(IE2) == 2) THEN ! tip edge
367 crklvset(pp3)%ENR0(1,iadc(nx1)) = abs(crklvset(pp2)%ENR0(1,iadc(nx1)))
368 crklvset(pp3)%ENR0(1,iadc(nx2)) = crklvset(pp2)%ENR0(1,iadc(nx2))
369 crklvset(pp3)%ENR0(1,iadc(nx3)) = crklvset(pp2)%ENR0(1,iadc(nx3))
370 crklvset(pp2)%ENR0(1,iadc(nx1)) = -crknodiad(iadc(nx1)) - knod2elc(nx(nx1))*(ilay-1)
371c
372 beta = crkedge(ilay)%RATIO(ie2)
373 nod1 = nodedge(1,ie2)
374 nod2 = nodedge(2,ie2)
375 k = nx3 ! starting edge node(local)
376 IF (nod1 == ixtg(k+1,i+nft) .and.
377 . nod2 == ixtg(dd(k)+1,i+nft)) THEN
378 p1 = k
379 p2 = dd(k)
380 sigbeta = ie2
381 ELSEIF (nod2 == ixtg(k+1,i+nft) .and.
382 . nod1 == ixtg(dd(k)+1,i+nft)) THEN
383 p1 = dd(k)
384 p2 = k
385 sigbeta = -ie2
386 ENDIF
387c
388c-- AREA factors
389 x1 = xxl(nx1,i)
390 y1 = yyl(nx1,i)
391 ied = crkedge(ilay)%IEDGETG(nx1,elcrktg)
392 x2 = xin(ied,i)
393 y2 = yin(ied,i)
394 ied = crkedge(ilay)%IEDGETG(nx3,elcrktg)
395 x3 = xin(ied,i)
396 y3 = yin(ied,i)
397 area1 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
398 area1 = area1 / area(i)
399 x1 = xxl(nx2,i)
400 y1 = yyl(nx2,i)
401 x2 = xxl(nx3,i)
402 y2 = yyl(nx3,i)
403 area2 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
404 area3 = one - area1 - area2
405c
406 ELSEIF (itri > 0) THEN ! sign ILEV3 = ILEV1 > 0
407c
408 ie1 = xedge3n(nx1,elcrktg) ! tip edge
409c IF (CRKEDGE(ILAY)%ICUTEDGE(IE1) == 2) THEN
410 crklvset(pp3)%ENR0(1,iadc(nx1)) = abs(crklvset(pp1)%ENR0(1,iadc(nx1)))
411 crklvset(pp3)%ENR0(1,iadc(nx2)) = crklvset(pp1)%ENR0(1,iadc(nx2))
412 crklvset(pp3)%ENR0(1,iadc(nx3)) = crklvset(pp1)%ENR0(1,iadc(nx3))
413 crklvset(pp1)%ENR0(1,iadc(nx1)) = -crknodiad(iadc(nx1)) - knod2elc(nx(nx1))*(ilay-1)
414c
415 beta = crkedge(ilay)%RATIO(ie1)
416 nod1 = nodedge(1,ie1)
417 nod2 = nodedge(2,ie1)
418 k = nx1 ! starting edge node (local)
419 IF (nod1 == ixtg(k+1,i+nft) .and.
420 . nod2 == ixtg(dd(k)+1,i+nft)) THEN
421 p1 = k
422 p2 = dd(k)
423 sigbeta = ie1
424 ELSEIF (nod2 == ixtg(k+1,i+nft) .and.
425 . nod1 == ixtg(dd(k)+1,i+nft)) THEN
426 p1 = dd(k)
427 p2 = k
428 sigbeta = -ie1
429 ENDIF
430c
431c-- AREA factors for each phantom
432 x1 = xxl(nx1,i)
433 y1 = yyl(nx1,i)
434 ied = crkedge(ilay)%IEDGETG(nx3,elcrktg)
435 x2 = xin(ied,i)
436 y2 = yin(ied,i)
437 ied = crkedge(ilay)%IEDGETG(nx1,elcrktg)
438 x3 = xin(ied,i)
439 y3 = yin(ied,i)
440 area1 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
441 area1 = area1 / area(i)
442 x1 = xxl(nx2,i)
443 y1 = yyl(nx2,i)
444 x3 = xxl(nx3,i)
445 y3 = yyl(nx3,i)
446 area2 = half*abs((x1-x3)*(y2-y1) - (x1-x2)*(y3-y1))
447 area3 = one - area1 - area2
448c
449 ENDIF ! ITRI
450c
451 crklvset(pp3)%ENR0(2,iadc(1)) = crklvset(pp3)%ENR0(1,iadc(1))
452 crklvset(pp3)%ENR0(2,iadc(2)) = crklvset(pp3)%ENR0(1,iadc(2))
453 crklvset(pp3)%ENR0(2,iadc(3)) = crklvset(pp3)%ENR0(1,iadc(3))
454c
455 crklvset(pp1)%AREA(elcrk) = area1
456 crklvset(pp2)%AREA(elcrk) = area2
457 crklvset(pp3)%AREA(elcrk) = area3
458C-------------------
459 ENDIF ! LAYCUT /= 0
460 ENDDO ! DO I=JFT,JLT
461 ENDDO ! DO ILAY=1,NXLAY
462C-------------
463 RETURN
clskew3
subroutine clskew3(jft, jlt, irep, rx, ry, rz, sx, sy, sz, e1x, e2x, e3x, e1y, e2y, e3y, e1z, e2z, e3z, det)
Definition clskew.F:34
my_real
#define my_real
Definition cppsort.cpp:32
lsint4
subroutine lsint4(y1, z1, y2, z2, y, z, fi)
Definition crklayer4n_adv.F:710
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
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::crklvset
type(xfem_lvset_), dimension(:), allocatable crklvset
Definition crackxfem_mod.F:146
arret
subroutine arret(nn)
Definition arret.F:87