OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
prony25c.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!|| prony25c ../engine/source/materials/mat/mat025/prony25c.F
25!||--- called by ------------------------------------------------------
26!|| sigeps25c ../engine/source/materials/mat/mat025/sigeps25c.F
27!||--- calls -----------------------------------------------------
28!|| roto_sig ../engine/source/airbag/roto.F
29!||====================================================================
30 SUBROUTINE prony25c(NEL ,NPRONY ,BETA ,KV ,
31 2 GV ,TIMESTEP,RHO0 ,OFF ,DIR ,
32 3 EPSPXX ,EPSPYY ,EPSPXY ,EPSPYZ ,EPSPZX ,
33 4 SIGVXX ,SIGVYY ,SIGVXY ,SIGVYZ ,SIGVZX ,
34 5 SIGV ,SOUNDSP ,VARI ,IGTYP )
35C-----------------------------------------------
36C I m p l i c i t T y p e s
37C-----------------------------------------------
38#include "implicit_f.inc"
39C-----------------------------------------------
40C G l o b a l P a r a m e t e r s
41C-----------------------------------------------
42#include "mvsiz_p.inc"
43C---------+---------+---+---+--------------------------------------------
44C VAR | SIZE |TYP| RW| DEFINITION
45C---------+---------+---+---+--------------------------------------------
46C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
47C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
48C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
49C---------+---------+---+---+--------------------------------------------
50C---------+---------+---+---+--------------------------------------------
51C TIME | 1 | F | R | CURRENT TIME
52C TIMESTEP| 1 | F | R | CURRENT TIME STEP
53C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
54C RHO0 | NEL | F | R | INITIAL DENSITY
55C EPSPXX | NEL | F | R | STRAIN RATE XX
56C EPSPYY | NEL | F | R | STRAIN RATE YY
57C ... | | | |
58C ... | | | |
59C SIGVXX | NEL | F | W | VISCOUS STRESS XX
60C SIGVYY | NEL | F | W | VISCOUS STRESS YY
61C ... | | | |
62C VISCMAX | NEL | F | W | MAXIMUM DAMPING MODULUS(NEEDED FOR TIME STEP)
63C---------+---------+---+---+--------------------------------------------
64C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
65C OFF | NEL | F |R/W| DELETED ELEMENT FLAG (=1. ON, =0. OFF)
66C---------+---------+---+---+--------------------------------------------
67C-----------------------------------------------
68C I N P U T A r g u m e n t s
69C-----------------------------------------------
70C
71 INTEGER NEL ,NPRONY,IGTYP
72 my_real
73 . TIME,TIMESTEP,KV,GV(NPRONY),BETA(NPRONY),
74 . RHO0(*),EPSPXX(*),EPSPYY(*),
75 . epspxy(*),epspyz(*),epspzx(*),sigv(nel,5),dir(nel,2)
76C-----------------------------------------------
77C O U T P U T A r g u m e n t s
78C-----------------------------------------------
79 my_real
80 . sigvxx(*),sigvyy(*),
81 . sigvxy(*),sigvyz(*),sigvzx(*),
82 . soundsp(*)
83C-----------------------------------------------
84C I N P U T O U T P U T A r g u m e n t s
85C-----------------------------------------------
86 my_real
87 . vari(nel,*) , off(*)
88C-----------------------------------------------
89C VARIABLES FOR FUNCTION INTERPOLATION
90C-----------------------------------------------
91C Y = FINTER(IFUNC(J),X,NPF,TF,DYDX)
92C Y : y = f(x)
93C X : x
94C DYDX : f'(x) = dy/dx
95C IFUNC(J): FUNCTION INDEX
96C J : FIRST(J=1), SECOND(J=2) .. FUNCTION USED FOR THIS LAW
97C NPF,TF : FUNCTION PARAMETER
98C-----------------------------------------------
99C L o c a l V a r i a b l e s
100C-----------------------------------------------
101 integer
102 . i,j,ii,ll
103 my_real
104 . cc,g1,dav,epxx(mvsiz),epyy(mvsiz),epzz(mvsiz),svxx,
105 . svyy,svzz,p(mvsiz),
106 . aa(mvsiz,nprony),bb(mvsiz,nprony),h0(6),epspzz(mvsiz),
107 . h(mvsiz,6,nprony),s(mvsiz,6),a1(mvsiz),a2(mvsiz),fac,h01,h02,h03,h04,
108 . h05,h06
109C-----------------------------------------------
110C USER VARIABLES INITIALIZATION
111C-----------------------------------------------
112C------------------------------------------
113C estiamte stress
114C------------------------------------------
115C
116 g1 = zero
117 DO i=1,nprony
118 g1 = g1 + gv(i)
119 ENDDO
120C
121 a1(1:nel) = zero
122 a2(1:nel) = zero
123 epspzz(1:nel) = zero
124C
125 DO j=1,nprony
126 DO i= 1,nel
127 aa(i,j)= exp(-beta(j)*timestep)
128 bb(i,j)=timestep*gv(j)*exp(-half*beta(j)*timestep)
129C
130C for computing epszz_dot (sigzz=0)
131C
132 h03 = vari(i, (j - 1)*6 + 3)
133 a1(i) = a1(i) + aa(i,j)*h03
134 a2(i) = a2(i) + bb(i,j)
135 ENDDO
136 ENDDO
137C
138C compute epszz_dot sig33= 0
139C
140 DO i= 1,nel
141 fac = one/max(em20,two_third*a2(i) + kv)
142 epspzz(i) = -a1(i)+(third*a2(i)-kv)*(epspxx(i)+epspyy(i))
143 epspzz(i)= fac*epspzz(i)
144 ENDDO
145C
146 DO i= 1,nel
147C spherique part
148 dav = third*(epspxx(i) + epspyy(i) + epspzz(i))
149 p(i) = -three*kv*dav
150c deviatorique part
151 epxx(i) = epspxx(i) - dav
152 epyy(i) = epspyy(i) - dav
153 epzz(i) = epspzz(i) - dav
154 ENDDO
155C
156 DO j= 1,nprony
157 DO i=1,nel
158C
159 h01 = vari(i, (j - 1)*6 + 1)
160 h02 = vari(i, (j - 1)*6 + 2)
161 h03 = vari(i, (j - 1)*6 + 3)
162 h04 = vari(i, (j - 1)*6 + 4)
163 h05 = vari(i, (j - 1)*6 + 5)
164 h06 = vari(i, (j - 1)*6 + 6)
165C
166
167 h(i,1,j) = aa(i,j)*h01 + bb(i,j)*epxx(i)
168 h(i,2,j) = aa(i,j)*h02 + bb(i,j)*epyy(i)
169 h(i,3,j) = aa(i,j)*h03 + bb(i,j)*epzz(i)
170 h(i,4,j) = aa(i,j)*h04 + half*bb(i,j)*epspxy(i)
171 h(i,5,j) = aa(i,j)*h05 + half*bb(i,j)*epspyz(i)
172 h(i,6,j) = aa(i,j)*h06 + half*bb(i,j)*epspzx(i)
173C
174
175 vari(i, (j - 1)*6 + 1) = h(i,1,j)
176 vari(i, (j - 1)*6 + 2) = h(i,2,j)
177 vari(i, (j - 1)*6 + 3) = h(i,3,j)
178 vari(i, (j - 1)*6 + 4) = h(i,4,j)
179 vari(i, (j - 1)*6 + 5) = h(i,5,j)
180 vari(i, (j - 1)*6 + 6) = h(i,6,j)
181 ENDDO
182 ENDDO
183C
184C comppute stress
185C
186
187 DO ii = 1,6
188 s(1:mvsiz,ii) = zero
189 ENDDO
190
191 DO j= 1,nprony
192 DO i=1,nel
193 s(i,1) = s(i,1) + h(i,1,j)
194 s(i,2) = s(i,2) + h(i,2,j)
195cc S(I,3) = S(I,3) + H(I,3,J)
196 s(i,4) = s(i,4) + h(i,4,j)
197 s(i,5) = s(i,5) + h(i,5,j)
198 s(i,6) = s(i,6) + h(i,6,j)
199 ENDDO
200 ENDDO
201 DO i=1,nel
202 sigvxx(i) = s(i,1) - p(i)
203 sigvyy(i) = s(i,2) - p(i)
204 sigvxy(i) = s(i,4)
205 sigvyz(i) = s(i,5)
206 sigvzx(i) = s(i,6)
207 ENDDO
208C
209 DO i=1,nel
210 sigvxx(i) = sigvxx(i)*off(i)
211 sigvyy(i) = sigvyy(i)*off(i)
212 sigvxy(i) = sigvxy(i)*off(i)
213 sigvyz(i) = sigvyz(i)*off(i)
214 sigvzx(i) = sigvzx(i)*off(i)
215C
216 sigv(i,1) = sigvxx(i)
217 sigv(i,2) = sigvyy(i)
218 sigv(i,3) = sigvxy(i)
219 sigv(i,4) = sigvyz(i)
220 sigv(i,5) = sigvzx(i)
221C
222 soundsp(i) = sqrt(soundsp(i)**2 + g1/rho0(i))
223 ENDDO
224
225C transformation dans le rep re d'orthotropie
226 ll = 1
227 CALL roto_sig(ll,nel,sigv, dir,nel)
228C
229 DO i=1,nel
230 sigvxx(i) = sigv(i,1)
231 sigvyy(i) = sigv(i,2)
232 sigvxy(i) = sigv(i,3)
233 sigvyz(i) = sigv(i,4)
234 sigvzx(i) = sigv(i,5)
235 ENDDO
236C
237 RETURN
238 END
max
#define max(a, b)
Definition macros.h:21
prony25c
subroutine prony25c(nel, nprony, beta, kv, gv, timestep, rho0, off, dir, epspxx, epspyy, epspxy, epspyz, epspzx, sigvxx, sigvyy, sigvxy, sigvyz, sigvzx, sigv, soundsp, vari, igtyp)
Definition prony25c.F:35
roto_sig
subroutine roto_sig(jft, jlt, sig, dir, nel)
Definition roto.F:126