46
47
48
49#include "implicit_f.inc"
50
51
52
53#include "mvsiz_p.inc"
54
55
56
57#include "scr05_c.inc"
58#include "impl1_c.inc"
59
60
61
62 INTEGER NEL, NUPARAM, NUVAR,ISMSTR,NGL(*),IHET,IEXPAN
64 . time , timestep , uparam(nuparam),
65 . rho(nel), rho0(nel), volume(nel), eint(nel),
66 . epspxx(nel), epspyy(nel), epspzz(nel),
67 . epspxy(nel), epspyz(nel), epspzx(nel),
68 . depsxx(nel), depsyy(nel), depszz(nel),
69 . depsxy(nel), depsyz(nel), depszx(nel),
70 . epsxx(nel), epsyy(nel), epszz(nel),
71 . epsxy(nel), epsyz(nel), epszx(nel),
72 . sigoxx(nel), sigoyy(nel), sigozz(nel),
73 . sigoxy(nel), sigoyz(nel), sigozx(nel),
74 . mfxx(nel) , mfxy(nel), mfxz(nel),
75 . mfyx(nel) , mfyy(nel), mfyz(nel),
76 . mfzx(nel) , mfzy(nel), mfzz(nel),offg(nel),
77 . epsth3(nel)
78
79
80
82 . signxx(nel), signyy(nel), signzz(nel),
83 . signxy(nel), signyz(nel), signzx(nel),
84 . sigvxx(nel), sigvyy(nel), sigvzz(nel),
85 . sigvxy(nel), sigvyz(nel), sigvzx(nel),
86 . soundsp(nel), viscmax(nel), et(nel)
87
88
89
91 . uvar(nel,nuvar), off(nel)
92
93
94
95 INTEGER NPF(*), NFUNC, IFUNC(NFUNC)
96 my_real finter,fintte,tf(*),fint2v
97 EXTERNAL finter,fintte
98
99
100
101 INTEGER I,J,K,II,NORDRE,NMAXW,JJ
103 . mu(100),al(100),s(mvsiz,3),dd,pp,gmax,
104 . d(mvsiz),av(mvsiz,6),rbulk,ev(mvsiz,3),
105 . evd(mvsiz,3),rv(mvsiz),evv(mvsiz,3),
106 . dirprv(mvsiz,3,3),pui,pui_tab(3),lam_al(3)
108 . a(3,3),dpra(3,3),eigen(3) ,amax ,eti,p,efac
109 my_real ,
DIMENSION(NEL,3) :: cii
110 my_real ,
DIMENSION(NEL) :: gtmax,rkmax
111
112
113
114 gmax = zero
115
116 nordre = nint(uparam(1))
117 DO i= 1,nordre
118 mu(i) = uparam(1 + i)
119 al(i) = uparam(1 + nordre + i)
120 d(i) = uparam(1 + 2*nordre + i)
121 gmax = gmax + mu(i)
122
123 ENDDO
124
125
126 rbulk = two/d(1)
127
128 IF(time==zero)THEN
129 DO j = 1,nuvar
130 DO i = 1, nel
131 uvar(i,j) = zero
132 ENDDO
133 ENDDO
134 ENDIF
135
136
137 DO i=1,nel
138 av(i,1)=epsxx(i)
139 av(i,2)=epsyy(i)
140 av(i,3)=epszz(i)
141 av(i,4)=epsxy(i)/2
142 av(i,5)=epsyz(i)/2
143 av(i,6)=epszx(i)/2
144 ENDDO
145
146
147 IF (iresp==1) THEN
149 ELSE
151 ENDIF
152
153
154
155
156 IF(ismstr==0.OR.ismstr==2.OR.ismstr==4) THEN
157 DO i=1,nel
158
159 ev(i,1)=exp(evv(i,1))
160 ev(i,2)=exp(evv(i,2))
161 ev(i,3)=exp(evv(i,3))
162 ENDDO
163 ELSEIF(ismstr==10.OR.ismstr==12) THEN
164 DO i =1,nel
165 IF(offg(i)<=one) THEN
166 ev(i,1)=sqrt(evv(i,1)+ one)
167 ev(i,2)=sqrt(evv(i,2)+ one)
168 ev(i,3)=sqrt(evv(i,3)+ one)
169 ELSE
170 ev(i,1)=evv(i,1)+ one
171 ev(i,2)=evv(i,2)+ one
172 ev(i,3)=evv(i,3)+ one
173 END IF
174 ENDDO
175 ELSE
176
177 DO i=1,nel
178 ev(i,1)=evv(i,1)+ one
179 ev(i,2)=evv(i,2)+ one
180 ev(i,3)=evv(i,3)+ one
181 ENDDO
182 ENDIF
183
184 IF (impl_s > 0 .OR. ihet > 1) THEN
185 DO j =1,3
186 DO i=1,nel
187 amax = ev(i,j)
188 eti=zero
189 IF(amax/=zero) THEN
190 DO k= 1,nordre
191 eti = eti + mu(k)*amax**(al(k)-one)
192 ENDDO
193 ENDIF
194 et(i)=
max(eti,et(i))
195 ENDDO
196 ENDDO
197
198 DO i=1,nel
199 et(i)=
max(one,et(i)/gmax)
200 ENDDO
201 ENDIF
202
203
204
205 DO i=1,nel
206
207 rv(i) = ev(i,1)*ev(i,2)*ev(i,3)
208 ENDDO
209
210 IF(iexpan > 0.AND.(ismstr==10.OR.ismstr==11.OR.ismstr==12)) THEN
211 DO i=1,nel
212 rv(i) = rv(i) -epsth3(i)
213 ENDDO
214 ENDIF
215
216 DO i=1,nel
217 IF(rv(i)/=zero) THEN
218 pui = rv(i)**(-third)
219 ELSE
220 pui = zero
221 ENDIF
222 evd(i,1) = ev(i,1)*pui
223 evd(i,2) = ev(i,2)*pui
224 evd(i,3) = ev(i,3)*pui
225 ENDDO
226 DO i=1,nel
227 s(i,1) = zero
228 s(i,2) = zero
229 s(i,3) = zero
230 p = zero
231 DO j=1,nordre
232 dd = two*mu(j)/al(j)
233 dd = dd /rv(i)
234 pp = zero
235 IF(al(j)/=zero) THEN
236 DO ii=1,3
237 IF(evd(i,ii)/=zero) THEN
238 pui_tab(ii) = evd(i,ii)**al(j)
239 ELSE
240 pui_tab(ii) = zero
241 ENDIF
242 ENDDO
243 ELSE
244 DO ii=1,3
245 IF(evd(i,ii)/=zero) THEN
246 pui_tab(ii) = one
247 ELSE
248 pui_tab(ii) = zero
249 ENDIF
250 ENDDO
251 ENDIF
252
253
254 IF(d(j) /= zero)pp = one/d(j)
255 s(i,1) = s(i,1) + dd*(two_third*pui_tab(1)
256 . - third*(pui_tab(2) + pui_tab(3)))
257
258 s(i,2) = s(i,2) + dd*(two_third*pui_tab(2)
259 . - third*(pui_tab(1) + pui_tab(3)))
260
261 s(i,3) = s(i,3) + dd*(two_third*pui_tab(3)
262 . - third*(pui_tab(2) + pui_tab(1)))
263
264 p = p + 2*j*pp*(rv(i) - one)**(2*j - 1)
265 ENDDO
266
267 s(i,1) = s(i,1) + p
268 s(i,2) = s(i,2) + p
269 s(i,3) = s(i,3) + p
270 ENDDO
271
272
273 gtmax(1:nel) = gmax
274 rkmax(1:nel) = rbulk
275
276 cii(1:nel,1:3) = zero
277 DO ii = 1,nordre
278 IF(mu(ii)/=zero) THEN
279 DO i=1,nel
280 lam_al(1:3) = exp(al(ii)*log(evd(i,1:3)))
281 amax = third*(lam_al(1)+lam_al(2)+lam_al(3))
282 cii(i,1:3) = cii(i,1:3) + mu(ii) * (lam_al(1:3)+amax)
283 ENDDO
284 ENDIF
285 ENDDO
286 DO ii = 2,nordre
287 IF (d(ii)==zero) cycle
288 pp = 2*ii*(2*ii-1)/d(ii)
289 jj = 2*ii - 2
290 DO i=1,nel
291 IF (abs(rv(i)- one)<em20) cycle
292 rkmax(i) = rkmax(i) + pp * (rv(i)- one)**jj
293 ENDDO
294 ENDDO
295
296 DO i=1,nel
297 amax= half*
max(cii(i,1),cii(i,2),cii(i,3))/gmax
298
299 eti =
max(one,amax*0.81)
300 gtmax(i) = gmax*eti
301 rkmax(i) =
max(rbulk,rkmax(i))
302
304 ENDDO
305
306 DO i=1,nel
307 signxx(i) = dirprv(i,1,1)*dirprv(i,1,1)*s(i,1)
308 . + dirprv(i,1,2)*dirprv(i,1,2)*s(i,2)
309 . + dirprv(i,1,3)*dirprv(i,1,3)*s(i,3)
310
311 signyy(i) = dirprv(i,2,2)*dirprv(i,2,2)*s(i,2)
312 . + dirprv(i,2,3)*dirprv(i,2,3)*s(i,3)
313 . + dirprv(i,2,1)*dirprv(i,2,1)*s(i,1)
314
315 signzz(i) = dirprv(i,3,3)*dirprv(i,3,3)*s(i,3)
316 . + dirprv(i,3,1)*dirprv(i,3,1)*s(i,1)
317 . + dirprv(i,3,2)*dirprv(i,3,2)*s(i,2)
318
319 signxy(i) = dirprv(i,1,1)*dirprv(i,2,1)*s(i,1)
320 . + dirprv(i,1,2)*dirprv(i,2,2)*s(i,2)
321 . + dirprv(i,1,3)*dirprv(i,2,3)*s(i,3)
322
323 signyz(i) = dirprv(i,2,2)*dirprv(i,3,2)*s(i,2)
324 . + dirprv(i,2,3)*dirprv(i,3,3)*s(i,3)
325 . + dirprv(i,2,1)*dirprv(i,3,1)*s(i,1)
326
327 signzx(i) = dirprv(i,3,3)*dirprv(i,1,3)*s(i,3)
328 . + dirprv(i,3,1)*dirprv(i,1,1)*s(i,1)
329 . + dirprv(i,3,2)*dirprv(i,1,2)*s(i,2)
330
331
332 soundsp(i)=sqrt((four_over_3*gtmax(i) + rkmax(i))/rho(i)
333
334 viscmax(i) = zero
335 ENDDO
336
337 RETURN
subroutine valpvecdp_v(sig, val, vec, nel)
subroutine valpvec_v(sig, val, vec, nel)