OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
m4law.F File Reference
#include "implicit_f.inc"
#include "com08_c.inc"
#include "param_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine m4law (pm, off, sig, epxe, mat, ssp, vol, d1, d2, d3, d4, d5, d6, rho0, dpdm, epd, ipla, sigy, defp, dpla, epsp, tstar, tempel, nel, jthe, fheat, jlag)

Function/Subroutine Documentation

◆ m4law()

subroutine m4law ( dimension(npropm,*), intent(in) pm,
dimension(nel) off,
sig,
dimension(nel) epxe,
integer, dimension(nel), intent(in) mat,
dimension(nel) ssp,
dimension(nel), intent(in) vol,
dimension(nel) d1,
dimension(nel) d2,
dimension(nel) d3,
dimension(nel) d4,
dimension(nel) d5,
dimension(nel) d6,
dimension(nel) rho0,
dimension(nel) dpdm,
dimension(nel) epd,
integer, intent(in) ipla,
dimension(nel) sigy,
dimension(nel) defp,
dimension(nel) dpla,
dimension(nel) epsp,
intent(in) tstar,
intent(inout) tempel,
integer, intent(in) nel,
integer, intent(in) jthe,
intent(inout) fheat,
integer, intent(in) jlag )

Definition at line 28 of file m4law.F.

