OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_wind_xfem.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!|| fail_wind_xfem ../engine/source/materials/fail/alter/fail_wind_xfem.F
25!||--- called by ------------------------------------------------------
26!|| mulawc ../engine/source/materials/mat_share/mulawc.F90
27!|| usermat_shell ../engine/source/materials/mat_share/usermat_shell.F
28!||====================================================================
29 SUBROUTINE fail_wind_xfem(
30 1 NEL ,NUPARAM ,NUVAR ,UPARAM ,UVAR ,
31 2 TIME ,TIMESTEP ,SSP ,ALDT ,CRKLEN ,
32 3 ELCRKINI ,TDEL ,OFF ,OFFLY ,
33 4 SIGNXX ,SIGNYY ,SIGNXY ,SIGNYZ ,SIGNZX ,
34 5 NGL ,IXEL ,ILAY ,IPT ,NPT ,
35 6 IXFEM ,IORTH ,DIR_A ,DIR1_CRK ,DIR2_CRK ,
36 7 CRKDIR )
37C-----------------------------------------------
38c Windshield failure model for XFEM (ref : PhD thesis Christian Alter 2018)
39c now : only for monolayer xfem
40C-----------------------------------------------
41C I m p l i c i t T y p e s
42C-----------------------------------------------
43#include "implicit_f.inc"
44C-----------------------------------------------
45C G l o b a l P a r a m e t e r s
46C-----------------------------------------------
47#include "units_c.inc"
48#include "comlock.inc"
49#include "com_xfem1.inc"
50#include "mvsiz_p.inc"
51C---------+---------+---+---+--------------------------------------------
52C VAR | SIZE |TYP| RW| DEFINITION
53C---------+---------+---+---+--------------------------------------------
54C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
55C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
56C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
57C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
58C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
59C---------+---------+---+---+--------------------------------------------
60C TIME | 1 | F | R | CURRENT TIME
61C TIMESTEP| 1 | F | R | CURRENT TIME STEP
62C---------+---------+---+---+--------------------------------------------
63C SIGNXX | NEL | F | R | STRESS XX
64C SIGNYY | NEL | F | R | STRESS YY
65C ... | | | |
66C---------+---------+---+---+--------------------------------------------
67C OFF | NEL | F |R/W| DELETED ELEMENT FLAG (=1. ON, =0. OFF)
68C---------+---------+---+---+--------------------------------------------
69C ILAY CURRENT LAYER
70C IPT CURRENT INTEGRATION POINT IN THE LAYER (FOR OUTPUT ONLY)
71C---------+---------+---+---+--------------------------------------------
72C I N P U T A r g u m e n t s
73C-----------------------------------------------
74 INTEGER NEL,NUPARAM,NUVAR,IXFEM,IORTH,IXEL,ILAY,IPT,NPT
75 INTEGER ELCRKINI(NXLAYMAX,NEL),NGL(NEL)
76 my_real TIME,TIMESTEP
77 my_real, DIMENSION(NEL) :: SSP,ALDT,TDEL,CRKLEN,
78 . SIGNXX,SIGNYY,SIGNZZ,SIGNXY,SIGNYZ,SIGNZX
79 my_real UPARAM(NUPARAM),DIR_A(NEL,2),CRKDIR(NEL,2),DIR1_CRK(NXLAYMAX,MVSIZ),
80 . DIR2_CRK(NXLAYMAX,NEL)
81 TARGET :: dir1_crk,dir2_crk
82C-----------------------------------------------
83C I N P U T O U T P U T A r g u m e n t s
84C-----------------------------------------------
85 INTEGER OFFLY(NEL)
86 my_real uvar(nel,nuvar),off(nel)
87C-----------------------------------------------
88C L o c a l V a r i a b l e s
89C-----------------------------------------------
90 INTEGER I,J,NINDX,NINDXD,IDEB,ISRATE
91 INTEGER , DIMENSION(NEL) :: INDX,INDXD,RFLAG
92 my_real sigmax,aa,bb,cc,s1,s2,s3,cr,k_ic,k_th,v0,vc,tfact,reflen,
93 . cosi,cosx,cosy,sinx,siny,cos2,sin2,dmg1,dmg3,
94 . sp1,sp2,sp3,geored,cr_foil,cr_air,cr_core,cr_edge,alpha,alphai,exp_n,exp_m
95 my_real, DIMENSION(NEL) :: dadv,tpropg,ai,formf,sigp_akt,dmg_scale,
96 . sigp1,sigp2,sigdt1,sigdt2,sig_dtf1,sig_dtf2,sigp_min,sigp_max
97 my_real, DIMENSION(:), POINTER :: dir1,dir2
98 CHARACTER (LEN=3) :: XCHAR
99c--------------------------------------------------------------------
100c! cm(6): initial crack depth in contact to foil
101c! cm(7): initial crack depth in contact to air
102c! cm(8): initial crack depth on edge
103c! cm(9): K_IC - fracture toughness from literature
104c! cm(10): K_TH - fatigue threshold from literature
105c--------------------------------------------------------------------
106c UVAR(1) = SIGP1 ! principal stress DIR1 ! uvar 17
107c UVAR(1) = SIGP2 ! principal stress DIR2 ! uvar 18
108c UVAR(3) = SIGDT1 ! principal stress rate DIR1
109c UVAR(4) = SIGDT2 ! principal stress rate DIR2
110c UVAR(5) = SIGP_MIN ! uvar 24
111c UVAR(6) = SIGP_MAX ! uvar 14
112c UVAR(7) = FORMF ! uvar 25
113c UVAR(8) = ai ! uvar 26
114c UVAR(9) = sigp_akt
115c UVAR(10) = edge element flag
116c UVAR(11) = reduction factor flag
117c UVAR(12) = DADV
118c UVAR(13) = SIG_DTF1
119c UVAR(14) = SIG_DTF2
120C=======================================================================
121 xchar = ' '
122 IF (ixel > 0) THEN ! testing phantom elements
123 IF (ixel == 1) THEN
124 xchar = '1st'
125 ELSEIF (ixel == 2) THEN
126 xchar = '2nd'
127 ELSEIF (ixel == 3) THEN
128 xchar = '3rd'
129 ENDIF
130 ELSE
131 xchar = 'standard'
132 ENDIF
133c
134 rflag(1:nel) = 0
135 indx(1:nel) = 0
136 nindx = 0
137 nindxd = 0
138c--------------------------------------------------------------------
139 exp_n = uparam(1)
140 cr_foil = uparam(2)
141 cr_air = uparam(3)
142 cr_core = uparam(4)
143 cr_edge = uparam(5)
144 k_ic = uparam(6)
145 k_th = uparam(7)
146 v0 = uparam(8)
147 vc = uparam(9)
148 alpha = uparam(10)
149 geored = uparam(11)
150 reflen = nint(uparam(14))
151 israte = nint(uparam(16))
152 ideb = nint(uparam(17))
153c
154 alphai = one - alpha
155 exp_m = one / (one + exp_n)
156c---
157c IF (IXEL == 0) THEN
158 DO i = 1,nel
159 dadv(i) = min(one, geored / sqrt(aldt(i)/reflen) ) ! reduction factor for advancement
160 tpropg(i) = crklen(i) / min(vc,ssp(i))
161 tpropg(i) = max(tpropg(i), timestep)
162 uvar(i,12)= dadv(i)
163 END DO
164c ENDIF
165c--------------------------------------------------------------------
166 IF (time == zero) then
167c Crack depth modification and definition of Y (formfactor)
168 DO i=1,nel
169 IF (uvar(i,10) == one) THEN ! edge-element
170 ai(i) = cr_edge
171 formf(i) = 1.12
172 ELSE
173 formf(i) = one
174 IF (ipt == 1) THEN
175 ai(i) = cr_foil ! crack depth at foil side in mm (unit system: length)
176 ELSEIF (ipt == npt) THEN
177 ai(i) = cr_air ! crack depth at top side in mm (unit system: length)
178 ELSE
179 ai(i) = cr_core ! crack depth inside glass
180 ENDIF
181 ENDIF
182 sigp_akt(i) = (two*(exp_n + one)*k_ic**exp_n) /
183 . ((exp_n - two)*v0*(formf(i)*sqrt(pi))**exp_n*
184 . ai(i)**((exp_n-two)/two))
185 sigp_akt(i) = sigp_akt(i)**exp_m
186c
187 sigp_min(i) = k_th / (formf(i)*sqrt(pi*ai(i)))
188 sigp_max(i) = k_ic / (formf(i)*sqrt(pi*ai(i)))
189 uvar(i,5) = sigp_min(i)
190 uvar(i,6) = sigp_max(i)
191 uvar(i,7) = formf(i)
192 uvar(i,8) = ai(i)
193 uvar(i,9) = sigp_akt(i)
194 ENDDO
195 ELSE
196 sigp_min(1:nel) = uvar(1:nel,5)
197 sigp_max(1:nel) = uvar(1:nel,6)
198 formf(1:nel) = uvar(1:nel,7)
199 ai(1:nel) = uvar(1:nel,8)
200 sigp_akt(1:nel) = uvar(1:nel,9)
201 ENDIF
202c--------------------------------------------------------------------
203c Principal stress
204c--------------------------------------------------------------------
205 DO i = 1,nel
206 aa = (signxx(i) + signyy(i))*half
207 bb = (signxx(i) - signyy(i))*half
208 cc = signxy(i)
209 cr = sqrt(bb**2 + cc**2)
210 s1 = aa + cr
211 s2 = aa - cr
212 sigp1(i) = s1
213 sigp2(i) = s2
214c
215c IF (S12 == ZERO .and. S2 >= ZERO) THEN
216c MAJANG = ZERO
217c ELSEIF (S12 == ZERO .and. S2 < ZERO) THEN
218c MAJANG = HALF*PI
219c ELSE
220c MAJANG = ATAN((SIGP1(I) - SIGNXX(I)) / S12)
221c ENDIF
222 END DO
223c--------------------------------------------------------------------
224c Calculate stress rate + filtering
225c--------------------------------------------------------------------
226 IF (israte == 0) THEN
227c exponential moving average with smoohting coefficient = alpha
228 DO i = 1,nel
229 sigdt1(i) = abs(sigp1(i) - uvar(i,1)) / max(timestep, em20)
230 sigdt2(i) = abs(sigp2(i) - uvar(i,2)) / max(timestep, em20)
231 sigdt1(i) = alpha * sigdt1(i) + alphai * uvar(i,3)
232 sigdt2(i) = alpha * sigdt2(i) + alphai * uvar(i,4)
233 END DO
234 ELSE
235c arythmetic moving average over 50 cycles
236 DO i = 1,nel
237 IF (off(i) == one) THEN
238 sigdt1(i) = abs(sigp1(i) - uvar(i,1)) / max(timestep, em20)
239 uvar(i,21) = sigdt1(i)
240 DO j = 1,49
241 sigdt1(i) = sigdt1(i) + uvar(i,21+j)
242 ENDDO
243 sigdt1(i) = sigdt1(i) / 50
244 DO j = 70,22,-1
245 uvar(i,j) = uvar(i,j-1)
246 ENDDO
247c
248 sigdt2(i) = abs(sigp2(i) - uvar(i,2)) / max(timestep, em20)
249 uvar(i,71) = sigdt2(i)
250 DO j = 1,49
251 sigdt2(i) = sigdt2(i) + uvar(i,71+j)
252 ENDDO
253 sigdt2(i) = sigdt2(i) / 50
254 DO j = 120,72,-1
255 uvar(i,j) = uvar(i,j-1)
256 ENDDO
257 ELSE
258 sigdt1(i) = zero
259 sigdt2(i) = zero
260 ENDIF
261 END DO
262 ENDIF
263c--------------------------------------------------------------------
264c Max strength calculation
265c--------------------------------------------------------------------
266 IF (ipt == npt) THEN ! Top layer integration point
267 DO i = 1,nel
268 sig_dtf1(i) = sigp_akt(i) *abs(sigdt1(i))**exp_m
269 sig_dtf2(i) = sigp_akt(i) *abs(sigdt2(i))**exp_m
270 sig_dtf1(i) = max(sig_dtf1(i),sigp_min(i))
271 sig_dtf1(i) = min(sig_dtf1(i),sigp_max(i))
272 sig_dtf2(i) = max(sig_dtf2(i),sigp_min(i))
273 sig_dtf2(i) = min(sig_dtf2(i),sigp_max(i))
274 END DO
275 ELSEIF (ipt < npt) THEN ! inner integration points
276 DO i = 1,nel
277 IF (uvar(i,10) == one) THEN ! edge element => stress rate dependent
278 sig_dtf1(i) = sigp_akt(i)*abs(sigdt1(i))**exp_m
279 sig_dtf2(i) = sigp_akt(i)*abs(sigdt2(i))**exp_m
280 sig_dtf1(i) = max(sig_dtf1(i),sigp_min(i))
281 sig_dtf1(i) = min(sig_dtf1(i),sigp_max(i))
282 sig_dtf2(i) = max(sig_dtf2(i),sigp_min(i))
283 sig_dtf2(i) = min(sig_dtf2(i),sigp_max(i))
284 ELSE
285 sig_dtf1(i) = sigp_max(i)
286 sig_dtf2(i) = sigp_max(i)
287 ENDIF
288 END DO
289 ENDIF
290c--------------------------
291 IF (ixel == 0) THEN ! original element
292 DO i = 1,nel
293 IF (off(i) == one) THEN
294 IF (elcrkini(ilay,i) == 0) THEN
295 ! no cracks in the neibourhood => initialization
296 IF (sigp1(i) > sig_dtf1(i) .and. tdel(i) == zero .and. offly(i) == 1) THEN
297 ! start damage
298 tdel(i) = time
299 offly(i) = 0
300 elcrkini(ilay,i) = -5 ! for crack length estimation
301 nindx = nindx+1
302 indx(nindx) = i
303 rflag(i) = 1
304 ENDIF
305 IF (tdel(i) > 0) THEN ! damage started already
306 tfact = (time - tdel(i)) / tpropg(i)
307 tfact = min(one, tfact)
308 dmg_scale(i) = one - tfact
309 nindxd = nindxd+1
310 indxd(nindxd) = i
311 IF (tfact >= one) THEN ! end of damage in normal direction
312 elcrkini(ilay,i) = -1 ! failure completed (by initiation)
313 off(i) = four_over_5
314 IF (iorth > 0) THEN
315 cosx = dir_a(i,1)
316 sinx = dir_a(i,2)
317 cosy = crkdir(i,1)
318 siny = crkdir(i,2)
319 dir1_crk(ilay,i) = cosx*cosy - sinx*siny
320 dir2_crk(ilay,i) = cosx*siny + sinx*cosy
321 ELSE
322 dir1_crk(ilay,i) = crkdir(i,1)
323 dir2_crk(ilay,i) = crkdir(i,2)
324 ENDIF
325 ELSE
326 ENDIF
327 ENDIF
328c
329 ELSEIF (elcrkini(ilay,i) == 2) THEN
330 ! neighbor element cut at common edge => advancement
331 uvar(i,11) = one ! reduction factor flag for output
332 sig_dtf1(i) = sig_dtf1(i)*dadv(i) ! apply reduction factor to criteria
333 sig_dtf2(i) = sig_dtf2(i)*dadv(i)
334 sig_dtf1(i) = max(sig_dtf1(i),sigp_min(i)*dadv(i))
335 sig_dtf1(i) = min(sig_dtf1(i),sigp_max(i)*dadv(i))
336 sig_dtf2(i) = max(sig_dtf2(i),sigp_min(i)*dadv(i))
337 sig_dtf2(i) = min(sig_dtf2(i),sigp_max(i)*dadv(i))
338 IF (sigp1(i) > sig_dtf1(i) .and. tdel(i) == zero .and. offly(i) == 1) THEN
339 tdel(i) = time
340 offly(i) = 0
341 elcrkini(ilay,i) = 5 ! tag for crack length estimation
342 nindx = nindx+1
343 indx(nindx) = i
344 rflag(i) =-1
345 ENDIF
346 IF (offly(i) == 0) THEN ! damage started already
347 tfact = (time - tdel(i)) / tpropg(i)
348 nindxd = nindxd+1
349 indxd(nindxd) = i
350 IF (tfact >= one) THEN
351 elcrkini(ilay,i) = 1 ! advance existing crack
352 dmg_scale(i) = zero
353 off(i) = four_over_5
354 IF (iorth > 0) THEN
355 cosx = dir_a(i,1)
356 sinx = dir_a(i,2)
357 cosy = crkdir(i,1)
358 siny = crkdir(i,2)
359 dir1_crk(ilay,i) = cosx*cosy - sinx*siny
360 dir2_crk(ilay,i) = cosx*siny + sinx*cosy
361 ELSE
362 dir1_crk(ilay,i) = crkdir(i,1)
363 dir2_crk(ilay,i) = crkdir(i,2)
364 ENDIF
365 ELSE
366 dmg_scale(i) = one - tfact
367 ENDIF
368 ENDIF
369 ENDIF
370 ENDIF ! OFF == 1
371 ENDDO
372c--------------------------
373c calculate and save principal stress direction when first failure occurs
374c--------------------------
375 IF (nindx > 0) THEN
376 DO j=1,nindx
377 i = indx(j)
378 s1 = signxy(i)
379 s2 = sigp1(i) - signxx(i)
380 cr = sqrt(s1**2 + s2**2)
381 IF (cr > zero) THEN
382 crkdir(i,1) = s1 / cr
383 crkdir(i,2) = s2 / cr
384 ELSEIF (s1 > s2) THEN
385 crkdir(i,1) = zero
386 crkdir(i,2) = one
387 ELSE
388 crkdir(i,1) = one
389 crkdir(i,2) = zero
390 ENDIF
391 ENDDO
392c
393 IF (iorth > 0) THEN ! stress is in orthotropic frame (DIR_A)
394 DO j=1,nindx
395 i = indx(j)
396 cosx = dir_a(i,1)
397 sinx = dir_a(i,2)
398 cosy = crkdir(i,1)
399 siny = crkdir(i,2)
400 dir1_crk(ilay,i) = cosx*cosy - sinx*siny
401 dir2_crk(ilay,i) = cosx*siny + sinx*cosy
402 ENDDO
403 ELSE ! stress is in local element frame
404 DO j=1,nindx
405 i = indx(j)
406 dir1_crk(ilay,i) = crkdir(i,1)
407 dir2_crk(ilay,i) = crkdir(i,2)
408 ENDDO
409 ENDIF
410 ENDIF
411c
412c apply progressive damage and turn the stress tensor back to the local system
413c
414 IF (nindxd > 0) THEN
415 DO j=1,nindxd
416 i = indxd(j)
417 s1 = signxx(i)
418 s2 = signyy(i)
419 s3 = signxy(i)
420 cosx = crkdir(i,1)
421 sinx = crkdir(i,2)
422 cos2 = cosx * cosx
423 sin2 = sinx * sinx
424 cosi = cosx * sinx
425c rotate stress to previously saved crack direction
426 sp1 = cos2*s1 + sin2*s2 + two*cosi*s3
427 sp2 = sin2*s1 + cos2*s2 - two*cosi*s3
428 sp3 = cosi*(s2 - s1) + (cos2 - sin2)*s3
429
430 dmg1 = dmg_scale(i)
431 dmg3 = min(one, 0.6 + 0.4 * dmg1)
432c stress reduction
433 sp1 = sp1 * dmg1
434 sp3 = sp3 * dmg3
435c rotate reduced stress back to current element coordinate system
436 signxx(i) = cos2*sp1 + sin2*sp2 - two*cosi*sp3
437 signyy(i) = sin2*sp1 + cos2*sp2 + two*cosi*sp3
438 signxy(i) = cosi*(sp1 - sp2) + (cos2 - sin2)*sp3
439 ENDDO
440c
441 if (ideb==1) then
442 DO j=1,nindxd
443 i = indxd(j)
444 write(iout,'(A,I10,3E16.9)') 'DAMAGE ELEMENT N, CRLEN,TPROPG,DMG_SCALE=',
445 . ngl(i),crklen(i),tpropg(i),dmg_scale(i)
446 ENDDO
447 endif
448 ENDIF
449c
450 ELSEIF (ixel > 0) THEN ! phantom element
451 DO i = 1,nel
452 IF (off(i) == one) THEN
453 IF (elcrkini(ilay,i) == 2) THEN ! second crack arrived => delete phantom element
454 nindx = nindx+1
455 indx(nindx) = i
456 rflag(i) = 3
457 off(i) = four_over_5
458 ELSE
459 sigp_min(i) = k_th / (sqrt(pi*cr_foil))
460 sigp_max(i) = k_ic / (sqrt(pi*cr_foil))
461 sig_dtf1(i) = sigp_max(i) * dadv(i)
462 IF (sigp1(i) > sig_dtf1(i)) THEN
463 ! phantom element reaches criteria
464c delete phantom but should transfer information to remaining neighbor element
465 nindx = nindx+1
466 indx(nindx) = i
467 rflag(i) = 3
468 off(i) = four_over_5
469 ENDIF
470 ENDIF
471 ENDIF
472 ENDDO
473 ENDIF
474c-----------------------------------------------
475 DO i=1,nel
476 uvar(i,1) = sigp1(i)
477 uvar(i,2) = sigp2(i)
478 uvar(i,3) = sigdt1(i)
479 uvar(i,4) = sigdt2(i)
480 END DO
481c-----------------------------------------------
482 IF (nindx > 0) THEN
483 DO j=1,nindx
484 i = indx(j)
485#include "lockon.inc"
486 WRITE(iout, 3000)ngl(i),ilay,ipt
487 WRITE(istdo,3100)ngl(i),ilay,ipt,time
488c initialization std element
489 IF (rflag(i) == 1) WRITE(iout ,4000)
490 IF (rflag(i) == 1) WRITE(istdo,4000)
491c advancement std element
492 IF (rflag(i) == -1) WRITE(iout, 4200)
493 IF (rflag(i) == -1) WRITE(istdo,4200)
494 if (ideb==1 .and. abs(rflag(i))==1) then
495 write(iout,'(A,5E16.9)') 'SIGP1,SIG_DTF1,SIGP_AKT,SIGP_MAX,SIGDT1=',
496 . sigp1(i),sig_dtf1(i),sigp_akt(i),sigp_max(i),sigdt1(i)
497 write(iout,'(A,5E16.9)') 'SIGP2,SIG_DTF2,SIGP_AKT,SIGP_MIN,SIGDT2=',
498 . sigp2(i),sig_dtf2(i),sigp_akt(i),sigp_min(i),sigdt2(i)
499 if (rflag(i)==-1) then
500 write(iout,'(A,E16.9)') 'advancement => crit reduction, DADV =',dadv(i)
501 endif
502 endif
503c delete phantom
504 IF (rflag(i) == 3) WRITE(iout, 4400) xchar,ngl(i)
505 IF (rflag(i) == 3) WRITE(istdo,4500) xchar,ngl(i),time
506#include "lockoff.inc"
507 ENDDO
508 ENDIF
509c-----------------------------------------------
510 3000 FORMAT(1x,'FAILURE IN SHELL',i10,1x,'LAYER',i2,1x,'INT POINT',i2)
511 3100 FORMAT(1x,'FAILURE IN SHELL',i10,1x,'LAYER',i2,1x,'INT POINT',i2,
512 . 1x,'AT TIME ',1pe12.4)
513 4000 FORMAT(10x,'CRACK INITIALIZATION')
514 4200 FORMAT(10x,'CRACK ADVANCEMENT')
515 4400 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10)
516 4500 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10,/
517 . 1x,'AT TIME :',1pe12.4)
518c-----------------------------------------------
519 RETURN
520 END
alpha
#define alpha
Definition eval.h:35
fail_wind_xfem
subroutine fail_wind_xfem(nel, nuparam, nuvar, uparam, uvar, time, timestep, ssp, aldt, crklen, elcrkini, tdel, off, offly, signxx, signyy, signxy, signyz, signzx, ngl, ixel, ilay, ipt, npt, ixfem, iorth, dir_a, dir1_crk, dir2_crk, crkdir)
Definition fail_wind_xfem.F:37
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21