OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
main_beam18.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!|| main_beam18 ../engine/source/elements/beam/main_beam18.F
25!||--- called by ------------------------------------------------------
26!|| pforc3 ../engine/source/elements/beam/pforc3.F
27!||--- calls -----------------------------------------------------
28!|| fail_beam18 ../engine/source/elements/beam/fail_beam18.F
29!|| m2lawpi ../engine/source/materials/mat/mat002/m2lawpi.F
30!|| mulaw_ib ../engine/source/elements/beam/mulaw_ib.F
31!||--- uses -----------------------------------------------------
32!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
33!|| mat_elem_mod ../common_source/modules/mat_elem/mat_elem_mod.F90
34!||====================================================================
35 SUBROUTINE main_beam18(ELBUF_STR,
36 1 NEL ,NPT ,MTN ,IMAT ,
37 2 PID ,NGL ,PM ,IPM ,
38 3 GEO ,OFF ,FOR ,MOM ,
39 4 EINT ,AL ,EPSD ,BUFMAT ,NPF ,
40 5 TF ,EXX ,EXY ,EXZ ,KXX ,
41 6 KYY ,KZZ ,F1 ,F2 ,F3 ,
42 7 M1 ,M2 ,M3 ,JTHE ,TEMPEL ,
43 8 IFAIL ,SBUFMAT ,SNPC ,STF ,NUMMAT ,
44 9 NUMGEO ,IOUT ,ISTDO ,NPROPMI ,NPROPM ,
45 A NPROPG ,TIME ,DTIME ,IDEL7NOK,ISIGI ,
46 B IMCONV ,ISMSTR ,MAT_PARAM ,
47 C NTABLE ,TABLE )
48
49C-----------------------------------------------
50C M o d u l e s
51C-----------------------------------------------
52 USE mat_elem_mod
53 USE elbufdef_mod
54C-----------------------------------------------
55C I m p l i c i t T y p e s
56C-----------------------------------------------
57#include "implicit_f.inc"
58#include "comlock.inc"
59C-----------------------------------------------
60C D u m m y A r g u m e n t s
61C-----------------------------------------------
62 INTEGER ,INTENT(IN) :: IMAT
63 INTEGER ,INTENT(IN) :: NEL,MTN,NPT,JTHE,IFAIL
64 INTEGER ,INTENT(IN) :: SBUFMAT
65 INTEGER ,INTENT(IN) :: SNPC
66 INTEGER ,INTENT(IN) :: STF
67 INTEGER ,INTENT(IN) :: NUMMAT
68 INTEGER ,INTENT(IN) :: NUMGEO
69 INTEGER ,INTENT(IN) :: NPROPMI
70 INTEGER ,INTENT(IN) :: NPROPM
71 INTEGER ,INTENT(IN) :: NPROPG
72 INTEGER ,INTENT(IN) :: IOUT
73 INTEGER ,INTENT(IN) :: ISTDO
74 INTEGER ,INTENT(IN) :: ISIGI
75 INTEGER ,INTENT(IN) :: IMCONV
76 INTEGER ,INTENT(IN) :: ISMSTR
77 INTEGER,INTENT(IN) :: NTABLE
78 TYPE(ttable), DIMENSION(NTABLE), INTENT(INOUT) :: TABLE ! TABLE DATA
79 INTEGER ,DIMENSION(NEL) ,INTENT(IN) :: PID
80 INTEGER ,DIMENSION(NEL) ,INTENT(IN) :: NGL
81 INTEGER ,DIMENSION(SNPC) ,INTENT(IN) :: NPF
82 INTEGER IPM(NPROPMI,NUMMAT)
83 INTEGER ,INTENT(INOUT) :: IDEL7NOK
84 my_real ,INTENT(IN) :: time
85 my_real ,INTENT(IN) :: dtime
86 my_real ,INTENT(IN) :: pm(npropm,nummat)
87 my_real ,INTENT(IN) :: geo(npropg,numgeo)
88 my_real ,DIMENSION(NEL) ,INTENT(IN) :: al
89 my_real ,DIMENSION(NEL) ,INTENT(OUT) :: epsd
90 my_real ,DIMENSION(SBUFMAT) ,INTENT(IN) :: bufmat
91 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: tempel
92 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: off
93 my_real ,DIMENSION(NEL,2),INTENT(INOUT) :: eint
94 my_real ,DIMENSION(NEL,3),INTENT(INOUT) :: for,mom
95 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: exx,exy,exz
96 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: kxx,kyy,kzz
97 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: f1,f2,f3
98 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: m1,m2,m3
99 my_real ,DIMENSION(STF) ,INTENT(IN) :: tf
100C
101 TYPE (ELBUF_STRUCT_) ,INTENT(INOUT) :: ELBUF_STR
102 TYPE (MATPARAM_STRUCT_) ,INTENT(IN) :: MAT_PARAM
103C-----------------------------------------------
104C L o c a l V a r i a b l e s
105C-----------------------------------------------
106 INTEGER :: I,J,IPT,IPID,IPLA,NPAR,IADBUF,NFUNC,ISRATE,IPY,IPZ,IPA,II(3)
107 my_real :: dtinv,asrate,epsdi,dmpm,dmpf,spe,spg,rho,fact,area,
108 . dmm,dmpfl,ixx,iyy,izz,g,e,damage_loc,dfxx,dfxy,dfxz,
109 . signxx,signxy,signxz,sigoxx,sigoxy,sigoxz
110 my_real ,DIMENSION(NEL ):: degmb,degsh,degfx,ypt,zpt,apt,off_old
111 my_real ,DIMENSION(:,:) ,ALLOCATABLE :: dpla,sigy,vm,vm0
112C=======================================================================
113 ipy = 200
114 ipz = 300
115 ipa = 400
116 ipid = pid(1)
117 area = geo(1 ,ipid)
118 DO i=1,3
119 ii(i) = nel*(i-1)
120 ENDDO
121C-------------------
122C STRAIN RATE
123C-------------------
124 israte = ipm(3,imat)
125 asrate = min(one, pm(9,imat)*dtime)
126c calculate total strain rate on neutral fiber position
127 DO i = 1,nel
128 epsdi = sqrt(exx(i)**2 + half*(exy(i)**2 + exz(i)**2))
129 IF (israte > 0) THEN
130 epsd(i)= asrate*epsdi + (one - asrate)*epsd(i)
131 ELSE
132 epsd(i)= epsdi
133 ENDIF
134 ENDDO
135c
136 ! Save old OFF value
137 off_old(1:nel) = off(1:nel)
138c-------------------
139c STRAIN
140c-------------------
141 DO i = 1,nel
142 exx(i) = exx(i)*dtime
143 exy(i) = exy(i)*dtime
144 exz(i) = exz(i)*dtime
145 kxx(i) = kxx(i)*dtime
146 kyy(i) = kyy(i)*dtime
147 kzz(i) = kzz(i)*dtime
148 ENDDO
149C
150 DO i = 1,nel
151 degmb(i) = for(i,1)*exx(i)
152 degsh(i) = for(i,2)*exy(i) + for(i,3)*exz(i)
153 degfx(i) = mom(i,1)*kxx(i) + mom(i,2)*kyy(i) + mom(i,3)*kzz(i)
154 ENDDO
155C
156 IF (isigi == 0 .OR. (isigi /= 0 .AND. time /= zero)) THEN
157 DO i = 1,nel
158 for(i,1) = zero
159 for(i,2) = zero
160 for(i,3) = zero
161 mom(i,1) = zero
162 mom(i,2) = zero
163 mom(i,3) = zero
164 ENDDO
165 ENDIF
166c
167 ! Plasticity specific tables
168 ipla = elbuf_str%BUFLY(1)%L_PLA
169 ALLOCATE (dpla(nel*ipla,npt*ipla))
170 ALLOCATE (sigy(nel*ipla,npt*ipla))
171 ALLOCATE (vm0(nel*ipla,npt*ipla))
172 ALLOCATE (vm(nel*ipla,npt*ipla))
173 IF (ipla > 0) THEN
174 DO j=1,npt
175 DO i = 1,nel
176 dpla(i,j) = elbuf_str%BUFLY(1)%LBUF(1,1,j)%PLA(i)
177 sigy(i,j) = 1.e30
178 sigoxx = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(1)+i)
179 sigoxy = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(2)+i)
180 sigoxz = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(3)+i)
181 vm0(i,j) = sqrt(sigoxx**2 + three*(sigoxy**2 + sigoxz**2))
182 ENDDO
183 ENDDO
184 ENDIF
185c---------------------------
186c material models
187c---------------------------
188 IF (mtn == 2) THEN
189 CALL m2lawpi(elbuf_str,
190 1 1 ,nel ,npt ,pm ,geo ,
191 2 eint ,off ,imat ,
192 3 pid ,epsd ,exx ,exy ,exz ,
193 4 kxx ,kyy ,kzz ,al ,nel ,
194 5 ipm ,asrate ,dtime ,nummat ,sigy)
195c
196 ELSE
197c
198 CALL mulaw_ib(elbuf_str,
199 1 nel ,npt ,mtn ,imat ,
200 2 pid ,ngl ,ipm ,
201 3 geo ,off ,
202 4 epsd ,bufmat ,npf ,tf ,
203 5 exx ,exy ,exz ,kxx ,
204 6 kyy ,kzz ,jthe ,tempel ,sigy )
205 ENDIF
206!
207c---------------------------
208c Computation of Plastic Work
209c---------------------------
210 IF (ipla > 0) THEN
211 DO j=1,npt
212 DO i = 1,nel
213 dpla(i,j) = elbuf_str%BUFLY(1)%LBUF(1,1,j)%PLA(i) - dpla(i,j)
214 signxx = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(1)+i)
215 signxy = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(2)+i)
216 signxz = elbuf_str%BUFLY(1)%LBUF(1,1,j)%SIG(ii(3)+i)
217 vm(i,j) = sqrt(signxx**2 + three*(signxy**2 + signxz**2))
218 IF (elbuf_str%GBUF%G_WPLA > 0) elbuf_str%GBUF%WPLA(i) = elbuf_str%GBUF%WPLA(i) +
219 . half*(vm(i,j)+vm0(i,j))*dpla(i,j)*al(i)*area/npt
220 ENDDO
221 ENDDO
222 ENDIF
223c---------------------------
224c failure models
225c---------------------------
226 IF (ifail > 0) THEN
227c
228 CALL fail_beam18(
229 . elbuf_str,mat_param%FAIL(1),nummat ,numgeo ,
230 . npropm ,npropg ,snpc ,stf ,
231 . nel ,npt ,imat ,ipid ,jthe ,
232 . tempel ,ngl ,pm ,geo ,
233 . off ,epsd ,npf ,tf ,
234 . dpla ,eint ,time ,iout ,istdo ,
235 . al ,ismstr ,exx ,exy ,exz ,
236 . kxx ,kyy ,kzz ,dtime ,
237 . ntable ,table ,sigy )
238 END IF
239C------------------------------------
240c resultant force and moment
241c------------------------------------
242 IF (elbuf_str%BUFLY(1)%L_DMGSCL > 0 ) THEN
243 DO ipt = 1,npt
244 DO i=1,nel
245 ypt(i) = geo(ipy+ipt,pid(i))
246 zpt(i) = geo(ipz+ipt,pid(i))
247 apt(i) = geo(ipa+ipt,pid(i))
248 damage_loc = elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%DMGSCL(i)
249 dfxx = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(1)+i)*damage_loc
250 dfxy = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(2)+i)*damage_loc
251 dfxz = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(3)+i)*damage_loc
252 for(i,1) = for(i,1) + dfxx
253 for(i,2) = for(i,2) + dfxy
254 for(i,3) = for(i,3) + dfxz
255 mom(i,1) = mom(i,1) + dfxy*zpt(i) - dfxz*ypt(i)
256 mom(i,2) = mom(i,2) + dfxx*zpt(i)
257 mom(i,3) = mom(i,3) - dfxx*ypt(i)
258 ENDDO
259 ENDDO
260 ELSE
261C-----------------------
262C FORCES ET MOMENTS
263C-----------------------
264 DO ipt = 1,npt
265 DO i=1,nel
266 ypt(i) = geo(ipy+ipt,pid(i))
267 zpt(i) = geo(ipz+ipt,pid(i))
268 apt(i) = geo(ipa+ipt,pid(i))
269 dfxx = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(1)+i)
270 dfxy = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(2)+i)
271 dfxz = apt(i)*elbuf_str%BUFLY(1)%LBUF(1,1,ipt)%SIG(ii(3)+i)
272 for(i,1) = for(i,1) + dfxx
273 for(i,2) = for(i,2) + dfxy
274 for(i,3) = for(i,3) + dfxz
275 mom(i,1) = mom(i,1) + dfxy*zpt(i) - dfxz*ypt(i)
276 mom(i,2) = mom(i,2) + dfxx*zpt(i)
277 mom(i,3) = mom(i,3) - dfxx*ypt(i)
278 ENDDO
279 ENDDO
280 ENDIF
281
282
283C-------------------------------------
284 DO i = 1,nel
285 for(i,1) = for(i,1)*off(i)
286 for(i,2) = for(i,2)*off(i)
287 for(i,3) = for(i,3)*off(i)
288 mom(i,1) = mom(i,1)*off(i)
289 mom(i,2) = mom(i,2)*off(i)
290 mom(i,3) = mom(i,3)*off(i)
291 f1(i) = for(i,1)
292 f2(i) = for(i,2)
293 f3(i) = for(i,3)
294 m1(i) = mom(i,1)
295 m2(i) = mom(i,2)
296 m3(i) = mom(i,3)
297 ENDDO
298c---------------------------
299c element damping is removed outside with a common call
300C---------------------------
301C Internal energy
302C---------------------------
303 DO i = 1,nel
304 degmb(i) = degmb(i) + for(i,1)*exx(i)
305 degsh(i) = degsh(i) + for(i,2)*exy(i) + for(i,3)*exz(i)
306 degfx(i) = degfx(i)
307 . + mom(i,1)*kxx(i) + mom(i,2)*kyy(i) + mom(i,3)*kzz(i)
308 fact = half*off(i)*al(i)
309 eint(i,1) = eint(i,1) + fact*(degmb(i)+degsh(i))
310 eint(i,2) = eint(i,2) + fact* degfx(i)
311 ENDDO
312c----------------------------
313c Check element erosion (IMCONV = 1 => convergence of implicit)
314c----------------------------
315 DO i = 1,nel
316 IF (off(i) == four_over_5 .AND. imconv == 1) THEN
317#include "lockon.inc"
318 WRITE(iout, 1000) ngl(i)
319 WRITE(istdo,1100) ngl(i),time
320#include "lockoff.inc"
321 ENDIF
322 ENDDO
323C
324 DO i = 1,nel
325 IF (off(i) < em01) off(i) = zero
326 IF (off(i) < one ) off(i) = off(i)*four_over_5
327 ENDDO
328c
329c--------------------------------------------------------
330c SHOOTING NODES ALGORITHM ACTIVATION
331c--------------------------------------------------------
332 DO i = 1,nel
333 IF ((off_old(i) > zero) .AND. (off(i) == zero)) THEN
334 idel7nok = 1
335 ENDIF
336 ENDDO
337c
338 DEALLOCATE (dpla)
339C------------------------------------------
340 1000 FORMAT(1x,'-- RUPTURE OF BEAM ELEMENT NUMBER ',i10)
341 1100 FORMAT(1x,'-- RUPTURE OF BEAM ELEMENT :',i10,' AT TIME :',g11.4)
342C------------------------------------------
343 RETURN
344 END SUBROUTINE main_beam18
my_real
#define my_real
Definition cppsort.cpp:32
fail_beam18
subroutine fail_beam18(elbuf_str, fail, nummat, numgeo, npropm, npropg, snpc, stf, nel, npt, imat, iprop, jthe, tempel, ngl, pm, geo, off, epsd, npf, tf, dpla, eint, time, iout, istdo, al, ismstr, exx, exy, exz, kxx, kyy, kzz, dtime, ntable, table, sigy)
Definition fail_beam18.F:54
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
m2lawpi
subroutine m2lawpi(elbuf_str, jft, jlt, npt, pm, geo, eint, off, imat, pid, epsp, exx, exy, exz, kxx, kyy, kzz, al, nel, ipm, asrate, timestep, nummat, sigy)
Definition m2lawpi.F:36
min
#define min(a, b)
Definition macros.h:20
main_beam18
subroutine main_beam18(elbuf_str, nel, npt, mtn, imat, pid, ngl, pm, ipm, geo, off, for, mom, eint, al, epsd, bufmat, npf, tf, exx, exy, exz, kxx, kyy, kzz, f1, f2, f3, m1, m2, m3, jthe, tempel, ifail, sbufmat, snpc, stf, nummat, numgeo, iout, istdo, npropmi, npropm, npropg, time, dtime, idel7nok, isigi, imconv, ismstr, mat_param, ntable, table)
Definition main_beam18.F:48
mulaw_ib
subroutine mulaw_ib(elbuf_str, nel, npt, mtn, imat, pid, ngl, ipm, geo, off, epsp, bufmat, npf, tf, exx, exy, exz, kxx, kyy, kzz, jthe, tempel, sigy)
Definition mulaw_ib.F:42
for
for(i8=*sizetab-1;i8 >=0;i8--)
Definition mumps_common.c:91