OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
conversion_mod Module Reference

Functions/Subroutines

subroutine conversion (zxx, zyy, zzz, zxy, zyz, zzx, cijkl, dfirst, drate, xscal, scal, nel, jac, slope, numtabl, table, nvartmp, vartmp, yld, ndim_table)

Function/Subroutine Documentation

◆ conversion()

subroutine conversion_mod::conversion ( dimension(nel), intent(in) zxx,
dimension(nel), intent(in) zyy,
dimension(nel), intent(in) zzz,
dimension(nel), intent(in) zxy,
dimension(nel), intent(in) zyz,
dimension(nel), intent(in) zzx,
dimension(nel,6,6), intent(out) cijkl,
dimension(nel,3,6), intent(out) dfirst,
dimension(nel,3), intent(inout) drate,
intent(in) xscal,
intent(in) scal,
integer, intent(in) nel,
dimension(nel), intent(in) jac,
dimension(nel,3), intent(out) slope,
integer, intent(in) numtabl,
type (table_4d_), dimension(numtabl), target table,
integer, intent(in) nvartmp,
integer, dimension(nel,nvartmp), intent(inout) vartmp,
dimension(nel), intent(out) yld,
integer, intent(in) ndim_table )

Definition at line 41 of file conversion.F.

46c computation of coefficients to convert stiffness
47c matrix from local to global reference system
48c===========================================================================
49C M o d u l e s
50C-----------------------------------------------
51 USE table4d_mod
53C-----------------------------------------------
54C I m p l i c i t T y p e s
55C-----------------------------------------------
56#include "implicit_f.inc"
57#include "comlock.inc"
58C-----------------------------------------------
59#include "scr05_c.inc"
60#include "impl1_c.inc"
61#include "com01_c.inc"
62#include "mvsiz_p.inc"
63c===========================================================================
64C-----------------------------------------------
65C I N T P U T A r g u m e n t s
66C-----------------------------------------------
67 INTEGER, INTENT(IN) :: NEL, NUMTABL,NVARTMP,NDIM_TABLE
68 my_real, INTENT(IN) :: scal,zxx(nel),zyy(nel),zzz(nel),zxy(nel),
69 . zzx(nel),zyz(nel), jac(nel), xscal
70 INTEGER, DIMENSION(NEL,NVARTMP), INTENT(INOUT) :: VARTMP
71C-----------------------------------------------
72C I N T P U T O U T P U T A r g u m e n t s
73C-----------------------------------------------
74 my_real, INTENT(OUT) :: dfirst(nel,3,6),cijkl(nel,6,6),slope(nel,3),yld(nel)
75 my_real, INTENT(INOUT) :: drate(nel,3)
76 TYPE (TABLE_4D_), DIMENSION(NUMTABL) ,TARGET :: TABLE
77c===========================================================================
78C-----------------------------------------------
79C L o c a l V a r i a b l e s
80C-----------------------------------------------
81 INTEGER I,J,K,L
82
83 my_real, DIMENSION(NEL) :: dydx1,volstr,
84 . fsig1,dsig1,fsig2,dsig2,
85 . sigener,dsigener,fsig3,dsig3
86 my_real, DIMENSION(NEL,3) :: alambda,epsilon,engirate, sigma, dynsigma,
87 . strainrate
88 my_real
89 . dik(nel,3), sik(nel,3),
90 . ac,as,z1(nel),z2(nel),z3(nel),dsecond(nel,3,6,6),
91 . ah1(nel),ah2(nel),ah3(nel),j2(nel),j3(nel),aa3(nel),aa2(nel),
92 . dzero(nel,3),dat2(nel,6),dat3(nel,6),da(nel,6),da2(nel,6,6),dat22(nel,6,6),
93 . dat32(nel,6,6),dx2(nel,6,6),
94 . dh1zxx,dh1zyy,dh1zxy,dh1zyz,dh1zzx,dh1zzz,
95 . dh2zxx,dh2zyy,dh2zxy,dh2zyz,dh2zzx,dh2zzz,
96 . dh3zxx,dh3zyy,dh3zxy,dh3zyz,dh3zzx,dh3zzz,
97 . dt2zxx,dt2zyy,dt2zxy,dt2zyz,dt2zzx,dt2zzz,
98 . dt3zxx,dt3zyy,dt3zxy,dt3zyz,dt3zzx,dt3zzz,
99 . temp1,temp2,temp3,temp4,temp5,temp11,temp12,
100 . dz1zxx,dz1zyy,dz1zxy,dz1zyz,dz1zzx,dz1zzz,
101 . dz2zxx,dz2zyy,dz2zxy,dz2zyz,dz2zzx,dz2zzz,
102 . dz3zxx,dz3zyy,dz3zxy,dz3zyz,dz3zzx,dz3zzz,
103 . da3zxx,da3zyy,da3zxy,da3zyz,da3zzx,da3zzz,
104 . dah2zxxzxx,dah2zxxzyy,dah2zxxzzz,dah2zxxzxy,dah2zxxzyz,dah2zxxzzx,
105 . dah2zyyzxx,dah2zyyzyy,dah2zyyzzz,dah2zyyzxy,dah2zyyzyz,dah2zyyzzx,
106 . dah2zzzzxx,dah2zzzzyy,dah2zzzzzz,dah2zzzzxy,dah2zzzzyz,dah2zzzzzx,
107 . dah2zxyzxx,dah2zxyzyy,dah2zxyzzz,dah2zxyzxy,dah2zxyzyz,dah2zxyzzx,
108 . dah2zyzzxx,dah2zyzzyy,dah2zyzzzz,dah2zyzzxy,dah2zyzzyz,dah2zyzzzx,
109 . dah2zzxzxx,dah2zzxzyy,dah2zzxzzz,dah2zzxzxy,dah2zzxzyz,dah2zzxzzx,
110 . dah3zxxzxx,dah3zxxzyy,dah3zxxzzz,dah3zxxzxy,dah3zxxzyz,dah3zxxzzx,
111 . dah3zyyzxx,dah3zyyzyy,dah3zyyzzz,dah3zyyzxy,dah3zyyzyz,dah3zyyzzx,
112 . dah3zzzzxx,dah3zzzzyy,dah3zzzzzz,dah3zzzzxy,dah3zzzzyz,dah3zzzzzx,
113 . dah3zxyzxx,dah3zxyzyy,dah3zxyzzz,dah3zxyzxy,dah3zxyzyz,dah3zxyzzx,
114 . dah3zyzzxx,dah3zyzzyy,dah3zyzzzz,dah3zyzzxy,dah3zyzzyz,dah3zyzzzx,
115 . dah3zzxzxx,dah3zzxzyy,dah3zzxzzz,dah3zzxzxy,dah3zzxzyz,dah3zzxzzx,
116 . x1,x2,x22,x3,x4,x5,x6,
117 . t1,t2,t3,t4,t5,t6,
118 . afac1,afac2,afac3,afac4,afac33,afac5,
119 . ahelp3,ahelp33,ahelp,ahelp4,teta,
120 . xcap,dxcap,three_teta,termdat32,dacos_daa3,d2acos_daa32,
121 . fac1_da2, fac2_da2,fac3_da2,t1da2,t2da2,t3da2,t4da2,t5da2,t6da2,
122 . daa3dt2,daa3dt3,d2a_dt22,d2a_dt2dt3,d2a_dt3dt2
123 ! TABLE
124 INTEGER IPOS(NEL,3),IPOS2(NEL,2),IPOS1(NEL,1)
125 my_real, DIMENSION(NEL,3) :: xvec
126 my_real, DIMENSION(NEL,2) :: xvec2
127 my_real, DIMENSION(NEL,1) :: xvec1
128 TYPE(TABLE_4D_), POINTER :: FUNC_SIG
129C===========================================================================
130C===========================================================================
131 func_sig => table(1)
132
133
134 DO i=1,nel
135C
136C PART 1 : PRINCIPAL STRAINS
137C
138C COMPUTE H1, H2 AND H3
139C
140 ah1(i) = zxx(i) + zyy(i) + zzz(i)
141 ah2(i) = zyy(i)*zzz(i)+zzz(i)*zxx(i)+zxx(i)*zyy(i)
142 . - zyz(i)**2 - zzx(i)**2 - zxy(i)**2
143 ah3(i) = zxx(i)*zyy(i)*zzz(i)+two*zxy(i)*zyz(i)*zzx(i)
144 . - zxx(i)*zyz(i)**2-zyy(i)*zzx(i)**2-zzz(i)*zxy(i)**2
145C
146C DERIVATIVES OF H1
147C
148 dh1zxx=one
149 dh1zyy=one
150 dh1zzz=one
151 dh1zxy=zero
152 dh1zyz=zero
153 dh1zzx=zero
154C
155C DERIVATIVES OF H2
156C
157 dh2zxx=zyy(i)+zzz(i)
158 dh2zyy=zxx(i)+zzz(i)
159 dh2zzz=zyy(i)+zxx(i)
160 dh2zxy=-two*zxy(i)
161 dh2zyz=-two*zyz(i)
162 dh2zzx=-two*zzx(i)
163C
164C DERIVATIVES OF H3
165C
166 dh3zxx = zyy(i)*zzz(i)-zyz(i)**2
167 dh3zyy = zxx(i)*zzz(i)-zzx(i)**2
168 dh3zzz = zyy(i)*zxx(i)-zxy(i)**2
169 dh3zxy = -two*zxy(i)*zzz(i)+two*zyz(i)*zzx(i)
170 dh3zyz = -two*zyz(i)*zxx(i)+two*zxy(i)*zzx(i)
171 dh3zzx = -two*zzx(i)*zyy(i)+two*zyz(i)*zxy(i)
172C
173C COMPUTE T2 AND T3
174C
175 j3(i) = (two/twenty7)*ah1(i)**3 - ah1(i)*ah2(i)*third + ah3(i)
176 j2(i) = third*ah1(i)**2 - ah2(i)
177 j2(i) = max(j2(i),em16)
178 ! J2 = -3*Q
179 ! J3 = 2*R
180 ! D = Q**3 + R**2 must be negative or nul to have real solutions and define teta
181C
182C COMPUTE COS(A)
183C
184 aa3(i) = j3(i)/two/sqrt(j2(i)**3/twenty7)
185 aa3(i) = min( one,aa3(i))
186 aa3(i) = max(-one,aa3(i))
187 aa2(i) = sqrt(j2(i)*third)
188 three_teta = acos(aa3(i))
189C
190C COMPUTE PRINCIPAL STRAINS
191C
192 z1(i) = two*aa2(i)*cos(third*three_teta) + ah1(i)*third
193 z2(i) = two*aa2(i)*cos(third*three_teta - two *pi*third) + ah1(i)*third
194 z3(i) = two*aa2(i)*cos(third*three_teta - four*pi*third) + ah1(i)*third
195C
196C PART 2
197C FIRST DERIVATIVES OF PRINCIPAL STRAINS
198C
199C DERIVATIVES OF T2
200C
201 dt2zxx = two*ah1(i)*third-(zyy(i)+zzz(i))
202 dt2zyy = two*ah1(i)*third-(zxx(i)+zzz(i))
203 dt2zzz = two*ah1(i)*third-(zyy(i)+zxx(i))
204 dt2zxy = two*zxy(i)
205 dt2zyz = two*zyz(i)
206 dt2zzx = two*zzx(i)
207C
208C STORE DERIVATIVES OF T2
209C
210 dat2(i,1)=dt2zxx
211 dat2(i,2)=dt2zyy
212 dat2(i,3)=dt2zzz
213 dat2(i,4)=dt2zxy
214 dat2(i,5)=dt2zyz
215 dat2(i,6)=dt2zzx
216C
217C DERIVATIVES OF T3
218C
219 dt3zxx = six * ah1(i)**2/twenty7 - ah2(i)*third
220 . - ah1(i)*(zyy(i)+zzz(i))*third
221 . +(zyy(i)*zzz(i)-zyz(i)**2)
222 dt3zyy = six * ah1(i)**2/twenty7-ah2(i)*third
223 . - ah1(i)*(zxx(i)+zzz(i))*third
224 . +(zxx(i)*zzz(i)-zzx(i)**2)
225 dt3zzz = six * ah1(i)**2/twenty7-ah2(i)*third
226 . - ah1(i)*(zyy(i)+zxx(i))*third
227 . +(zyy(i)*zxx(i)-zxy(i)**2)
228 dt3zxy = two*ah1(i)*zxy(i)*third
229 . +two*(zyz(i)*zzx(i)-zxy(i)*zzz(i))
230 dt3zyz = two*ah1(i)*zyz(i)*third
231 . +two*(zxy(i)*zzx(i)-zyz(i)*zxx(i))
232 dt3zzx = two*ah1(i)*zzx(i)*third
233 . +two*(zyz(i)*zxy(i)-zzx(i)*zyy(i))
234C
235C STORE DERIVATIVES OF T3
236C
237 dat3(i,1)=dt3zxx
238 dat3(i,2)=dt3zyy
239 dat3(i,3)=dt3zzz
240 dat3(i,4)=dt3zxy
241 dat3(i,5)=dt3zyz
242 dat3(i,6)=dt3zzx
243C
244C AUXILIARY COEFFICIENTS
245C
246
247C
248C 'IF' TEST NEEDED TO GET CORRECT LIMIT VALUE
249C IN THE UNIAXIAL CASES
250C BIAXIAL CASES SEEM FINE WITHOUT THE TEST
251C
252
253C
254C DERIVATIVES OF A
255C
256 dacos_daa3 = -one - (aa3(i)**2)/two - (aa3(i)**4)/eight -
257 . (aa3(i)**6)/48.0d0 - (aa3(i)**8)/384.0d0 -
258 . (aa3(i)**10)/(3840.0d0)
259 daa3dt3 = one/(two*sqrt((j2(i)**3)/twenty7))
260 daa3dt2 = -(j3(i)/36.0d0)*(j2(i)**2)/(sqrt(((j2(i)**3)/twenty7)**3))
261
262 da3zxx = dacos_daa3*(daa3dt3*dt3zxx + daa3dt2*dt2zxx)
263 da3zyy = dacos_daa3*(daa3dt3*dt3zyy + daa3dt2*dt2zyy)
264 da3zzz = dacos_daa3*(daa3dt3*dt3zzz + daa3dt2*dt2zzz)
265 da3zxy = dacos_daa3*(daa3dt3*dt3zxy + daa3dt2*dt2zxy)
266 da3zyz = dacos_daa3*(daa3dt3*dt3zyz + daa3dt2*dt2zyz)
267 da3zzx = dacos_daa3*(daa3dt3*dt3zzx + daa3dt2*dt2zzx)
268C
269C STORE DERIVATIVES OF A
270C
271 da(i,1)=da3zxx
272 da(i,2)=da3zyy
273 da(i,3)=da3zzz
274 da(i,4)=da3zxy
275 da(i,5)=da3zyz
276 da(i,6)=da3zzx
277C
278 temp4 = one*third/sqrt(j2(i)*third)
279 temp5 = two*sqrt(j2(i)*third)*third
280 teta = three_teta * third
281C
282C DERIVATIVES OF PRINCIPAL STRAINS
283C
284 dz1zxx =temp4*dt2zxx*cos(-teta) +temp5*sin(-teta)*da3zxx
285 dz2zxx =temp4*dt2zxx*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zxx
286 dz3zxx =temp4*dt2zxx*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zxx
287
288 dz1zyy =temp4*dt2zyy*cos(-teta) +temp5*sin(-teta)*da3zyy
289 dz2zyy =temp4*dt2zyy*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zyy
290 dz3zyy =temp4*dt2zyy*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zyy
291
292 dz1zzz =temp4*dt2zzz*cos(-teta) +temp5*sin(-teta)*da3zzz
293 dz2zzz =temp4*dt2zzz*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zzz
294 dz3zzz =temp4*dt2zzz*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zzz
295
296 dz1zxy =temp4*dt2zxy*cos(-teta) +temp5*sin(-teta)*da3zxy
297 dz2zxy =temp4*dt2zxy*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zxy
298 dz3zxy =temp4*dt2zxy*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zxy
299
300 dz1zyz =temp4*dt2zyz*cos(-teta) +temp5*sin(-teta)*da3zyz
301 dz2zyz =temp4*dt2zyz*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zyz
302 dz3zyz =temp4*dt2zyz*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zyz
303
304 dz1zzx =temp4*dt2zzx*cos(-teta) +temp5*sin(-teta)*da3zzx
305 dz2zzx =temp4*dt2zzx*cos(-teta+two*pi*third) +temp5*sin(-teta+two*pi*third)*da3zzx
306 dz3zzx =temp4*dt2zzx*cos(-teta+four*pi*third) +temp5*sin(-teta+four*pi*third)*da3zzx
307c
308c store first derivatives
309c store zero'd derivatives
310c
311 dzero(i,1) =z1(i)
312 dzero(i,2) =z2(i)
313 dzero(i,3) =z3(i)
314 dfirst(i,1,1)=dz1zxx + third
315 dfirst(i,1,2)=dz1zyy + third
316 dfirst(i,1,3)=dz1zzz + third
317 dfirst(i,1,4)=dz1zxy
318 dfirst(i,1,5)=dz1zyz
319 dfirst(i,1,6)=dz1zzx
320 dfirst(i,2,1)=dz2zxx + third
321 dfirst(i,2,2)=dz2zyy + third
322 dfirst(i,2,3)=dz2zzz + third
323 dfirst(i,2,4)=dz2zxy
324 dfirst(i,2,5)=dz2zyz
325 dfirst(i,2,6)=dz2zzx
326 dfirst(i,3,1)=dz3zxx + third
327 dfirst(i,3,2)=dz3zyy + third
328 dfirst(i,3,3)=dz3zzz + third
329 dfirst(i,3,4)=dz3zxy
330 dfirst(i,3,5)=dz3zyz
331 dfirst(i,3,6)=dz3zzx
332
333C
334C PART 3
335C SECOND DERIVATIVES OF PRINCIPAL STRAINS
336C SECOND DERIVATIVES OF H1 ARE ALL ZERO
337C SECOND DERIVATIVES OF H2
338C
339 dah2zxxzxx = zero
340 dah2zxxzyy = one
341 dah2zxxzzz = one
342 dah2zxxzxy = zero
343 dah2zxxzyz = zero
344 dah2zxxzzx = zero
345C
346 dah2zyyzxx = one
347 dah2zyyzyy = zero
348 dah2zyyzzz = one
349 dah2zyyzxy = zero
350 dah2zyyzyz = zero
351 dah2zyyzzx = zero
352C
353 dah2zzzzxx = one
354 dah2zzzzyy = one
355 dah2zzzzzz = zero
356 dah2zzzzxy = zero
357 dah2zzzzyz = zero
358 dah2zzzzzx = zero
359C
360 dah2zxyzxx = zero
361 dah2zxyzyy = zero
362 dah2zxyzzz = zero
363 dah2zxyzxy = -two
364 dah2zxyzyz = zero
365 dah2zxyzzx = zero
366C
367 dah2zyzzxx = zero
368 dah2zyzzyy = zero
369 dah2zyzzzz = zero
370 dah2zyzzxy = zero
371 dah2zyzzyz = -two
372 dah2zyzzzx = zero
373C
374 dah2zzxzxx = zero
375 dah2zzxzyy = zero
376 dah2zzxzzz = zero
377 dah2zzxzxy = zero
378 dah2zzxzyz = zero
379 dah2zzxzzx = -two
380C
381C SECOND DERIVATIVES OF H3
382C
383 dah3zxxzxx=zero
384 dah3zxxzyy=zzz(i)
385 dah3zxxzzz=zyy(i)
386 dah3zxxzxy=zero
387 dah3zxxzyz=-two*zyz(i)
388 dah3zxxzzx=zero
389C
390 dah3zyyzxx=zzz(i)
391 dah3zyyzyy=zero
392 dah3zyyzzz=zxx(i)
393 dah3zyyzxy=zero
394 dah3zyyzyz=zero
395 dah3zyyzzx=-two*zzx(i)
396C
397 dah3zzzzxx=zyy(i)
398 dah3zzzzyy=zxx(i)
399 dah3zzzzzz=zero
400 dah3zzzzxy=-two*zxy(i)
401 dah3zzzzyz=zero
402 dah3zzzzzx=zero
403C
404 dah3zxyzxx=zero
405 dah3zxyzyy=zero
406 dah3zxyzzz=-two*zxy(i)
407 dah3zxyzxy=-two*zzz(i)
408 dah3zxyzyz=two*zzx(i)
409 dah3zxyzzx=two*zyz(i)
410C
411 dah3zyzzxx=-two*zyz(i)
412 dah3zyzzyy=zero
413 dah3zyzzzz=zero
414 dah3zyzzxy=two*zzx(i)
415 dah3zyzzyz=-two*zxx(i)
416 dah3zyzzzx=two*zxy(i)
417C
418 dah3zzxzxx=zero
419 dah3zzxzyy=-two*zzx(i)
420 dah3zzxzzz=zero
421 dah3zzxzxy=two*zyz(i)
422 dah3zzxzyz=two*zxy(i)
423 dah3zzxzzx=-two*zyy(i)
424C
425C SECOND DERIVATIVES OF T2
426C
427 dat22(i,1,1)=-dah2zxxzxx+(two_third)
428 dat22(i,1,2)=-dah2zxxzyy+(two_third)
429 dat22(i,1,3)=-dah2zxxzzz+(two_third)
430 dat22(i,1,4)=-dah2zxxzxy
431 dat22(i,1,5)=-dah2zxxzyz
432 dat22(i,1,6)=-dah2zxxzzx
433C
434 dat22(i,2,1)=-dah2zyyzxx+(two_third)
435 dat22(i,2,2)=-dah2zyyzyy+(two_third)
436 dat22(i,2,3)=-dah2zyyzzz+(two_third)
437 dat22(i,2,4)=-dah2zyyzxy
438 dat22(i,2,5)=-dah2zyyzyz
439 dat22(i,2,6)=-dah2zyyzzx
440C
441 dat22(i,3,1)=-dah2zzzzxx+(two_third)
442 dat22(i,3,2)=-dah2zzzzyy+(two_third)
443 dat22(i,3,3)=-dah2zzzzzz+(two_third)
444 dat22(i,3,4)=-dah2zzzzxy
445 dat22(i,3,5)=-dah2zzzzyz
446 dat22(i,3,6)=-dah2zzzzzx
447C
448 dat22(i,4,1)=-dah2zxyzxx
449 dat22(i,4,2)=-dah2zxyzyy
450 dat22(i,4,3)=-dah2zxyzzz
451 dat22(i,4,4)=-dah2zxyzxy
452 dat22(i,4,5)=-dah2zxyzyz
453 dat22(i,4,6)=-dah2zxyzzx
454C
455 dat22(i,5,1)=-dah2zyzzxx
456 dat22(i,5,2)=-dah2zyzzyy
457 dat22(i,5,3)=-dah2zyzzzz
458 dat22(i,5,4)=-dah2zyzzxy
459 dat22(i,5,5)=-dah2zyzzyz
460 dat22(i,5,6)=-dah2zyzzzx
461C
462 dat22(i,6,1)=-dah2zzxzxx
463 dat22(i,6,2)=-dah2zzxzyy
464 dat22(i,6,3)=-dah2zzxzzz
465 dat22(i,6,4)=-dah2zzxzxy
466 dat22(i,6,5)=-dah2zzxzyz
467 dat22(i,6,6)=-dah2zzxzzx
468C
469C SECOND DERIVATIVES OF T3
470C
471 dat32(i,1,1) = dah3zxxzxx - (ah1(i)*third)*dah2zxxzxx
472 dat32(i,1,2) = dah3zxxzyy - (ah1(i)*third)*dah2zxxzyy
473 dat32(i,1,3) = dah3zxxzzz - (ah1(i)*third)*dah2zxxzzz
474 dat32(i,1,4) = dah3zxxzxy - (ah1(i)*third)*dah2zxxzxy
475 dat32(i,1,5) = dah3zxxzyz - (ah1(i)*third)*dah2zxxzyz
476 dat32(i,1,6) = dah3zxxzzx - (ah1(i)*third)*dah2zxxzzx
477C
478 dat32(i,2,1) = dah3zyyzxx - (ah1(i)*third)*dah2zyyzxx
479 dat32(i,2,2) = dah3zyyzyy - (ah1(i)*third)*dah2zyyzyy
480 dat32(i,2,3) = dah3zyyzzz - (ah1(i)*third)*dah2zyyzzz
481 dat32(i,2,4) = dah3zyyzxy - (ah1(i)*third)*dah2zyyzxy
482 dat32(i,2,5) = dah3zyyzyz - (ah1(i)*third)*dah2zyyzyz
483 dat32(i,2,6) = dah3zyyzzx - (ah1(i)*third)*dah2zyyzzx
484C
485 dat32(i,3,1) = dah3zzzzxx - (ah1(i)*third)*dah2zzzzxx
486 dat32(i,3,2) = dah3zzzzyy - (ah1(i)*third)*dah2zzzzyy
487 dat32(i,3,3) = dah3zzzzzz - (ah1(i)*third)*dah2zzzzzz
488 dat32(i,3,4) = dah3zzzzxy - (ah1(i)*third)*dah2zzzzxy
489 dat32(i,3,5) = dah3zzzzyz - (ah1(i)*third)*dah2zzzzyz
490 dat32(i,3,6) = dah3zzzzzx - (ah1(i)*third)*dah2zzzzzx
491C
492 dat32(i,4,1) = dah3zxyzxx - (ah1(i)*third)*dah2zxyzxx
493 dat32(i,4,2) = dah3zxyzyy - (ah1(i)*third)*dah2zxyzyy
494 dat32(i,4,3) = dah3zxyzzz - (ah1(i)*third)*dah2zxyzzz
495 dat32(i,4,4) = dah3zxyzxy - (ah1(i)*third)*dah2zxyzxy
496 dat32(i,4,5) = dah3zxyzyz - (ah1(i)*third)*dah2zxyzyz
497 dat32(i,4,6) = dah3zxyzzx - (ah1(i)*third)*dah2zxyzzx
498C
499 dat32(i,5,1) = dah3zyzzxx - (ah1(i)*third)*dah2zyzzxx
500 dat32(i,5,2) = dah3zyzzyy - (ah1(i)*third)*dah2zyzzyy
501 dat32(i,5,3) = dah3zyzzzz - (ah1(i)*third)*dah2zyzzzz
502 dat32(i,5,4) = dah3zyzzxy - (ah1(i)*third)*dah2zyzzxy
503 dat32(i,5,5) = dah3zyzzyz - (ah1(i)*third)*dah2zyzzyz
504 dat32(i,5,6) = dah3zyzzzx - (ah1(i)*third)*dah2zyzzzx
505
506 dat32(i,6,1) = dah3zzxzxx - (ah1(i)*third)*dah2zzxzxx
507 dat32(i,6,2) = dah3zzxzyy - (ah1(i)*third)*dah2zzxzyy
508 dat32(i,6,3) = dah3zzxzzz - (ah1(i)*third)*dah2zzxzzz
509 dat32(i,6,4) = dah3zzxzxy - (ah1(i)*third)*dah2zzxzxy
510 dat32(i,6,5) = dah3zzxzyz - (ah1(i)*third)*dah2zzxzyz
511 dat32(i,6,6) = dah3zzxzzx - (ah1(i)*third)*dah2zzxzzx
512C
513C ADD FIRST DERIVATIVE TERMS
514C
515 termdat32 = (ten+two) * ah1(i)/twenty7
516 dat32(i,1,1)=dat32(i,1,1)+termdat32*dh1zxx*dh1zxx
517 . -(third)*dh1zxx*dh2zxx-(third)*dh1zxx*dh2zxx
518 dat32(i,1,2)=dat32(i,1,2)+termdat32*dh1zyy*dh1zxx
519 . -(third)*dh1zxx*dh2zyy-(third)*dh1zyy*dh2zxx
520 dat32(i,1,3)=dat32(i,1,3)+termdat32*dh1zzz*dh1zxx
521 . -(third)*dh1zxx*dh2zzz-(third)*dh1zzz*dh2zxx
522 dat32(i,1,4)=dat32(i,1,4)+termdat32*dh1zxy*dh1zxx
523 . -(third)*dh1zxx*dh2zxy-(third)*dh1zxy*dh2zxx
524 dat32(i,1,5)=dat32(i,1,5)+termdat32*dh1zyz*dh1zxx
525 . -(third)*dh1zxx*dh2zyz-(third)*dh1zyz*dh2zxx
526 dat32(i,1,6)=dat32(i,1,6)+termdat32*dh1zzx*dh1zxx
527 . -(third)*dh1zxx*dh2zzx-(third)*dh1zzx*dh2zxx
528C
529 dat32(i,2,1)=dat32(i,2,1)+termdat32*dh1zxx*dh1zyy
530 . -(third)*dh1zyy*dh2zxx-(third)*dh1zxx*dh2zyy
531 dat32(i,2,2)=dat32(i,2,2)+termdat32*dh1zyy*dh1zyy
532 . -(third)*dh1zyy*dh2zyy-(third)*dh1zyy*dh2zyy
533 dat32(i,2,3)=dat32(i,2,3)+termdat32*dh1zzz*dh1zyy
534 . -(third)*dh1zyy*dh2zzz-(third)*dh1zzz*dh2zyy
535 dat32(i,2,4)=dat32(i,2,4)+termdat32*dh1zxy*dh1zyy
536 . -(third)*dh1zyy*dh2zxy-(third)*dh1zxy*dh2zyy
537 dat32(i,2,5)=dat32(i,2,5)+termdat32*dh1zyz*dh1zyy
538 . -(third)*dh1zyy*dh2zyz-(third)*dh1zyz*dh2zyy
539 dat32(i,2,6)=dat32(i,2,6)+termdat32*dh1zzx*dh1zyy
540 . -(third)*dh1zyy*dh2zzx-(third)*dh1zzx*dh2zyy
541C
542 dat32(i,3,1)=dat32(i,3,1)+termdat32*dh1zxx*dh1zzz
543 . -(third)*dh1zzz*dh2zxx-(third)*dh1zxx*dh2zzz
544 dat32(i,3,2)=dat32(i,3,2)+termdat32*dh1zyy*dh1zzz
545 . -(third)*dh1zzz*dh2zyy-(third)*dh1zyy*dh2zzz
546 dat32(i,3,3)=dat32(i,3,3)+termdat32*dh1zzz*dh1zzz
547 . -(third)*dh1zzz*dh2zzz-(third)*dh1zzz*dh2zzz
548 dat32(i,3,4)=dat32(i,3,4)+termdat32*dh1zxy*dh1zzz
549 . -(third)*dh1zzz*dh2zxy-(third)*dh1zxy*dh2zzz
550 dat32(i,3,5)=dat32(i,3,5)+termdat32*dh1zyz*dh1zzz
551 . -(third)*dh1zzz*dh2zyz-(third)*dh1zyz*dh2zzz
552 dat32(i,3,6)=dat32(i,3,6)+termdat32*dh1zzx*dh1zzz
553 . -(third)*dh1zzz*dh2zzx-(third)*dh1zzx*dh2zzz
554C
555 dat32(i,4,1)=dat32(i,4,1)+termdat32*dh1zxx*dh1zxy
556 . -(third)*dh1zxy*dh2zxx-(third)*dh1zxx*dh2zxy
557 dat32(i,4,2)=dat32(i,4,2)+termdat32*dh1zyy*dh1zxy
558 . -(third)*dh1zxy*dh2zyy-(third)*dh1zyy*dh2zxy
559 dat32(i,4,3)=dat32(i,4,3)+termdat32*dh1zzz*dh1zxy
560 . -(third)*dh1zxy*dh2zzz-(third)*dh1zzz*dh2zxy
561 dat32(i,4,4)=dat32(i,4,4)+termdat32*dh1zxy*dh1zxy
562 . -(third)*dh1zxy*dh2zxy-(third)*dh1zxy*dh2zxy
563 dat32(i,4,5)=dat32(i,4,5)+termdat32*dh1zyz*dh1zxy
564 . -(third)*dh1zxy*dh2zyz-(third)*dh1zyz*dh2zxy
565 dat32(i,4,6)=dat32(i,4,6)+termdat32*dh1zzx*dh1zxy
566 . -(third)*dh1zxy*dh2zzx-(third)*dh1zzx*dh2zxy
567C
568 dat32(i,5,1)=dat32(i,5,1)+termdat32*dh1zxx*dh1zyz
569 . -(third)*dh1zyz*dh2zxx-(third)*dh1zxx*dh2zyz
570 dat32(i,5,2)=dat32(i,5,2)+termdat32*dh1zyy*dh1zyz
571 . -(third)*dh1zyz*dh2zyy-(third)*dh1zyy*dh2zyz
572 dat32(i,5,3)=dat32(i,5,3)+termdat32*dh1zzz*dh1zyz
573 . -(third)*dh1zyz*dh2zzz-(third)*dh1zzz*dh2zyz
574 dat32(i,5,4)=dat32(i,5,4)+termdat32*dh1zxy*dh1zyz
575 . -(third)*dh1zyz*dh2zxy-(third)*dh1zxy*dh2zyz
576 dat32(i,5,5)=dat32(i,5,5)+termdat32*dh1zyz*dh1zyz
577 . -(third)*dh1zyz*dh2zyz-(third)*dh1zyz*dh2zyz
578 dat32(i,5,6)=dat32(i,5,6)+termdat32*dh1zzx*dh1zyz
579 . -(third)*dh1zyz*dh2zzx-(third)*dh1zzx*dh2zyz
580C
581 dat32(i,6,1)=dat32(i,6,1)+termdat32*dh1zxx*dh1zzx
582 . -(third)*dh1zzx*dh2zxx-(third)*dh1zxx*dh2zzx
583 dat32(i,6,2)=dat32(i,6,2)+termdat32*dh1zyy*dh1zzx
584 . -(third)*dh1zzx*dh2zyy-(third)*dh1zyy*dh2zzx
585 dat32(i,6,3)=dat32(i,6,3)+termdat32*dh1zzz*dh1zzx
586 . -(third)*dh1zzx*dh2zzz-(third)*dh1zzz*dh2zzx
587 dat32(i,6,4)=dat32(i,6,4)+termdat32*dh1zxy*dh1zzx
588 . -(third)*dh1zzx*dh2zxy-(third)*dh1zxy*dh2zzx
589 dat32(i,6,5)=dat32(i,6,5)+termdat32*dh1zyz*dh1zzx
590 . -(third)*dh1zzx*dh2zyz-(third)*dh1zyz*dh2zzx
591 dat32(i,6,6)=dat32(i,6,6)+termdat32*dh1zzx*dh1zzx
592 . -(third)*dh1zzx*dh2zzx-(third)*dh1zzx*dh2zzx
593C
594C SECOND DERIVATIVES OF A
595C
596 d2acos_daa32 = - aa3(i) - (aa3(i)**3)/two - (aa3(i)**5)/eight -
597 . (aa3(i)**7)/48.0d0 - (aa3(i)**9)/384.0d0
598 d2a_dt22 = (j3(i)*j2(i)/36.0d0)*(one/sqrt(((j2(i)**3)/twenty7)**3))*(five/two)
599 d2a_dt3dt2 = -(j2(i)**2)/(36.0d0*sqrt(((j2(i)**3)/twenty7)**3))
600 d2a_dt2dt3 = -((j2(i)**2)/36.0d0)*(one/sqrt(((j2(i)**3)/twenty7)**3))
601 DO j=1,6
602 DO k=1,6
603 da2(i,j,k) = d2acos_daa32*((daa3dt3*dat3(i,k) + daa3dt2*dat2(i,k))*
604 . (daa3dt3*dat3(i,j) + daa3dt2*dat2(i,j))) +
605 . dacos_daa3*(d2a_dt3dt2*dat2(i,j)*dat3(i,k) +
606 . d2a_dt2dt3*dat3(i,j)*dat2(i,k) +
607 . d2a_dt22*dat2(i,j)*dat2(i,k))
608 ENDDO
609 ENDDO
610C
611C SIMILAR ZEROING OUT AS IN THE FIRST DERIVATIVES
612C UNLIKE WITH THE FIRST DERIVATIVES, HERE IT MAKES A DIFFERENCE
613
614C SECOND DERIVATIVES OF PRINCIPAL STRAINS
615C MAKES USE OF FIRST DERIVATIVES OF T2 : DAT2
616C MAKES USE OF FIRST DERIVATIVES OF A : DA
617C MAKES USE OF SECOND DERIVATIVES OF T2 : DAT22
618C MAKES USE OF SECOND DERIVATIVES OF A : DA2
619C
620
621 DO l=1,3
622 ac = cos(-teta +(l-1)*two*pi*third)
623 as = sin(-teta +(l-1)*two*pi*third)
624 aa2(i)= sqrt(j2(i)*third)
625 DO j=1,6
626 DO k=1,6
627 t1=(-one/six/aa2(i)**3*third)* ac * dat2(i,j)*dat2(i,k)
628 t2=( third/aa2(i)) * ac * dat22(i,j,k)
629 t3=(one/nine/aa2(i)) * as * dat2(i,k) * da(i,j)
630 t4=(one/nine/aa2(i)) * as * dat2(i,j) * da(i,k)
631 t5=(-two/nine) *aa2(i) * ac * da(i,j) * da(i,k)
632 t6=(two*third) *aa2(i) * as * da2(i,j,k) ! IF J2==ZERO T6=ZERO
633 dsecond(i,l,j,k) = t1 + t2 + t3 + t4 + t5 + t6
634 ENDDO
635 ENDDO
636 ENDDO
637 ENDDO
638C
639CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
640C END OF CALCULATION OF
641C PRINCIPAL STRAIN DERIVATIVES
642CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
643C
644C
645C
646CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
647C COMPUTE STATIC PRINCIPAL 2PK STRESSES
648C COMPUTE CURRENT TANGENT TO STRESS-STRAIN CURVE
649C IN GENERAL : USE THE LOAD CURVE WITH LOWEST STRAIN RATE
650CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
651C
652C
653C
654 strainrate(1:nel,1)= zero
655 strainrate(1:nel,2)= zero
656 strainrate(1:nel,3)= zero
657C
658 DO i=1,nel
659C
660c VOLUMETRIC STRAIN FROM JACOBIAN
661c POSITIVE IN COMPRESSION, FOR INSTANCE 0.8 MEANS 80% COMPRESSION
662c TO BE USED AS ABCISSA FOR THE TABLE3D
663 volstr(i) = one-jac(i)
664
665c PRINCIPAL STRETCHES FROM PRINCIPAL GL STRAINS
666 alambda(i,1)=sqrt(two*z1(i)+one)
667 alambda(i,2)=sqrt(two*z2(i)+one)
668 alambda(i,3)=sqrt(two*z3(i)+one)
669c
670c PRINCIPAL ENGINEERING STRAIN POSITIVE IN COMPRESSION
671c TO BE USED AS ABCISSA FOR THE LOAD CURVES
672C
673 epsilon(i,1) = one - alambda(i,1)
674 epsilon(i,2) = one - alambda(i,2)
675 epsilon(i,3) = one - alambda(i,3)
676C
677C ENGINEERING STRAIN RATE, COMPUTED FROM TRUE STRAIN RATE
678C AND WE TAKE THE ABSOLUTE VALUE ( SO LINEAR STRAIN RATE DEFINITION )
679C TO BE USED AS ABCISSA FOR THE TABLE3D (SECOND ABCISSA)
680C
681 strainrate(i,1)= abs(drate(i,1))
682 strainrate(i,2)= abs(drate(i,2))
683 strainrate(i,3)= abs(drate(i,3))
684C
685 ENDDO !I=1,NEL
686C
687C 2 STRESS LOOKUPS ARE DONE NOW, STATIC AND DYNAMIC
688C BOTH USE THE TABLE2D THAT CORRESPONDS TO THE CURRENT VALUE OF VOLSTR
689C
690
691C PRINCIPAL ENGINEERING STRESS POSITIVE IN COMPRESSION
692C THIS IS FROM THE QUASISTATIC LOAD CURVE ( LOWEST STRAIN RATE )
693
694
695 !check table dimension
696 IF (ndim_table == 1) THEN
697 xvec1(1:nel,1) = epsilon(1:nel,1)
698 ipos1(1:nel,1) = vartmp(1:nel,1)
699 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig1,dsig1)
700 vartmp(1:nel,1) = ipos1(1:nel,1)
701 fsig1(1:nel) = fsig1(1:nel)*scal
702 dsig1(1:nel) = dsig1(1:nel)*scal
703
704 xvec1(1:nel,1) = epsilon(1:nel,2)
705 ipos1(1:nel,1) = vartmp(1:nel,3)
706 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig2,dsig2)
707 vartmp(1:nel,3) = ipos1(1:nel,1)
708 fsig2(1:nel) = fsig2(1:nel)*scal
709 dsig2(1:nel) = dsig2(1:nel)*scal
710
711 xvec1(1:nel,1) = epsilon(1:nel,3)
712 ipos1(1:nel,1) = vartmp(1:nel,5)
713 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig3,dsig3)
714 vartmp(1:nel,5) = ipos1(1:nel,1)
715 fsig3(1:nel) = fsig3(1:nel)*scal
716 dsig3(1:nel) = dsig3(1:nel)*scal
717c
718 !if 2D table is used
719 ELSEIF (ndim_table == 2) THEN
720 xvec2(1:nel,1) = epsilon(1:nel,1)
721 xvec2(1:nel,2) = zero ! QUASISTATIC
722 ipos2(1:nel,1) = vartmp(1:nel,1)
723 ipos2(1:nel,2) = 1
724 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig1,dsig1)
725 vartmp(1:nel,1) = ipos2(1:nel,1)
726 fsig1(1:nel) = fsig1(1:nel)*scal
727 dsig1(1:nel) = dsig1(1:nel)*scal
728
729 xvec2(1:nel,1) = epsilon(1:nel,2)
730 xvec2(1:nel,2) = zero ! QUASISTATIC
731 ipos2(1:nel,1) = vartmp(1:nel,3)
732 ipos2(1:nel,2) = 1
733 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig2,dsig2)
734 vartmp(1:nel,3) = ipos2(1:nel,1)
735 fsig2(1:nel) = fsig2(1:nel)*scal
736 dsig2(1:nel) = dsig2(1:nel)*scal
737
738 xvec2(1:nel,1) = epsilon(1:nel,3)
739 xvec2(1:nel,2) = zero ! QUASISTATIC
740 ipos2(1:nel,1) = vartmp(1:nel,5)
741 ipos2(1:nel,2) = 1
742 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig3,dsig3)
743 vartmp(1:nel,5) = ipos2(1:nel,1)
744 fsig3(1:nel) = fsig3(1:nel)*scal
745 dsig3(1:nel) = dsig3(1:nel)*scal
746
747 !if 3D table is used
748 ELSEIF (ndim_table == 3) THEN
749
750 xvec(1:nel,1) = epsilon(1:nel,1)
751 xvec(1:nel,2) = zero ! QUASISTATIC
752 xvec(1:nel,3) = volstr(1:nel)
753 ipos(1:nel,1) = vartmp(1:nel,1)
754 ipos(1:nel,2) = 1
755 ipos(1:nel,3) = vartmp(1:nel,2)
756 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig1,dsig1)
757 vartmp(1:nel,1) = ipos(1:nel,1)
758 vartmp(1:nel,2) = ipos(1:nel,3)
759 fsig1(1:nel) = fsig1(1:nel)*scal
760 dsig1(1:nel) = dsig1(1:nel)*scal
761c
762 xvec(1:nel,1) = epsilon(1:nel,2)
763 xvec(1:nel,2) = zero ! QUASISTATIC
764 xvec(1:nel,3) = volstr(1:nel)
765 ipos(1:nel,1) = vartmp(1:nel,3)
766 ipos(1:nel,2) = 1
767 ipos(1:nel,3) = vartmp(1:nel,4)
768 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig2,dsig2)
769 vartmp(1:nel,3) = ipos(1:nel,1)
770 vartmp(1:nel,4) = ipos(1:nel,3)
771 fsig2(1:nel) = fsig2(1:nel)*scal
772 dsig2(1:nel) = dsig2(1:nel)*scal
773c
774 xvec(1:nel,1) = epsilon(1:nel,3)
775 xvec(1:nel,2) = zero ! QUASISTATIC
776 xvec(1:nel,3) = volstr(1:nel)
777 ipos(1:nel,1) = vartmp(1:nel,5)
778 ipos(1:nel,2) = 1
779 ipos(1:nel,3) = vartmp(1:nel,6)
780 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig3,dsig3)
781 vartmp(1:nel,5) = ipos(1:nel,1)
782 vartmp(1:nel,6) = ipos(1:nel,3)
783 fsig3(1:nel) = fsig3(1:nel)*scal
784 dsig3(1:nel) = dsig3(1:nel)*scal
785 ENDIF
786c
787 DO i=1,nel
788C
789 sigma(i,1) = fsig1(i)
790 slope(i,1) = dsig1(i)
791
792 sigma(i,2) = fsig2(i)
793 slope(i,2) = dsig2(i)
794
795 sigma(i,3) = fsig3(i)
796 slope(i,3) = dsig3(i)
797 ENDDO !I=1,NEL
798C
799C PRINCIPAL ENGINEERING STRESS POSITIVE IN COMPRESSION
800C THIS IS FROM THE COMPLETE RATE DEPENDENT TABLE3D
801C SO 'ENGIRATE' HAS TO COME IN AS AN ADDITIONAL ARGUMENT
802 !check table dimension
803 IF (ndim_table == 1) THEN
804C DIRECTION 1
805 xvec1(1:nel,1) = epsilon(1:nel,1)
806 ipos1(1:nel,1) = vartmp(1:nel,7)
807 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig1,dsig1)
808 vartmp(1:nel,7) = ipos1(1:nel,1)
809 fsig1(1:nel) = fsig1(1:nel)*scal
810 dsig1(1:nel) = dsig1(1:nel)*scal
811
812C DIRECTION 2
813
814 xvec1(1:nel,1) = epsilon(1:nel,2)
815 ipos1(1:nel,1) = vartmp(1:nel,10)
816 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig2,dsig2)
817 vartmp(1:nel,10) = ipos1(1:nel,1)
818 fsig2(1:nel) = fsig2(1:nel)*scal
819 dsig2(1:nel) = dsig2(1:nel)*scal
820
821C DIRECTION 3
822
823 xvec1(1:nel,1) = epsilon(1:nel,3)
824 ipos1(1:nel,1) = vartmp(1:nel,13)
825 CALL table_mat_vinterp(func_sig,nel,nel,ipos1,xvec1,fsig3,dsig3)
826 vartmp(1:nel,13) = ipos1(1:nel,1)
827 fsig3(1:nel) = fsig3(1:nel)*scal
828 dsig3(1:nel) = dsig3(1:nel)*scal
829
830 !if 2D table is used
831 ELSEIF (ndim_table == 2) THEN
832C DIRECTION 1
833 xvec2(1:nel,1) = epsilon(1:nel,1)
834 xvec2(1:nel,2) = strainrate(1:nel,1) / xscal
835 ipos2(1:nel,1) = vartmp(1:nel,7)
836 ipos2(1:nel,2) = vartmp(1:nel,8)
837 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig1,dsig1)
838 vartmp(1:nel,7) = ipos2(1:nel,1)
839 vartmp(1:nel,8) = ipos2(1:nel,2)
840 fsig1(1:nel) = fsig1(1:nel)*scal
841 dsig1(1:nel) = dsig1(1:nel)*scal
842C DIRECTION 2
843 xvec2(1:nel,1) = epsilon(1:nel,2)
844 xvec2(1:nel,2) = strainrate(1:nel,2) / xscal
845 ipos2(1:nel,1) = vartmp(1:nel,10)
846 ipos2(1:nel,2) = vartmp(1:nel,11)
847 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig2,dsig2)
848 vartmp(1:nel,10) = ipos2(1:nel,1)
849 vartmp(1:nel,11) = ipos2(1:nel,2)
850 fsig2(1:nel) = fsig2(1:nel)*scal
851 dsig2(1:nel) = dsig2(1:nel)*scal
852C DIRECTION 3
853 xvec2(1:nel,1) = epsilon(1:nel,3)
854 xvec2(1:nel,2) = strainrate(1:nel,3) / xscal
855 ipos2(1:nel,1) = vartmp(1:nel,13)
856 ipos2(1:nel,2) = vartmp(1:nel,14)
857 CALL table_mat_vinterp(func_sig,nel,nel,ipos2,xvec2,fsig3,dsig3)
858 vartmp(1:nel,13) = ipos2(1:nel,1)
859 vartmp(1:nel,14) = ipos2(1:nel,2)
860 fsig3(1:nel) = fsig3(1:nel)*scal
861 dsig3(1:nel) = dsig3(1:nel)*scal
862c
863 !if 3D table is used
864 ELSEIF (ndim_table == 3) THEN
865C DIRECTION 1
866 xvec(1:nel,1) = epsilon(1:nel,1)
867 xvec(1:nel,2) = strainrate(1:nel,1) / xscal
868 xvec(1:nel,3) = volstr(1:nel)
869!
870 ipos(1:nel,1) = vartmp(1:nel,7)
871 ipos(1:nel,2) = vartmp(1:nel,8)
872 ipos(1:nel,3) = vartmp(1:nel,9)
873 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig1,dsig1)
874 vartmp(1:nel,7) = ipos(1:nel,1)
875 vartmp(1:nel,8) = ipos(1:nel,2)
876 vartmp(1:nel,9) = ipos(1:nel,3)
877 fsig1(1:nel) = fsig1(1:nel)*scal
878 dsig1(1:nel) = dsig1(1:nel)*scal
879
880C DIRECTION 2
881
882 xvec(1:nel,1) = epsilon(1:nel,2)
883 xvec(1:nel,2) = strainrate(1:nel,2) / xscal
884 xvec(1:nel,3) = volstr(1:nel)
885 ipos(1:nel,1) = vartmp(1:nel,10)
886 ipos(1:nel,2) = vartmp(1:nel,11)
887 ipos(1:nel,3) = vartmp(1:nel,12)
888 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig2,dsig2)
889 vartmp(1:nel,10) = ipos(1:nel,1)
890 vartmp(1:nel,11) = ipos(1:nel,2)
891 vartmp(1:nel,12) = ipos(1:nel,3)
892 fsig2(1:nel) = fsig2(1:nel)*scal
893 dsig2(1:nel) = dsig2(1:nel)*scal
894
895C DIRECTION 3
896
897 xvec(1:nel,1) = epsilon(1:nel,3)
898 xvec(1:nel,2) = strainrate(1:nel,3) / xscal
899 xvec(1:nel,3) = volstr(1:nel)
900 ipos(1:nel,1) = vartmp(1:nel,13)
901 ipos(1:nel,2) = vartmp(1:nel,14)
902 ipos(1:nel,3) = vartmp(1:nel,15)
903 CALL table_mat_vinterp(func_sig,nel,nel,ipos,xvec,fsig3,dsig3)
904 vartmp(1:nel,13) = ipos(1:nel,1)
905 vartmp(1:nel,14) = ipos(1:nel,2)
906 vartmp(1:nel,15) = ipos(1:nel,3)
907 fsig3(1:nel) = fsig3(1:nel)*scal
908 dsig3(1:nel) = dsig3(1:nel)*scal
909
910 ENDIF !NDIM_TABLE
911
912
913
914 DO i=1,nel
915
916 dynsigma(i,1) = fsig1(i)
917 dynsigma(i,2) = fsig2(i)
918 dynsigma(i,3) = fsig3(i)
919 yld(i) = sqrt(dynsigma(i,1)**2 + dynsigma(i,2)**2 + dynsigma(i,3)**2 )
920C
921C PRINCIPAL 2PK STRESS POSITIVE IN TENSION
922C
923 sik(i,1)=-sigma(i,1)/max(alambda(i,1),em20)
924 sik(i,2)=-sigma(i,2)/max(alambda(i,2),em20)
925 sik(i,3)=-sigma(i,3)/max(alambda(i,3),em20)
926C
927C INCREMENTAL PRINCIPAL STIFFNESSES
928C
929 dik(i,1) = slope(i,1)/max(alambda(i,1)**2,em20)+sigma(i,1)/max(alambda(i,1)**3,em20)
930 dik(i,2) = slope(i,2)/max(alambda(i,2)**2,em20)+sigma(i,2)/max(alambda(i,2)**3,em20)
931 dik(i,3) = slope(i,3)/max(alambda(i,3)**2,em20)+sigma(i,3)/max(alambda(i,3)**3,em20)
932
933C
934C VISCOUS OVERSTRESS IN TERMS OF 2PK
935C POSITIVE IN TENSION
936C
937 drate(i,1)=-(dynsigma(i,1)-sigma(i,1))/max(alambda(i,1),em20)
938 drate(i,2)=-(dynsigma(i,2)-sigma(i,2))/max(alambda(i,2),em20)
939 drate(i,3)=-(dynsigma(i,3)-sigma(i,3))/max(alambda(i,3),em20)
940 ENDDO
941C
942CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
943C COMPUTE INCREMENTAL STIFFNESS MATRIX
944CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
945 DO l=1,3
946 DO k=1,6
947 DO j=1,6
948 DO i=1,nel
949 cijkl(i,j,k)=zero
950 ENDDO
951 ENDDO
952 ENDDO
953 ENDDO
954C
955 DO k=1,6
956 DO j=1,6
957 DO i=1,nel
958 DO l=1,3
959 cijkl(i,j,k)=cijkl(i,j,k)
960 . + dfirst(i,l,j)*dfirst(i,l,k)*dik(i,l)
961 . + dsecond(i,l,j,k)* sik(i,l)
962 ENDDO
963 ENDDO
964 ENDDO
965 ENDDO
966C
967C FACTOR 1/4 FOR VOIGT NOTATION WITH GAMMAS ( 2 GAMMAS)
968C FACTOR 1/2 FOR VOIGT NOTATION WITH GAMMAS ( 1 GAMMA)
969C
970 DO k=4,6
971 DO j=4,6
972 DO i=1,nel
973 cijkl(i,j,k)=cijkl(i,j,k)/four
974 ENDDO
975 ENDDO
976 ENDDO
977C
978 DO k=1,3
979 DO j=4,6
980 DO i=1,nel
981 cijkl(i,j,k)=cijkl(i,j,k)*half
982 ENDDO
983 ENDDO
984 ENDDO
985 DO k=4,6
986 DO j=1,3
987 DO i=1,nel
988 cijkl(i,j,k)=cijkl(i,j,k)*half
989 ENDDO
990 ENDDO
991 ENDDO
992C-------------------
993 RETURN
#define my_real
Definition cppsort.cpp:32
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine table_mat_vinterp(table, dimx, nel, ipos, xx, yy, dydx)