54 s11,s22,s33,s12,s23,s13,numnp,nstart,nend,numcstet,numat_vol, &
79 REAL*8,
DIMENSION(1:3,1:numnp) :: coor
81 REAL*8,
DIMENSION(1:9,1:numat_vol) :: ci
83 REAL*8,
DIMENSION(1:3*numnp) :: r_in
85 REAL*8,
DIMENSION(1:3*numnp) ::
d
87 REAL*8,
DIMENSION(1:numcstet) :: s11, s22, s33, s12, s23, s13
89 REAL*8 :: u1,u2,u3,u4,v1,v2,v3,v4,w1,w2,w3,w4
91 REAL*8 :: x1,x2,x3,x4,y1,y2,y3,y4,z1,z2,z3,z4
93 REAL*8 :: x14, x24, x34, y14, y24, y34, z14, z24, z34
95 REAL*8 :: c11, c21, c31
97 REAL*8 :: vx6, vx6inv, vol
99 REAL*8 :: b1,b2,b3,b4,b5,b6,b7,b8,b9,b10,b11,b12
101 REAL*8 :: dudx,dvdy,dwdz,dudy,dvdx,dvdz,dwdy,dudz,dwdx
103 REAL*8 :: e11,e22,e33,e12,e23,e13
105 INTEGER,
DIMENSION(1:4,1:numcstet) :: lmcstet
107 INTEGER,
DIMENSION(1:numcstet) :: matcstet
109 INTEGER :: n1,n2,n3,n4
111 INTEGER ::
i,
j,nstart,nend
112 INTEGER :: k1n1,k1n2,k1n3,k1n4,k2n1,k2n2,k2n3,k2n4
113 INTEGER :: k3n1,k3n2,k3n3,k3n4
115 REAL*8 :: f11, f12, f13, f21, f22, f23, f31, f32, f33
116 REAL*8 :: j1,j2,j2inv
117 REAL*8 :: c11, c12, c13, c21, c22, c23, c31, c32, c33
118 REAL*8 :: xmu, xlambda
180 c11 = y24*z34 - z24*y34
181 c21 = -( x24*z34 - z24*x34 )
182 c31 = x24*y34 - y24*x34
184 vx6 = -( x14*c11 + y14*c21 + z14*c31 )
198 b1 = -(y3*z4 - y4*z3 - y2*z4 + y2*z3 + z2*y4 - z2*y3) * vx6inv
199 b2 = (x3*z4 - x4*z3 - x2*z4 + x2*z3 + z2*x4 - z2*x3) * vx6inv
200 b3 = -(x3*y4 - x4*y3 - x2*y4 + x2*y3 + y2*x4 - y2*x3) * vx6inv
201 b4 = (y3*z4 - y4*z3 - y1*z4 + y1*z3 + z1*y4 - z1*y3) * vx6inv
202 b5 = -(x3*z4 - x4*z3 - x1*z4 + x1*z3 + z1*x4 - z1*x3) * vx6inv
203 b6 = (x3*y4 - x4*y3 - x1*y4 + x1*y3 + y1*x4 - y1*x3) * vx6inv
204 b7 = -(y2*z4 - z2*y4 - y1*z4 + y1*z2 + z1*y4 - z1*y2) * vx6inv
205 b8 = (x2*z4 - z2*x4 - x1*z4 + x1*z2 + z1*x4 - z1*x2) * vx6inv
206 b9 = -(x2*y4 - y2*x4 - x1*y4 + x1*y2 + y1*x4 - y1*x2) * vx6inv
207 b10 = (y2*z3 - z2*y3 - y1*z3 + y1*z2 + z1*y3 - z1*y2) * vx6inv
208 b11 = -(x2*z3 - z2*x3 - x1*z3 + x1*z2 + z1*x3 - z1*x2) * vx6inv
209 b12 = (x2*y3 - y2*x3 - x1*y3 + x1*y2 + y1*x3 - y1*x2) * vx6inv
213 dudx = b1*u1 + b4*u2 + b7*u3 + b10*u4
214 dvdy = b2*v1 + b5*v2 + b8*v3 + b11*v4
215 dwdz = b3*w1 + b6*w2 + b9*w3 + b12*w4
216 dudy = b2*u1 + b5*u2 + b8*u3 + b11*u4
217 dvdx = b1*v1 + b4*v2 + b7*v3 + b10*v4
218 dvdz = b3*v1 + b6*v2 + b9*v3 + b12*v4
219 dwdy = b2*w1 + b5*w2 + b8*w3 + b11*w4
220 dudz = b3*u1 + b6*u2 + b9*u3 + b12*u4
221 dwdx = b1*w1 + b4*w2 + b7*w3 + b10*w4
235 j1 = f11*(f22*f33-f23*f32)+f12*(f31*f23-f21*f33)+f13*(-f31*f22+f21*f32)
237 WRITE(*,*)
'area has become zero for element',
i
249 c11 = (f11*f11+f21*f21+f31*f31)*j2inv
250 c12 = (f11*f12+f21*f22+f31*f32)*j2inv
251 c13 = (f11*f13+f21*f23+f31*f33)*j2inv
253 c22 = (f12*f12+f22*f22+f32*f32)*j2inv
254 c23 = (f12*f13+f22*f23+f32*f33)*j2inv
257 c33 = (f13*f13+f23*f23+f33*f33)*j2inv
263 s11(
i) = xmu*(1.d0-(c22*c33-c23*c32)) + xlambda*log(j1)*(c22*c33-c23*c32)
264 s22(
i) = xmu*(1.d0-(c11*c33-c31*c13)) + xlambda*log(j1)*(c11*c33-c31*c13)
265 s33(
i) = xmu*(1.d0-(c11*c22-c12*c21)) + xlambda*log(j1)*(c11*c22-c12*c21)
266 s12(
i) = (-c12*c33+c13*c32)*(xmu - xlambda*log(j1))
267 s23(
i) = (-c11*c23+c13*c21)*(xmu - xlambda*log(j1))
268 s13(
i) = (c12*c23-c13*c22)*(xmu - xlambda*log(j1))
278 r_in(k1n1) = r_in(k1n1) - vol* &
279 ( s11(
i)*b1*(1.d0+dudx) + s22(
i)*b2*dudy + s33(
i)*b3*dudz &
280 + s12(
i)*( b2*(1.d0+dudx) + b1*dudy ) &
281 + s23(
i)*( b3*dudy + b2*dudz ) &
282 + s13(
i)*( b3*(1.d0+dudx) + b1*dudz ) )
283 r_in(k2n1) = r_in(k2n1) - vol* &
284 ( s11(
i)*b1*dvdx + s22(
i)*b2*(1.d0+dvdy) + s33(
i)*b3*dvdz &
285 + s12(
i)*( b1*(1.d0+dvdy) + b2*dvdx ) &
286 + s23(
i)*( b3*(1.d0+dvdy) + b2*dvdz ) &
287 + s13(
i)*( b3*dvdx + b1*dvdz ) )
288 r_in(k3n1) = r_in(k3n1) - vol* &
289 ( s11(
i)*b1*dwdx + s22(
i)*b2*dwdy + s33(
i)*b3*(1.d0+dwdz) &
290 + s12(
i)*( b2*dwdx + b1*dwdy ) &
291 + s23(
i)*( b3*dwdy + b2*(1.d0 + dwdz) ) &
292 + s13(
i)*( b3*dwdx + b1*(1.d0 + dwdz) ) )
294 r_in(k1n2) = r_in(k1n2) - vol* &
295 ( s11(
i)*b4*(1.d0+dudx) + s22(
i)*b5*dudy + s33(
i)*b6*dudz &
296 + s12(
i)*( b5*(1.d0+dudx) + b4*dudy ) &
297 + s23(
i)*( b6*dudy + b5*dudz ) &
298 + s13(
i)*( b6*(1.d0+dudx) + b4*dudz ) )
299 r_in(k2n2) = r_in(k2n2) - vol* &
300 ( s11(
i)*b4*dvdx + s22(
i)*b5*(1.d0+dvdy) + s33(
i)*b6*dvdz &
301 + s12(
i)*( b4*(1.d0+dvdy) + b5*dvdx ) &
302 + s23(
i)*( b6*(1.d0+dvdy) + b5*dvdz ) &
303 + s13(
i)*( b6*dvdx + b4*dvdz ) )
304 r_in(k3n2) = r_in(k3n2) - vol* &
305 ( s11(
i)*b4*dwdx + s22(
i)*b5*dwdy + s33(
i)*b6*(1.d0+dwdz) &
306 + s12(
i)*( b5*dwdx + b4*dwdy ) &
307 + s23(
i)*( b6*dwdy + b5*(1.d0 + dwdz) ) &
308 + s13(
i)*( b6*dwdx + b4*(1.d0 + dwdz) ) )
310 r_in(k1n3) = r_in(k1n3) - vol* &
311 ( s11(
i)*b7*(1.d0+dudx) + s22(
i)*b8*dudy + s33(
i)*b9*dudz &
312 + s12(
i)*( b8*(1.d0+dudx) + b7*dudy ) &
313 + s23(
i)*( b9*dudy + b8*dudz ) &
314 + s13(
i)*( b9*(1.d0+dudx) + b7*dudz ) )
315 r_in(k2n3) = r_in(k2n3) - vol* &
316 ( s11(
i)*b7*dvdx + s22(
i)*b8*(1.d0+dvdy) + s33(
i)*b9*dvdz &
317 + s12(
i)*( b7*(1.d0+dvdy) + b8*dvdx ) &
318 + s23(
i)*( b9*(1.d0+dvdy) + b8*dvdz ) &
319 + s13(
i)*( b9*dvdx + b7*dvdz ) )
320 r_in(k3n3) = r_in(k3n3) - vol* &
321 ( s11(
i)*b7*dwdx + s22(
i)*b8*dwdy + s33(
i)*b9*(1.d0+dwdz) &
322 + s12(
i)*( b8*dwdx + b7*dwdy ) &
323 + s23(
i)*( b9*dwdy + b8*(1.d0 + dwdz) ) &
324 + s13(
i)*( b9*dwdx + b7*(1.d0 + dwdz) ) )
326 r_in(k1n4) = r_in(k1n4) - vol* &
327 ( s11(
i)*b10*(1.d0+dudx) + s22(
i)*b11*dudy+s33(
i)*b12*dudz &
328 + s12(
i)*( b11*(1.d0+dudx) + b10*dudy ) &
329 + s23(
i)*( b12*dudy + b11*dudz ) &
330 + s13(
i)*( b12*(1.d0+dudx) + b10*dudz ) )
331 r_in(k2n4) = r_in(k2n4) - vol* &
332 ( s11(
i)*b10*dvdx + s22(
i)*b11*(1.d0+dvdy)+s33(
i)*b12*dvdz &
333 + s12(
i)*( b10*(1.d0+dvdy) + b11*dvdx ) &
334 + s23(
i)*( b12*(1.d0+dvdy) + b11*dvdz ) &
335 + s13(
i)*( b12*dvdx + b10*dvdz ) )
336 r_in(k3n4) = r_in(k3n4) - vol* &
337 ( s11(
i)*b10*dwdx + s22(
i)*b11*dwdy+s33(
i)*b12*(1.d0+dwdz) &
338 + s12(
i)*( b11*dwdx + b10*dwdy ) &
339 + s23(
i)*( b12*dwdy + b11*(1.d0 + dwdz) ) &
340 + s13(
i)*( b12*dwdx + b10*(1.d0 + dwdz) ) )
subroutine v3d4_neohookeancompress(coor, matcstet, lmcstet, R_in, d, ci, S11, S22, S33, S12, S23, S13, numnp, nstart, nend, numcstet, numat_vol, xmu, xlamda)