OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
m27plas.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!|| m27plas ../engine/source/materials/mat/mat027/m27plas.F
25!||--- called by ------------------------------------------------------
26!|| sigeps27c ../engine/source/materials/mat/mat027/sigeps27c.F
27!||====================================================================
28 SUBROUTINE m27plas(JFT ,JLT ,PM ,OFF ,SIG ,
29 2 PLA ,IMAT ,DT1 ,IPLA ,EZZ ,
30 3 EPSD_PG,ISRATE ,YLD ,ETSE ,DPLA ,
31 4 DEPSXX ,DEPSYY ,DEPSXY ,DEPSYZ ,DEPSZX,
32 5 NEL ,EPSD ,ASRATE )
33C-----------------------------------------------
34C I m p l i c i t T y p e s
35C-----------------------------------------------
36#include "implicit_f.inc"
37C-----------------------------------------------
38C G l o b a l P a r a m e t e r s
39C-----------------------------------------------
40#include "mvsiz_p.inc"
41C-----------------------------------------------
42C C o m m o n B l o c k s
43C-----------------------------------------------
44#include "param_c.inc"
45C-----------------------------------------------
46C D u m m y A r g u m e n t s
47C-----------------------------------------------
48 INTEGER JFT,JLT,IPLA,IMAT,ISRATE,NEL
49C REAL
50 my_real, INTENT(IN) :: ASRATE
51 my_real, INTENT(IN) :: DT1
52 my_real, INTENT(INOUT) :: sig(nel,5)
53 my_real
54 . pm(npropm,*),off(*),pla(*),ezz(*),dpla(*),
55 . depsxx(mvsiz),depsyy(mvsiz),depsxy(mvsiz),depsyz(mvsiz),
56 . depszx(mvsiz)
57 my_real
58 . yld(*),etse(*)
59 my_real, DIMENSION(NEL) ,INTENT(IN) :: epsd_pg
60 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: epsd
61C-----------------------------------------------
62C L o c a l V a r i a b l e s
63C-----------------------------------------------
64 INTEGER ICC(MVSIZ)
65 INTEGER I,J,N,NINDX,INDEX(MVSIZ),NMAX
66 my_real
67 . CA(MVSIZ),CB(MVSIZ),CN(MVSIZ),YMAX(MVSIZ),
68 . cc(mvsiz),epdr(mvsiz),svm(mvsiz),a(mvsiz),b(mvsiz),
69 . dpla_i(mvsiz),dpla_j(mvsiz),dr(mvsiz),h(mvsiz),
70 . nu1(mvsiz),nu2(mvsiz),p(mvsiz),q(mvsiz),nu(mvsiz),
71 . young(mvsiz),g(mvsiz),err,f,df,pla_i,p2,q2,r,s1,s2,s3,
72 . yld_i,nnu1,nnu2,nu3,nu4,nu5,nu6,sigz,pp,aa,bb,c,s11,
73 . s22,s12,s1s2,s122,umr,vm2,qq,small
74 my_real, DIMENSION(NEL) :: epsp
75 DATA nmax/3/
76C-----------------------------------------------
77 small = em7
78C
79#include "vectorize.inc"
80 DO i=jft,jlt
81 young(i)= pm(20,imat)
82 nu(i) = pm(21,imat)
83 g(i) = pm(22,imat)
84 ca(i) = pm(38,imat)
85 cb(i) = pm(39,imat)
86 cn(i) = pm(40,imat)
87 ymax(i) = pm(42,imat)
88 cc(i) = pm(43,imat)
89 epdr(i) = max(em20,pm(44,imat)*dt1)
90 icc(i) = nint(pm(49,imat))
91 etse(i) = one
92 ENDDO
93C-----------------------------------
94C VITESSE DE DEFORMATION & YIELD
95C-----------------------------------
96 DO i=jft,jlt
97 IF (israte == 0) THEN
98 epsp(i) = max(abs(depsxx(i)),abs(depsyy(i)),half*abs(depsxy(i)))
99 ELSE
100 epsd(i) = asrate*epsd_pg(i) + (one-asrate)*epsd(i)
101 epsp(i) = epsd(i) * dt1
102 END IF
103 epsp(i) = max(epsp(i),epdr(i))
104 q(i) = (one + cc(i) * log(epsp(i)/epdr(i)))
105 ENDDO
106C
107 DO i=jft,jlt
108 IF (icc(i) == 1) ymax(i) = ymax(i) * q(i)
109 ENDDO
110C
111 DO i=jft,jlt
112 dpla(i) = zero
113 yld(i) = (ca(i)+cb(i)*pla(i)**cn(i)) * q(i)
114 yld(i) = min(yld(i),ymax(i))
115 ca(i) = ca(i)*q(i)
116 cb(i) = cb(i)*q(i)
117 ENDDO
118C------------------------------------------
119C CONTRAINTES PLASTIQUEMENT ADMISSIBLES
120C------------------------------------------
121 IF (ipla == 0) THEN
122C-------------------
123C CONTRAINTE VM
124C-------------------
125 DO i=jft,jlt
126 svm(i) = sqrt(sig(i,1)*sig(i,1)
127 . +sig(i,2)*sig(i,2)
128 . -sig(i,1)*sig(i,2)
129 . +three*sig(i,3)*sig(i,3))
130 ENDDO
131C
132C projection radiale
133C
134#include "vectorize.inc"
135 DO i=jft,jlt
136 r = min(one,yld(i)/(svm(i)+ em15))
137 sig(i,1) = sig(i,1)*r
138 sig(i,2) = sig(i,2)*r
139 sig(i,3) = sig(i,3)*r
140 dpla(i) = off(i) * max(zero,(svm(i)-yld(i))/young(i))
141 h(i) = cn(i)*cb(i)*exp((cn(i)-one)*log(pla(i)+small))
142 IF (yld(i) >= ymax(i)) h(i) = zero
143 etse(i)= h(i)/(h(i)+young(i))
144 s1 = half*(sig(i,1)+sig(i,2))
145 ezz(i) = dpla(i) * s1 /yld(i)
146 pla(i) = pla(i) + dpla(i)
147 ENDDO
148C-------------------------
149 ELSEIF (ipla == 1) THEN
150C-------------------------
151C CRITERE DE VON MISES
152C-------------------------
153 DO i=jft,jlt
154 s1 = sig(i,1) + sig(i,2)
155 s2 = sig(i,1) - sig(i,2)
156 s3 = sig(i,3)
157 a(i) = fourth*s1*s1
158 b(i) = three_over_4*s2*s2+three*s3*s3
159 svm(i) = sqrt(a(i)+b(i))
160 ENDDO
161C-------------------------
162C GATHER PLASTIC FLOW
163C-------------------------
164 nindx=0
165 DO i=jft,jlt
166 IF (svm(i) > yld(i) .AND. off(i) == one) THEN
167 nindx = nindx + 1
168 index(nindx) = i
169 ENDIF
170 ENDDO
171 IF (nindx == 0) RETURN
172C---------------------------
173C DEP EN CONTRAINTE PLANE
174C---------------------------
175#include "vectorize.inc"
176 DO j=1,nindx
177 i = index(j)
178 nu1(i) = one/(one-nu(i))
179 nu2(i) = one/(one+nu(i))
180 h(i) = cn(i)*cb(i)*exp((cn(i)-one)*log(pla(i)+small))
181 IF (yld(i) >= ymax(i)) h(i) = zero
182 dpla_j(i) = (svm(i)-yld(i))/(three*g(i)+h(i))
183 etse(i) = h(i)/(h(i)+young(i))
184 ENDDO
185C
186 DO n=1,nmax
187#include "vectorize.inc"
188 DO j=1,nindx
189 i = index(j)
190 dpla_i(i) = dpla_j(i)
191 dpla(i) = dpla_i(i)
192 pla_i = pla(i)+dpla_i(i)
193 yld_i = min(ymax(i),ca(i)+cb(i)*pla_i**cn(i))
194 dr(i) = half*young(i)*dpla_i(i)/yld_i
195 p(i) = one/(one+dr(i)*nu1(i))
196 q(i) = one/(one + three*dr(i)*nu2(i))
197 p2 = p(i)*p(i)
198 q2 = q(i)*q(i)
199 f = a(i)*p2+b(i)*q2-yld_i*yld_i
200 df = -(a(i)*nu1(i)*p2*p(i)+3.*b(i)*nu2(i)*q2*q(i))
201 . *(young(i)-two*dr(i)*h(i))/yld_i
202 . -two*h(i)*yld_i
203 IF (dpla_i(i) > zero) THEN
204 dpla_j(i) = max(zero,dpla_i(i)-f/df)
205 ELSE
206 dpla_j(i) = zero
207 ENDIF
208 ENDDO
209 ENDDO
210C------------------------------------------
211C CONTRAINTES PLASTIQUEMENT ADMISSIBLES
212C------------------------------------------
213#include "vectorize.inc"
214 DO j=1,nindx
215 i = index(j)
216 pla(i) = pla(i) + dpla_i(i)
217 s1 = (sig(i,1)+sig(i,2))*p(i)
218 s2 = (sig(i,1)-sig(i,2))*q(i)
219 sig(i,1) = half*(s1+s2)
220 sig(i,2) = half*(s1-s2)
221 sig(i,3) = sig(i,3)*q(i)
222 ezz(i) = dr(i)*s1/young(i)
223 ENDDO
224C-------------------------
225 ELSEIF (ipla == 2) THEN
226C-------------------
227C CONTRAINTE VM
228C-------------------
229 DO i=jft,jlt
230 pp = -(sig(i,1)+sig(i,2))*third
231 s11 = sig(i,1)+pp
232 s22 = sig(i,2)+pp
233 s12 = sig(i,3)
234 p2 = pp*pp
235 s1s2 = s11*s22
236 s122 = s12*s12
237 nnu1 = nu(i) / (one - nu(i))
238 nnu2 = nnu1*nnu1
239 nu4 = one + nnu2 + nnu1
240 nu6 = half - nnu2 + half*nnu1
241 qq = (one-nnu1)*pp
242 aa = p2*nu4 + three*(s122 - s1s2)
243 bb = p2*nu6
244 c = qq*qq
245 vm2= aa+bb+bb+c
246 c = c - yld(i)*yld(i)
247C
248 r = min(one,(-bb+ sqrt(max(zero,bb*bb-aa*c)))/max(aa ,em20))
249C
250 umr = one - r
251 qq = qq*umr
252 sig(i,1) = sig(i,1)*r - qq
253 sig(i,2) = sig(i,2)*r - qq
254 sig(i,3) = s12*r
255 dpla(i) = off(i)*sqrt(vm2)*umr/(three*g(i))
256 h(i) = cn(i)*cb(i)*exp((cn(i)-one)*log(pla(i)+small))
257 IF (yld(i) >= ymax(i)) h(i) = zero
258 etse(i) = h(i)/(h(i)+young(i))
259 s1 = half*(sig(i,1)+sig(i,2))
260 ezz(i) = dpla(i) * s1 / yld(i)
261 pla(i) = pla(i) + dpla(i)
262 ENDDO
263 ENDIF
264C
265 RETURN
266 END
m27plas
subroutine m27plas(jft, jlt, pm, off, sig, pla, imat, dt1, ipla, ezz, epsd_pg, israte, yld, etse, dpla, depsxx, depsyy, depsxy, depsyz, depszx, nel, epsd, asrate)
Definition m27plas.F:33
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21