OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_setoff_c.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_setoff_c ../engine/source/materials/fail/fail_setoff_c.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!||--- uses -----------------------------------------------------
29!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
30!|| failwave_mod ../common_source/modules/failwave_mod.F
31!|| mat_elem_mod ../common_source/modules/mat_elem/mat_elem_mod.F90
32!|| stack_mod ../engine/share/modules/stack_mod.F
33!||====================================================================
34 SUBROUTINE fail_setoff_c(ELBUF_STR,MAT_ELEM ,GEO ,PID ,
35 . NGL ,NEL ,NLAY ,NPTTOT ,
36 . THK_LY ,THKLY ,OFF ,STACK ,
37 . ISUBSTACK,IGTYP ,FAILWAVE ,FWAVE_EL ,
38 . NLAY_MAX ,LAYNPT_MAX,NUMGEO ,NUMSTACK ,
39 . IGEO ,PRINT_FAIL)
40C-----------------------------------------------
41C M o d u l e s
42C-----------------------------------------------
43 USE mat_elem_mod
44 USE stack_mod
45 USE failwave_mod
46 USE stack_mod
47 USE elbufdef_mod
48C-----------------------------------------------
49C I m p l i c i t T y p e s
50C-----------------------------------------------
51#include "implicit_f.inc"
52C-----------------------------------------------
53C C o m m o n B l o c k s
54C-----------------------------------------------
55#include "param_c.inc"
56#include "com08_c.inc"
57#include "units_c.inc"
58#include "comlock.inc"
59C-----------------------------------------------
60C D u m m y A r g u m e n t s
61C-----------------------------------------------
62 TYPE(elbuf_struct_), INTENT(INOUT), TARGET :: ELBUF_STR
63 my_real, DIMENSION(NPROPG,NUMGEO), INTENT(IN) :: GEO
64 INTEGER, DIMENSION(NPROPGI,NUMGEO),INTENT(IN) :: IGEO
65 INTEGER,INTENT(IN) :: PID,NEL,NLAY,NPTTOT,IGTYP,ISUBSTACK,
66 . NLAY_MAX,LAYNPT_MAX,NUMGEO,NUMSTACK
67 INTEGER, DIMENSION(NEL), INTENT(IN) :: NGL
68 my_real, DIMENSION(NEL,NLAY_MAX*LAYNPT_MAX), INTENT(IN) :: thk_ly
69 my_real, DIMENSION(NPTTOT*NEL), INTENT(IN) :: thkly
70 my_real, DIMENSION(NEL), INTENT(INOUT) :: off
71 TYPE (STACK_PLY), INTENT(IN) :: STACK
72 TYPE (FAILWAVE_STR_), INTENT(IN), TARGET :: FAILWAVE
73 INTEGER, DIMENSION(NEL), INTENT(INOUT) :: FWAVE_EL
74 TYPE (MAT_ELEM_) ,INTENT(INOUT) :: MAT_ELEM
75 LOGICAL, DIMENSION(NEL), INTENT(INOUT) :: PRINT_FAIL
76C-----------------------------------------------
77C L o c a l V a r i a b l e s
78C-----------------------------------------------
79 INTEGER I,II,IEL,IPOS,IL,IFL,IPT,NPTT,
80 . nindxly,nfail,ipweight,ipthkly,imat,id_ply
81 my_real
82 . p_thickg,fail_exp,dfail
83 INTEGER, DIMENSION(NEL) :: INDXLY
84 INTEGER, DIMENSION(:), POINTER :: FOFF,LAY_OFF
85 my_real, DIMENSION(NLAY) :: WEIGHT,P_THKLY
86 my_real, DIMENSION(NLAY,100) :: pthkf
87 my_real, DIMENSION(NEL) :: thfact,norm,npfail
88 TYPE(l_bufel_) ,POINTER :: LBUF
89 CHARACTER(LEN=NCHARTITLE), DIMENSION(NEL) :: FAIL_NAME
90c-----------------------------------------------------------------------
91c NPTT NUMBER OF INTEGRATION POINTS IN CURRENT LAYER
92c NPTTOT NUMBER OF INTEGRATION POINTS IN ALL LAYERS (TOTAL)
93c THK_LY Ratio of layer thickness / element thickness
94c THK Total element thickness
95C=======================================================================
96c
97 !=================================================================
98 ! RECOVER PARAMETERS AND INITIALIZATION
99 !=================================================================
100 p_thickg = geo(42,pid) ! General PTHICKFAIL
101 fail_exp = geo(43,pid) ! Failure exponent for TYPE 51
102 IF (igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
103 ipthkly = 1 + 4*nlay ! Address of PTHKLY in STACK%GEO table
104 ipweight = ipthkly + nlay ! Address of WEIGHT in STACK%GEO table
105 ELSE
106 ipthkly = 700 ! Address of PTHKLY in GEO table
107 ipweight = 900 ! Address of WEIGHT in GEO table
108 ENDIF
109c
110 DO il=1,nlay
111 nfail = elbuf_str%BUFLY(il)%NFAIL
112 imat = elbuf_str%BUFLY(il)%IMAT
113 DO ifl = 1,nfail
114 pthkf(il,ifl) = mat_elem%MAT_PARAM(imat)%FAIL(ifl)%PTHK
115 END DO
116 END DO
117 !=================================================================
118 ! 1 LAYER PROPERTIES - TYPE1/TYPE9
119 !=================================================================
120 IF (nlay == 1) THEN
121c
122 ! Only 1 layer with several integration points
123 il = 1
124c
125 ! Material internal identifier number
126 imat = elbuf_str%BUFLY(il)%IMAT
127c
128 ! Check PTHICKFAIL coming from failure criteria
129 nfail = elbuf_str%BUFLY(il)%NFAIL
130 DO ifl = 1,nfail
131 ! -> Percentage of broken thickness
132 IF (pthkf(il,ifl) > zero) THEN
133 pthkf(il,ifl) = min(pthkf(il,ifl),abs(p_thickg))
134 pthkf(il,ifl) = max(min(pthkf(il,ifl),one-em06),em06)
135 ! -> Ratio of broken integration points
136 ELSEIF (pthkf(il,ifl) < zero) THEN
137 pthkf(il,ifl) = max(pthkf(il,ifl),-abs(p_thickg))
138 pthkf(il,ifl) = min(max(pthkf(il,ifl),-one+em6),-em06)
139 ! -> If not defined in the failure criterion, the value of the property is used by default
140 ELSE
141 pthkf(il,ifl) = p_thickg
142 ENDIF
143 ENDDO ! |-> IFL
144c
145 ! Check element failure
146 DO ifl = 1,nfail
147 thfact(1:nel) = zero
148 npfail(1:nel) = zero
149 nptt = elbuf_str%BUFLY(il)%NPTT
150 ! Computation of broken fraction of thickness /
151 ! ratio of intg. points
152 DO ipt=1,nptt
153 foff => elbuf_str%BUFLY(il)%FAIL(1,1,ipt)%FLOC(ifl)%OFF
154 DO iel=1,nel
155 IF (off(iel) == one) THEN
156 IF (foff(iel) < one) THEN
157 ipos = (ipt-1)*nel + iel
158 thfact(iel) = thfact(iel) + thkly(ipos) ! -> Percentage of broken thickness (accounting for integration weight)
159 npfail(iel) = npfail(iel) + one/nptt ! -> Ratio of broken integration points
160 ENDIF
161 ENDIF
162 ENDDO ! |-> IEL
163 ENDDO ! |-> IPT
164 ! Comparison with critical value PTHICKFAIL
165 DO iel=1,nel
166 IF (off(iel) == one) THEN
167 IF (((thfact(iel) >= pthkf(il,ifl)).AND.(pthkf(il,ifl) > zero)).OR.
168 . ((npfail(iel) >= abs(pthkf(il,ifl))).AND.(pthkf(il,ifl) < zero))) THEN
169 off(iel) = four_over_5
170 print_fail(iel) = .false.
171 fail_name(iel) = mat_elem%MAT_PARAM(imat)%FAIL(ifl)%KEYWORD
172#include "lockon.inc"
173 WRITE(iout, 1000) trim(fail_name(iel)),ngl(iel)
174 WRITE(istdo,1100) trim(fail_name(iel)),ngl(iel),tt
175#include "lockoff.inc"
176 IF (failwave%WAVE_MOD > 0) fwave_el(iel) = -1
177 ENDIF
178 ENDIF
179 ENDDO ! |-> IEL
180 ENDDO ! |-> ifl
181c
182 !=================================================================
183 ! MULTI LAYER PROPERTIES / 1 INTG POINT - TYPE10/11/16/17/51/52
184 !=================================================================
185 ELSEIF (nlay == npttot) THEN
186c
187 ! Only one integration points in each layer
188 ipt = 1
189c
190 ! Check layer failure
191 DO il=1,nlay
192 nindxly = 0
193 nfail = elbuf_str%BUFLY(il)%NFAIL
194 lay_off => elbuf_str%BUFLY(il)%OFF
195 imat = elbuf_str%BUFLY(il)%IMAT
196 DO ifl = 1,nfail
197 foff => elbuf_str%BUFLY(il)%FAIL(1,1,ipt)%FLOC(ifl)%OFF
198 DO iel=1,nel
199 IF (off(iel) == one .AND. lay_off(iel) == 1) THEN
200 IF (foff(iel) < 1) THEN
201 nindxly = nindxly + 1
202 indxly(nindxly) = iel
203 lay_off(iel) = 0
204 fail_name(iel) = mat_elem%MAT_PARAM(imat)%FAIL(ifl)%KEYWORD
205 ENDIF
206 ENDIF
207 ENDDO ! |-> IEL
208 ENDDO ! |-> IFL
209 ! Printing out layer/ply failure message
210 IF (nindxly > 0) THEN
211 ! -> Print out ply failure
212 IF (igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
213 IF (igtyp == 17 .OR. igtyp == 51) THEN
214 id_ply = igeo(1,stack%IGEO(2+il,isubstack))
215 ELSE
216 id_ply = ply_info(1,stack%IGEO(2+il,isubstack)-numstack)
217 ENDIF
218 DO i = 1,nindxly
219#include "lockon.inc"
220 WRITE(iout, 3000) trim(fail_name(indxly(i))),id_ply,ngl(indxly(i))
221 WRITE(istdo,3100) trim(fail_name(indxly(i))),id_ply,ngl(indxly(i)),tt
222#include "lockoff.inc"
223 ENDDO ! |-> I
224 ! -> Print out layer failure
225 ELSE
226 DO i = 1,nindxly
227#include "lockon.inc"
228 WRITE(iout, 2000) trim(fail_name(indxly(i))),il,ngl(indxly(i))
229 WRITE(istdo,2100) trim(fail_name(indxly(i))),il,ngl(indxly(i)),tt
230#include "lockoff.inc"
231 ENDDO ! |-> I
232 ENDIF
233 ENDIF
234 ENDDO ! |-> IL
235c
236 ! Check element failure
237 thfact(1:nel) = zero
238 norm(1:nel) = zero
239 npfail(1:nel) = zero
240 ! Computation of broken fraction of thickness / ratio of layers
241 DO il=1,nlay
242 IF (igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
243 weight(il) = stack%GEO(ipweight+ il,isubstack)
244 ELSE
245 weight(il) = geo(ipweight + il,pid)
246 ENDIF
247 lay_off => elbuf_str%BUFLY(il)%OFF
248 ii = (il-1)*nel
249 DO iel=1,nel
250 IF (off(iel) == one) THEN
251 ipos = ii + iel
252 dfail = thkly(ipos)*weight(il)
253 norm(iel) = norm(iel) + dfail
254 IF (lay_off(iel) == 0) THEN
255 thfact(iel) = thfact(iel) + thkly(ipos)*weight(il) ! -> Percentage of broken thickness
256 npfail(iel) = npfail(iel) + one/nlay ! -> Ratio of broken layers
257 ENDIF
258 ENDIF
259 ENDDO ! |-> IEL
260 ENDDO ! |-> IL
261 ! Comparison with critical value PTHICKFAIL
262 DO iel=1,nel
263 IF (off(iel) == one) THEN
264 IF (((thfact(iel) >= p_thickg*norm(iel)).AND.(p_thickg > zero)).OR.
265 . ((npfail(iel) >= abs(p_thickg)).AND.(p_thickg < zero))) THEN
266 off(iel) = four_over_5
267 IF (failwave%WAVE_MOD > 0) fwave_el(iel) = -1
268 ENDIF
269 ENDIF
270 ENDDO ! |-> IEL
271c
272 !=======================================================================
273 ! MULTI LAYER PROPERTIES / SEVERAL INTG. POINTS - TYPE51/TYPE52
274 !=======================================================================
275 ELSE
276c
277 ! Check PTHICKFAIL coming from failure criteria
278 DO il = 1,nlay
279 nfail = elbuf_str%BUFLY(il)%NFAIL
280 p_thkly(il) = stack%GEO(ipthkly + il,isubstack)
281 DO ifl = 1,nfail
282 ! -> Percentage of broken thickness
283 IF (pthkf(il,ifl) > zero) THEN
284 pthkf(il,ifl) = min(pthkf(il,ifl),abs(p_thkly(il)))
285 pthkf(il,ifl) = max(min(pthkf(il,ifl),one-em06),em06)
286 ! -> Ratio of broken integration points
287 ELSEIF (pthkf(il,ifl) < zero) THEN
288 pthkf(il,ifl) = max(pthkf(il,ifl),-abs(p_thkly(il)))
289 pthkf(il,ifl) = min(max(pthkf(il,ifl),-one+em6),-em06)
290 ! -> If not defined in the failure criterion, the value of the property is used by default
291 ELSE
292 pthkf(il,ifl) = p_thkly(il)
293 ENDIF
294 ENDDO ! |-> IFL
295 ENDDO ! |-> IL
296c
297 ! Check layer failure
298 DO il=1,nlay
299 nptt = elbuf_str%BUFLY(il)%NPTT
300 nindxly = 0
301 nfail = elbuf_str%BUFLY(il)%NFAIL
302 lay_off => elbuf_str%BUFLY(il)%OFF
303 imat = elbuf_str%BUFLY(il)%IMAT
304 ii = (il-1)*nel
305 ! Computation of broken fraction of thickness / ratio of integration points
306 DO ifl = 1,nfail
307 thfact(1:nel) = zero
308 npfail(1:nel) = zero
309 DO ipt = 1,nptt
310 foff => elbuf_str%BUFLY(il)%FAIL(1,1,ipt)%FLOC(ifl)%OFF
311 DO iel=1,nel
312 IF (off(iel) == one) THEN
313 IF (lay_off(iel) == 1) THEN
314 IF (foff(iel) < one) THEN
315 ipos = ii + iel
316 thfact(iel) = thfact(iel) + thkly(ipos)/thk_ly(iel,il) ! -> Percentage of broken thickness
317 npfail(iel) = npfail(iel) + one/nptt ! -> Ratio of broken integration points
318 ENDIF
319 ENDIF
320 ENDIF
321 ENDDO ! |-> IEL
322 ENDDO ! |-> IPT
323 ! Comparison with critical value PTHICKFAIL
324 DO iel=1,nel
325 IF (off(iel) == one) THEN
326 IF (((thfact(iel) >= pthkf(il,ifl)).AND.(pthkf(il,ifl) > zero)).OR.
327 . ((npfail(iel) >= abs(pthkf(il,ifl))).AND.(pthkf(il,ifl) < zero))) THEN
328 nindxly = nindxly + 1
329 indxly(nindxly) = iel
330 lay_off(iel) = 0
331 fail_name(iel) = mat_elem%MAT_PARAM(imat)%FAIL(ifl)%KEYWORD
332 DO ipt=1,nptt
333 foff => elbuf_str%BUFLY(il)%FAIL(1,1,ipt)%FLOC(ifl)%OFF
334 foff(iel) = 0
335 ENDDO
336 ENDIF
337 ENDIF
338 ENDDO ! |-> IEL
339 ENDDO ! |-> ifl
340 ! printing out ply failure message
341 IF (nindxly > 0) THEN
342 IF (igtyp == 51) THEN
343 id_ply = igeo(1,stack%IGEO(2+il,isubstack))
344 ELSE
345 id_ply = ply_info(1,stack%IGEO(2+il,isubstack)-numstack)
346 ENDIF
347 DO i = 1,nindxly
348#include "lockon.inc"
349 WRITE(iout, 3000) trim(fail_name(indxly(i))),id_ply,ngl(indxly(i))
350 WRITE(istdo,3100) trim(fail_name(indxly(i))),id_ply,ngl(indxly(i)),tt
351#include "lockoff.inc"
352 ENDDO ! |-> I
353 ENDIF
354 ENDDO ! |-> IL
355c
356 ! Check element failure
357 thfact(1:nel) = zero
358 npfail(1:nel) = zero
359 norm(1:nel) = zero
360 ! Computation of broken fraction of thickness / ratio of layers
361 DO il=1,nlay
362 weight(il) = stack%GEO(ipweight + il,isubstack)
363 lay_off => elbuf_str%BUFLY(il)%OFF
364 DO iel=1,nel
365 IF (off(iel) == one) THEN
366 dfail = (thk_ly(iel,il)*weight(il))**fail_exp
367 norm(iel) = norm(iel) + dfail
368 IF (lay_off(iel) == 0) THEN
369 thfact(iel) = thfact(iel) + dfail ! -> Percentage of broken thickness (among layers)
370 npfail(iel) = npfail(iel) + one/nlay ! -> Ratio of broken layers
371 ENDIF
372 ENDIF
373 ENDDO ! |-> IEL
374 ENDDO ! |-> IL
375 ! Comparison with critical value PTHICKFAIL
376 DO iel=1,nel
377 IF (off(iel) == one) THEN
378 thfact(iel) = thfact(iel)**(one/fail_exp)
379 norm(iel) = norm(iel)**(one/fail_exp)
380 IF (((thfact(iel) >= p_thickg*norm(iel)).AND.(p_thickg > zero)).OR.
381 . ((npfail(iel) >= abs(p_thickg)).AND.(p_thickg < zero))) THEN
382 off(iel) = four_over_5
383 IF (failwave%WAVE_MOD > 0) fwave_el(iel) = -1
384 ENDIF
385 ENDIF
386 ENDDO ! |-> IEL
387c
388 ENDIF ! IGTYP PROPERTY TYPE
389 !=======================================================================
390c
391 !=======================================================================
392 ! PRINTING OUT FORMATS FOR LAYERS/PLYS FAILURE
393 !=======================================================================
394 1000 FORMAT(1x,'-- RUPTURE (',a,') OF SHELL ELEMENT NUMBER ',i10)
395 1100 FORMAT(1x,'-- RUPTURE (',a,') OF SHELL ELEMENT :',i10,' AT TIME :',g11.4)
396 2000 FORMAT(1x,'-- FAILURE (',a,') OF LAYER',i3, ' ,SHELL ELEMENT NUMBER ',i10)
397 2100 FORMAT(1x,'-- FAILURE (',a,') OF LAYER',i3, ' ,SHELL ELEMENT NUMBER ',i10,
398 . 1x,'AT TIME :',g11.4)
399 3000 FORMAT(1x,'-- FAILURE (',a,') OF PLY ID ',i10, ' ,SHELL ELEMENT NUMBER ',i10)
400 3100 FORMAT(1x,'-- FAILURE (',a,') OF PLY ID ',i10, ' ,SHELL ELEMENT NUMBER ',i10,
401 . 1x,'AT TIME :',g11.4)
402 END
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
subroutine fail_setoff_c(elbuf_str, mat_elem, geo, pid, ngl, nel, nlay, npttot, thk_ly, thkly, off, stack, isubstack, igtyp, failwave, fwave_el, nlay_max, laynpt_max, numgeo, numstack, igeo, print_fail)
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
integer, dimension(:,:), allocatable ply_info
Definition stack_mod.F:133