36C-----------------------------------------------
37C I m p l i c i t T y p e s
38C-----------------------------------------------
39#include "implicit_f.inc"
40C-----------------------------------------------
41C C o m m o n B l o c k s
42C-----------------------------------------------
43#include "com08_c.inc"
44#include "param_c.inc"
45C-----------------------------------------------
46C D u m m y A r g u m e n t s
47C-----------------------------------------------
48 INTEGER, INTENT(IN) :: JTHE
49 INTEGER, INTENT(IN) :: NEL
50 INTEGER, INTENT(IN) :: IPLA
51 INTEGER, INTENT(IN) :: MAT(NEL)
52 INTEGER, INTENT(IN) :: JLAG
53 my_real, INTENT(IN) :: pm(npropm,*)
54 my_real, INTENT(IN) :: vol(nel)
55 my_real :: sig(nel,6)
56 my_real ,DIMENSION(NEL) :: off,epxe,ssp,d1,d2,d3,d4,d5,d6,rho0,dpdm,epd,sigy,defp,dpla,epsp
57 my_real ,DIMENSION(NEL) ,INTENT(IN) :: tstar
58 my_real ,DIMENSION(NEl) ,INTENT(INOUT) :: tempel
59 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: fheat ! Heat due to plastic work for Heat Equation with lagrangian framework
60C-----------------------------------------------
61C L o c a l V a r i a b l e s
62C-----------------------------------------------
63 INTEGER :: I,MX
64 my_real :: rhocp,tmax,ct, ce, ch
65 my_real g(nel) ,g1(nel) ,g2(nel) ,qs(nel) ,ak(nel),
66 . qh(nel) ,tmelt(nel),aj2(nel) ,dav(nel),p(nel) ,
67 . epmx(nel) ,ca(nel) ,cb(nel) ,cc(nel) ,
68 . cn(nel) ,epxo(nel) ,epdr(nel),cmx(nel),sigmx(nel),
69 . scale(nel),t0(nel)
70C-----------------------------------------------
71C B o d y
72C-----------------------------------------------
73 mx = mat(1)
74 tmax = pm(47,mx)
75 rhocp = pm(69,mx)
76 DO i=1,nel
77 g(i) =pm(22,mx)
78 ca(i) =pm(38,mx)
79 cb(i) =pm(39,mx)
80 cn(i) =pm(40,mx)
81 epmx(i) =pm(41,mx)
82 sigmx(i)=pm(42,mx)
83 cc(i) =pm(43,mx)
84 epdr(i) =pm(44,mx)
85 cmx(i) =pm(45,mx)
86 tmelt(i)=pm(46,mx)
87 t0(i) =pm(79,mx)
88 ENDDO
89C
90 DO i=1,nel
91 p(i) =-third*(sig(i,1)+sig(i,2)+sig(i,3))
92 dav(i)=-third*(d1(i)+d2(i)+d3(i))
93 g1(i)=dt1*g(i)*off(i)
94 g2(i)=two*g1(i)*off(i)
95 ENDDO
96C-----------------------
97C SOUND SPEED
98C-----------------------
99 DO i=1,nel
100 dpdm(i) = dpdm(i) + onep333*g(i)
101 ssp(i)=sqrt(abs(dpdm(i))/rho0(i))
102 ENDDO
103C-------------------------------
104C CONTRAINTES DEVIATORIQUES
105C-------------------------------
106 DO i=1,nel
107 sig(i,1)=sig(i,1)+p(i)+g2(i)*(d1(i)+dav(i))
108 sig(i,2)=sig(i,2)+p(i)+g2(i)*(d2(i)+dav(i))
109 sig(i,3)=sig(i,3)+p(i)+g2(i)*(d3(i)+dav(i))
110 sig(i,4)=sig(i,4)+g1(i)*d4(i)
111 sig(i,5)=sig(i,5)+g1(i)*d5(i)
112 sig(i,6)=sig(i,6)+g1(i)*d6(i)
113 ENDDO
114C
115 DO i=1,nel
116 epd(i)=off(i)*max( abs(d1(i)), abs(d2(i)), abs(d3(i)),
117 . half*abs(d4(i)),half*abs(d5(i)),half*abs(d6(i)))
118C
119 epsp(i) = epd(i)
120 aj2(i)=half*(sig(i,1)**2+sig(i,2)**2+sig(i,3)**2)
121 . +sig(i,4)**2+sig(i,5)**2+sig(i,6)**2
122 aj2(i)=sqrt(three*aj2(i))
123 ENDDO
124C-------------------------------
125C CRITERE
126C-------------------------------
127 DO i=1,nel
128 ct = one
129 IF (tempel(i) >= tmelt(i)) THEN
130 aj2(i)=zero
131 qh(i)=zero
132 ak(i)=zero
133 scale(i)=zero
134 cycle
135 ELSEIF (tempel(i) > t0(i)) THEN
136 IF (tempel(i) > tmax) cmx(i)=one
137 ct = one - tstar(i)**cmx(i)
138 ENDIF
139C
140 IF(epd(i)<=epdr(i)) THEN
141 ce=one
142 ELSE
143 ce=one + cc(i) * log(epd(i)/epdr(i))
144 ENDIF
145C
146 IF(epxe(i)<=zero) THEN
147 ch=ca(i)
148 ELSEIF(epxe(i)>epmx(i)) THEN
149 ch=zero
150 ELSE
151 ch=ca(i)+cb(i)*epxe(i)**cn(i)
152 ENDIF
153C
154 ak(i) = min(sigmx(i),ch)*ce*ct
155C-----------------------
156C MODULE ECROUISSAGE
157C-----------------------
158 IF(cn(i)>=one) THEN
159 qh(i)= (cb(i)*cn(i)*epxe(i)**(cn(i) - one))*ce*ct
160 ELSEIF(epxe(i)>zero)THEN
161 qh(i)= (cb(i)*cn(i)/epxe(i)**(one -cn(i)))*ce*ct
162 ELSE
163 qh(i)=zero
164 ENDIF
165C-------------------------------
166C PROJECTION
167C-------------------------------
168 IF(aj2(i)<=ak(i)) THEN
169 scale(i)=1.
170 ELSEIF(aj2(i)/=zero) THEN
171 scale(i)=ak(i)/aj2(i)
172 ELSE
173 scale(i)=zero
174 ENDIF
175 ENDDO
176C
177 IF(ipla==0)THEN
178 DO i=1,nel
179 sig(i,1)=scale(i)*sig(i,1)
180 sig(i,2)=scale(i)*sig(i,2)
181 sig(i,3)=scale(i)*sig(i,3)
182 sig(i,4)=scale(i)*sig(i,4)
183 sig(i,5)=scale(i)*sig(i,5)
184 sig(i,6)=scale(i)*sig(i,6)
185 dpla(i) =(one -scale(i))*aj2(i)/(three*g(i)+qh(i))
186 epxe(i)=epxe(i)+ dpla(i)
187 epxe(i)=epxe(i)*off(i)
188 ENDDO
189C
190 ELSEIF(ipla==2)THEN
191 DO i=1,nel
192 sig(i,1)=scale(i)*sig(i,1)
193 sig(i,2)=scale(i)*sig(i,2)
194 sig(i,3)=scale(i)*sig(i,3)
195 sig(i,4)=scale(i)*sig(i,4)
196 sig(i,5)=scale(i)*sig(i,5)
197 sig(i,6)=scale(i)*sig(i,6)
198 dpla(i) =(one -scale(i))*aj2(i)/(three*g(i))
199 epxe(i)=epxe(i)+dpla(i)
200 epxe(i)=epxe(i)*off(i)
201 ENDDO
202C
203 ELSEIF(ipla==1)THEN
204 DO i=1,nel
205C plastic strain increment.
206 dpla(i)=(one -scale(i))*aj2(i)/(three*g(i)+qh(i))
207C actual yield stress.
208 ak(i)=ak(i)+dpla(i)*qh(i)
209 scale(i)= min(one,ak(i)/ max(aj2(i),em15))
210 sig(i,1)=scale(i)*sig(i,1)
211 sig(i,2)=scale(i)*sig(i,2)
212 sig(i,3)=scale(i)*sig(i,3)
213 sig(i,4)=scale(i)*sig(i,4)
214 sig(i,5)=scale(i)*sig(i,5)
215 sig(i,6)=scale(i)*sig(i,6)
216 epxe(i)=epxe(i)+dpla(i)
217 epxe(i)=epxe(i)*off(i)
218 ENDDO
219 ENDIF
220C
221 DO i=1,nel
222 sigy(i)=ak(i)
223 defp(i)=epxe(i)
224 ENDDO
225C----------------------------------------------
226C TEMPERATURE (Heating due to plastic work)
227C----------------------------------------------
228 IF (jthe /= 0 .AND. jlag /= 0) THEN
229 DO i=1,nel
230 fheat(i) = fheat(i) + sigy(i)*dpla(i)*vol(i)
231 ENDDO
232 ELSEIF(rhocp > zero)THEN
233 DO i=1,nel
234 tempel(i) = tempel(i) + sigy(i)*dpla(i) / rhocp
235 ! temperature and internal energy must be incremented consistantly
236 ! internal energy is incremented later in parent subroutine (mmain)
237 ! with total energy deformation which already includes plastic work
238 ! Edef = 0.5 * VOL * sum ( sig.eps_dot , i=1..6)
239 ! so internal energy and temperature remain consistant
240 ENDDO
241 END IF
242c----------
243 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21