Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
v3d4n_nl_arrudaboyce.f90
Go to the documentation of this file.
1 !*********************************************************************
2 !* Illinois Open Source License *
3 !* *
4 !* University of Illinois/NCSA *
5 !* Open Source License *
6 !* *
7 !* Copyright@2008, University of Illinois. All rights reserved. *
8 !* *
9 !* Developed by: *
10 !* *
11 !* Center for Simulation of Advanced Rockets *
12 !* *
13 !* University of Illinois *
14 !* *
15 !* www.csar.uiuc.edu *
16 !* *
17 !* Permission is hereby granted, free of charge, to any person *
18 !* obtaining a copy of this software and associated documentation *
19 !* files (the "Software"), to deal with the Software without *
20 !* restriction, including without limitation the rights to use, *
21 !* copy, modify, merge, publish, distribute, sublicense, and/or *
22 !* sell copies of the Software, and to permit persons to whom the *
23 !* Software is furnished to do so, subject to the following *
24 !* conditions: *
25 !* *
26 !* *
27 !* @ Redistributions of source code must retain the above copyright *
28 !* notice, this list of conditions and the following disclaimers. *
29 !* *
30 !* @ Redistributions in binary form must reproduce the above *
31 !* copyright notice, this list of conditions and the following *
32 !* disclaimers in the documentation and/or other materials *
33 !* provided with the distribution. *
34 !* *
35 !* @ Neither the names of the Center for Simulation of Advanced *
36 !* Rockets, the University of Illinois, nor the names of its *
37 !* contributors may be used to endorse or promote products derived *
38 !* from this Software without specific prior written permission. *
39 !* *
40 !* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, *
41 !* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES *
42 !* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND *
43 !* NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR *
44 !* COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
45 !* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
46 !* ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
47 !* USE OR OTHER DEALINGS WITH THE SOFTWARE. *
48 !*********************************************************************
49 !* Please acknowledge The University of Illinois Center for *
50 !* Simulation of Advanced Rockets in works and publications *
51 !* resulting from this software or its derivatives. *
52 !*********************************************************************
53 SUBROUTINE v3d4n_nl_arruda_boyce(numnp, numel, coor, disp, nodes, Rnet, &
54  s11ev, s22ev, s33ev, s12ev, s23ev, s13ev, &
55  numelneigh,elconn,alpha, ahat,&
56  xmu,xkappa,nummatvol, &
57  nprocs,totnumndcomm,totnumneighprocs,neighproclist,numndcomm,neigh_lst,&
58  mpi_status_size,mpi_comm_rocfrac,mpi_double_precision, &
59  reqrcv,reqsnd,statrcv, statsnd)
60 
61  USE rocstar_rocfraccomm
62 
63  IMPLICIT NONE
64 
65  INTEGER :: k1, j1, k2
66 
67  INTEGER :: totnumneighprocs, nprocs,totnumndcomm
68  INTEGER, DIMENSION(1:TotNumNeighProcs) :: neighproclist
69  INTEGER, DIMENSION(1:TotNumNeighProcs) ::numndcomm
70 
71  TYPE(rcv_buf), pointer, DIMENSION(:) :: recvdatafrm
72  TYPE(send_buf), DIMENSION(1:TotNumNeighProcs) :: neigh_lst
73 
74  integer :: nummatvol
75  INTEGER :: numnp, numel
76  REAL*8, DIMENSION(1:numnp) :: ahat
77  REAL*8, DIMENSION(1:3*numnp) :: disp, rnet
78  REAL*8, DIMENSION(1:3,1:numnp) :: coor
79  integer, DIMENSION(1:numnp) ::numelneigh
80  INTEGER, DIMENSION(1:numnp,1:40) :: elconn ! fix 40 should not be hard coded
81  INTEGER,DIMENSION(1:4,1:numel) :: nodes
82  real*8, dimension(1:4,1:numel) :: alpha
83  real*8, dimension(1:NumMatVol) :: xmu, xkappa
84  REAL*8, DIMENSION(1:numel) :: s11ev, s22ev, s33ev, s12ev, s23ev, s13ev
85 
86  INTEGER :: j,i,k
87  REAL*8 :: aaa
88 !-- coordinate holding variable
89  REAL*8 :: x1,x2,x3,x4,y1,y2,y3,y4,z1,z2,z3,z4
90  REAL*8 :: x1d,x2d,x3d,x4d,y1d,y2d,y3d,y4d,z1d,z2d,z3d,z4d
91  INTEGER :: k1n1,k1n2,k1n3,k1n4,k2n1,k2n2,k2n3,k2n4
92  INTEGER :: k3n1,k3n2,k3n3,k3n4
93 
94  REAL*8 :: ahatinv
95  INTEGER :: ielnum
96  REAL*8 :: sixinv
97 !-- x, y, and z displacements of nodes
98  REAL*8 :: u1,u2,u3,u4,v1,v2,v3,v4,w1,w2,w3,w4
99 !-- partial derivatives of the displacement
100  REAL*8 :: dudx,dvdy,dwdz,dudy,dvdx,dvdz,dwdy,dudz,dwdx
101 
102 !-- Coordinate subtractions
103  REAL*8 :: x14, x24, x34, y14, y24, y34, z14, z24, z34
104 !-- Added these to speed up B calculation
105 
106  REAL*8 :: z12, z13,x12,x13,y12,y13
107  REAL*8 :: c11, c12, c13, c21, c22, c23, c31, c32, c33
108  REAL*8 :: sh1, sh2, sh3, sh4, sh5, sh6, sh7, sh8, sh9, sh10, sh11, sh12
109  REAL*8 :: sh13, sh14, sh15
110 ! -- 6*volume and the volume
111  REAL*8 :: vx6
112 !-- spacial derivatives
113  REAL*8 :: b1,b2,b3,b4,b5,b6,b7,b8,b9,b10,b11,b12
114 !-- strains
115  REAL*8 :: e11,e22,e33,e12,e23,e13
116  REAL*8 :: b13, b14, b15
117  REAL*8 :: sums11n,sums22n,sums33n,sums12n,sums23n,sums13n
118  REAL*8 :: sume11n,sume22n,sume33n,sume12n,sume23n,sume13n
119  REAL*8, DIMENSION(1:numnp) :: s11n, s22n, s33n,s12n, s23n, s13n
120 
121  REAL*8 :: e11n, e22n, e33n, e12n, e23n, e13n
122  INTEGER :: k3i,k2i,k1i
123  INTEGER :: ix
124  INTEGER :: n1, n2, n3, n4
125  REAL*8 :: volel, volel0,vx6inv,vol
126  REAL*8 :: s11e, s22e, s33e, s12e, s23e, s13e
127  REAL*8 :: vainv
128 
129  REAL*8 :: coeff1, coeff2
130  REAL*8 :: jacob
131  REAL*8 :: f11, f12, f13, f21, f22, f23, f31, f32, f33
132  REAL*8,DIMENSION(1:numnp) :: f11v, f12v, f13v, f21v, f22v, f23v, f31v, f32v, f33v, voldefv
133  REAL*8 :: vx6def, vx6invdef, voldef
134 
135 
136  INTEGER :: myid
137  INTEGER :: ierr
138 
139  REAL*8, pointer, dimension(:) :: buf
140 
141  integer :: mpi_status_size,mpi_comm_rocfrac,mpi_double_precision
142 
143 !-- Non-block receive, Non-block send request arrays
144  INTEGER, DIMENSION(1:TotNumNeighProcs) :: reqrcv, reqsnd
145  INTEGER, DIMENSION(1:MPI_STATUS_SIZE,1:TotNumNeighProcs) :: statrcv, statsnd
146 
147  sixinv = 1.d0/6.d0
148 
149  DO i = 1, numnp ! for each node
150 
151 ! /* Undeformed Volume of node */
152 
153  vainv = 1.d0/ahat(i)
154 
155  ! /* Initialize to 0, Va, Va0, Fa */
156 
157  f11 = 0.d0
158  f22 = 0.d0
159  f33 = 0.d0
160  f12 = 0.d0
161  f13 = 0.d0
162  f21 = 0.d0
163  f23 = 0.d0
164  f31 = 0.d0
165  f32 = 0.d0
166 
167  voldef = 0.d0
168 
169  DO j = 1, numelneigh(i) ! for each element associated with node i
170 
171  ielnum = elconn(i,j)
172 
173  n1 = nodes(1,ielnum)
174  n2 = nodes(2,ielnum)
175  n3 = nodes(3,ielnum)
176  n4 = nodes(4,ielnum)
177 
178  IF(i.EQ.n1) THEN
179  ix = 1
180  ELSE IF(i.EQ.n2)THEN
181  ix = 2
182  ELSE IF(i.EQ.n3)THEN
183  ix = 3
184  ELSE IF(i.EQ.n4)THEN
185  ix = 4
186  ENDIF
187 
188  k3n1 = 3*n1
189  k3n2 = 3*n2
190  k3n3 = 3*n3
191  k3n4 = 3*n4
192 
193  k2n1 = k3n1 - 1
194  k2n2 = k3n2 - 1
195  k2n3 = k3n3 - 1
196  k2n4 = k3n4 - 1
197 
198  k1n1 = k3n1 - 2
199  k1n2 = k3n2 - 2
200  k1n3 = k3n3 - 2
201  k1n4 = k3n4 - 2
202 
203  ! k#n# dummy variables replaces:
204  u1 = disp(k1n1) ! (3*n1-2)
205  u2 = disp(k1n2) ! (3*n2-2)
206  u3 = disp(k1n3) ! (3*n3-2)
207  u4 = disp(k1n4) ! (3*n4-2)
208  v1 = disp(k2n1) ! (3*n1-1)
209  v2 = disp(k2n2) ! (3*n2-1)
210  v3 = disp(k2n3) ! (3*n3-1)
211  v4 = disp(k2n4) ! (3*n4-1)
212  w1 = disp(k3n1) ! (3*n1)
213  w2 = disp(k3n2) ! (3*n2)
214  w3 = disp(k3n3) ! (3*n3)
215  w4 = disp(k3n4) ! (3*n4)
216 
217  x1 = coor(1,n1)
218  x2 = coor(1,n2)
219  x3 = coor(1,n3)
220  x4 = coor(1,n4)
221  y1 = coor(2,n1)
222  y2 = coor(2,n2)
223  y3 = coor(2,n3)
224  y4 = coor(2,n4)
225  z1 = coor(3,n1)
226  z2 = coor(3,n2)
227  z3 = coor(3,n3)
228  z4 = coor(3,n4)
229 
230  x12 = x1 - x2 ! not used in vol. calc
231  x13 = x1 - x3 ! not used in vol. calc
232  x14 = x1 - x4
233  x24 = x2 - x4
234  x34 = x3 - x4
235  y12 = y1 - y2 ! not used in vol. calc
236  y13 = y1 - y3 ! not used in vol. calc
237  y14 = y1 - y4
238  y24 = y2 - y4
239  y34 = y3 - y4
240  z12 = z1 - z2 ! not used in vol. calc
241  z13 = z1 - z3 ! not used in vol. calc
242  z14 = z1 - z4
243  z24 = z2 - z4
244  z34 = z3 - z4
245 
246  c11 = y24*z34 - z24*y34
247  c21 = -( x24*z34 - z24*x34 )
248  c31 = x24*y34 - y24*x34
249 
250  vx6 = -( x14*c11 + y14*c21 + z14*c31 )
251 
252  vx6inv = 1.d0/vx6
253 
254  volel0 = vx6/6.d0 ! undeformed volume of element (Ve_O)
255 
256 ! Compute the Shape functions
257 ! NOTE: Factored for a more equivalent/compact form then maple's
258 
259  b1 = (y34*z24 - y24*z34) * vx6inv
260  b2 = (z34*x24 - z24*x34) * vx6inv
261  b3 = (x34*y24 - x24*y34) * vx6inv
262  b4 = (y13*z14 - y14*z13) * vx6inv
263  b5 = (z13*x14 - z14*x13) * vx6inv
264  b6 = (x13*y14 - x14*y13) * vx6inv
265  b7 = (y14*z12 - y12*z14) * vx6inv
266  b8 = (z14*x12 - z12*x14) * vx6inv
267  b9 = (x14*y12 - x12*y14) * vx6inv
268  b10 = (y12*z13 - y13*z12) * vx6inv
269  b11 = (z12*x13 - z13*x12) * vx6inv
270  b12 = (x12*y13 - x13*y12) * vx6inv
271 !
272 ! deformation gradients F
273 !
274  f11 = f11 + alpha(ix,ielnum)*volel0*(1.d0 + ( b1*u1 + b4*u2 + b7*u3 + b10*u4 )) ! 1 + ( dudx )
275  f22 = f22 + alpha(ix,ielnum)*volel0*(1.d0 + ( b2*v1 + b5*v2 + b8*v3 + b11*v4 )) ! 1 + ( dvdy )
276  f33 = f33 + alpha(ix,ielnum)*volel0*(1.d0 + ( b3*w1 + b6*w2 + b9*w3 + b12*w4 )) ! 1 + ( dwdz )
277  f12 = f12 + alpha(ix,ielnum)*volel0*(b2*u1 + b5*u2 + b8*u3 + b11*u4) ! dudy
278  f21 = f21 + alpha(ix,ielnum)*volel0*(b1*v1 + b4*v2 + b7*v3 + b10*v4) ! dvdx
279  f23 = f23 + alpha(ix,ielnum)*volel0*(b3*v1 + b6*v2 + b9*v3 + b12*v4) ! dvdz
280  f32 = f32 + alpha(ix,ielnum)*volel0*(b2*w1 + b5*w2 + b8*w3 + b11*w4) ! dwdy
281  f13 = f13 + alpha(ix,ielnum)*volel0*(b3*u1 + b6*u2 + b9*u3 + b12*u4) ! dudz
282  f31 = f31 + alpha(ix,ielnum)*volel0*(b1*w1 + b4*w2 + b7*w3 + b10*w4) ! dwdx
283 
284  ENDDO
285 
286  f11v(i) = f11*vainv
287  f22v(i) = f22*vainv
288  f33v(i) = f33*vainv
289  f12v(i) = f12*vainv
290  f13v(i) = f13*vainv
291  f21v(i) = f21*vainv
292  f23v(i) = f23*vainv
293  f31v(i) = f31*vainv
294  f32v(i) = f32*vainv
295 
296  voldefv(i) = voldef
297 
298  ENDDO
299 
300 
301 
302  ALLOCATE(recvdatafrm(0:nprocs-1))
303 !
304 !----- FORM THE BUFFER CONTAINING COMMUNICATED STRESS MATRIX NODAL VALUES
305 
306  ALLOCATE(buf(1:totnumndcomm/3*10)) ! is really (TotNumNdComm*3 * 3 + 1)
307  k1 = 1
308  DO j1 = 1, totnumneighprocs
309  k = neighproclist(j1)
310  ALLOCATE(recvdatafrm(k)%rcvbuf(1:numndcomm(j1)*10))
311  DO j = 1, numndcomm(j1)
312  k2 = neigh_lst(j1)%NdId(j) !NdCommList(k)%NdId(j)
313  buf(k1 ) = f11v(k2)
314  buf(k1+1) = f12v(k2)
315  buf(k1+2) = f13v(k2)
316  buf(k1+3) = f21v(k2)
317  buf(k1+4) = f22v(k2)
318  buf(k1+5) = f23v(k2)
319  buf(k1+6) = f31v(k2)
320  buf(k1+7) = f32v(k2)
321  buf(k1+8) = f33v(k2)
322  buf(k1+9) = voldefv(k2)
323  k1 = k1 + 10
324  ENDDO
325  ENDDO
326 !
327 !-MPI- RECEIVE THE RECIPRICAL MASS MATRIX DIAGONAL FROM THE NEIGHBORS
328 !
329  DO j1 = 1, totnumneighprocs
330  k = neighproclist(j1)
331  CALL mpi_irecv(recvdatafrm(k)%rcvbuf(1),numndcomm(j1)*10, &
332  mpi_double_precision, k, 10, mpi_comm_rocfrac,reqrcv(j1),ierr)
333  ENDDO
334 !
335 !-MPI- SEND THE RECIPRICAL MASS MATRIX DIAGONAL TO THE NEIGHBORS
336 !
337  k2 = 1
338  DO j1 = 1, totnumneighprocs
339  k = neighproclist(j1)
340 
341  CALL mpi_isend(buf(k2),numndcomm(j1)*10,&
342  mpi_double_precision,k,10,mpi_comm_rocfrac,reqsnd(j1),ierr)
343  k2 = k2 + numndcomm(j1)*10
344  ENDDO
345 !
346 !-MPI- WAIT FOR INTERNAL FORCE VECTOR COMMUNICATION TO COMPLETE
347 !
348  IF(totnumneighprocs.GT.0)THEN
349  CALL mpi_waitall(totnumneighprocs,reqrcv,statrcv,ierr)
350  CALL mpi_waitall(totnumneighprocs,reqsnd,statsnd,ierr)
351  ENDIF
352  DEALLOCATE(buf)
353 
354 !
355 !----- ADD NEIGHBOR'S CONTRIBUTION TO THE INTERNAL FORCE VECTOR
356 !
357  DO j1 = 1, totnumneighprocs
358  k = neighproclist(j1)
359  k1 = 1
360  DO j = 1, numndcomm(j1)
361  k2 = neigh_lst(j1)%NdId(j)
362  f11v(k2) = f11v(k2) + recvdatafrm(k)%rcvbuf(k1)
363  f12v(k2) = f12v(k2) + recvdatafrm(k)%rcvbuf(k1+1)
364  f13v(k2) = f13v(k2) + recvdatafrm(k)%rcvbuf(k1+2)
365  f21v(k2) = f21v(k2) + recvdatafrm(k)%rcvbuf(k1+3)
366  f22v(k2) = f22v(k2) + recvdatafrm(k)%rcvbuf(k1+4)
367  f23v(k2) = f23v(k2) + recvdatafrm(k)%rcvbuf(k1+5)
368  f31v(k2) = f31v(k2) + recvdatafrm(k)%rcvbuf(k1+6)
369  f32v(k2) = f32v(k2) + recvdatafrm(k)%rcvbuf(k1+7)
370  f33v(k2) = f33v(k2) + recvdatafrm(k)%rcvbuf(k1+8)
371  voldefv(k2) = voldefv(k2) + recvdatafrm(k)%rcvbuf(k1+9)
372  k1 = k1 + 10
373  ENDDO
374  ENDDO
375 
376  DEALLOCATE(recvdatafrm)
377  DO i = 1, numnp ! for each node
378 
379  f11 = f11v(i)
380  f22 = f22v(i)
381  f33 = f33v(i)
382  f12 = f12v(i)
383  f13 = f13v(i)
384  f21 = f21v(i)
385  f23 = f23v(i)
386  f31 = f31v(i)
387  f32 = f32v(i)
388 
389  voldef = voldefv(i)
390 
391 ! WRONG FOR MORE THEN ONE MATERIAL
392 
393  CALL arruda_boyce_cauchy(f11,f12,f13,f21,f22,f23,f31,f32,f33, &
394  s11n(i),s22n(i),s33n(i),s12n(i),s13n(i),s23n(i),i,xmu(1),xkappa(1))
395  ENDDO
396 
397 
398  DO ielnum = 1, numel ! For each element
399 
400  n1 = nodes(1,ielnum)
401  n2 = nodes(2,ielnum)
402  n3 = nodes(3,ielnum)
403  n4 = nodes(4,ielnum)
404 
405  k3n1 = 3*n1
406  k3n2 = 3*n2
407  k3n3 = 3*n3
408  k3n4 = 3*n4
409 
410  k2n1 = k3n1 - 1
411  k2n2 = k3n2 - 1
412  k2n3 = k3n3 - 1
413  k2n4 = k3n4 - 1
414 
415  k1n1 = k3n1 - 2
416  k1n2 = k3n2 - 2
417  k1n3 = k3n3 - 2
418  k1n4 = k3n4 - 2
419 
420  ! k#n# dummy variables replaces:
421  u1 = disp(k1n1) ! (3*n1-2)
422  u2 = disp(k1n2) ! (3*n2-2)
423  u3 = disp(k1n3) ! (3*n3-2)
424  u4 = disp(k1n4) ! (3*n4-2)
425  v1 = disp(k2n1) ! (3*n1-1)
426  v2 = disp(k2n2) ! (3*n2-1)
427  v3 = disp(k2n3) ! (3*n3-1)
428  v4 = disp(k2n4) ! (3*n4-1)
429  w1 = disp(k3n1) ! (3*n1)
430  w2 = disp(k3n2) ! (3*n2)
431  w3 = disp(k3n3) ! (3*n3)
432  w4 = disp(k3n4) ! (3*n4)
433 
434  x1 = coor(1,n1)
435  x2 = coor(1,n2)
436  x3 = coor(1,n3)
437  x4 = coor(1,n4)
438  y1 = coor(2,n1)
439  y2 = coor(2,n2)
440  y3 = coor(2,n3)
441  y4 = coor(2,n4)
442  z1 = coor(3,n1)
443  z2 = coor(3,n2)
444  z3 = coor(3,n3)
445  z4 = coor(3,n4)
446 
447  x1 = x1 + u1
448  x2 = x2 + u2
449  x3 = x3 + u3
450  x4 = x4 + u4
451  y1 = y1 + v1
452  y2 = y2 + v2
453  y3 = y3 + v3
454  y4 = y4 + v4
455  z1 = z1 + w1
456  z2 = z2 + w2
457  z3 = z3 + w3
458  z4 = z4 + w4
459 
460  x12 = x1 - x2 ! not used in vol. calc
461  x13 = x1 - x3 ! not used in vol. calc
462  x14 = x1 - x4
463  x24 = x2 - x4
464  x34 = x3 - x4
465  y12 = y1 - y2 ! not used in vol. calc
466  y13 = y1 - y3 ! not used in vol. calc
467  y14 = y1 - y4
468  y24 = y2 - y4
469  y34 = y3 - y4
470  z12 = z1 - z2 ! not used in vol. calc
471  z13 = z1 - z3 ! not used in vol. calc
472  z14 = z1 - z4
473  z24 = z2 - z4
474  z34 = z3 - z4
475 
476  c11 = y24*z34 - z24*y34
477  c21 = -( x24*z34 - z24*x34 )
478  c31 = x24*y34 - y24*x34
479 
480  vx6def = -( x14*c11 + y14*c21 + z14*c31 )
481 
482  vx6invdef = 1.d0 / vx6def
483 
484 ! calculate the volume
485  voldef = vx6def/6.d0
486 
487  IF(vx6def.LE.0.d0) THEN
488  WRITE(*,100) i
489  stop
490  ENDIF
491 
492  ! equation (14)
493 
494  s11e = 0.d0
495  s22e = 0.d0
496  s33e = 0.d0
497  s12e = 0.d0
498  s23e = 0.d0
499  s13e = 0.d0
500 
501  DO k = 1, 4
502  s11e = s11e + alpha(k,ielnum)*s11n(nodes(k,ielnum))
503  s22e = s22e + alpha(k,ielnum)*s22n(nodes(k,ielnum))
504  s33e = s33e + alpha(k,ielnum)*s33n(nodes(k,ielnum))
505  s12e = s12e + alpha(k,ielnum)*s12n(nodes(k,ielnum))
506  s23e = s23e + alpha(k,ielnum)*s23n(nodes(k,ielnum))
507  s13e = s13e + alpha(k,ielnum)*s13n(nodes(k,ielnum))
508  ENDDO
509 
510  s11ev(ielnum) = s11e
511  s22ev(ielnum) = s22e
512  s33ev(ielnum) = s33e
513  s12ev(ielnum) = s12e
514  s23ev(ielnum) = s23e
515  s13ev(ielnum) = s13e
516 
517  sh1 = (y34*z24 - y24*z34) * vx6invdef
518  sh2 = (z34*x24 - z24*x34) * vx6invdef
519  sh3 = (x34*y24 - x24*y34) * vx6invdef
520  sh4 = (y13*z14 - y14*z13) * vx6invdef
521  sh5 = (z13*x14 - z14*x13) * vx6invdef
522  sh6 = (x13*y14 - x14*y13) * vx6invdef
523  sh7 = (y14*z12 - y12*z14) * vx6invdef
524  sh8 = (z14*x12 - z12*x14) * vx6invdef
525  sh9 = (x14*y12 - x12*y14) * vx6invdef
526  sh10 = (y12*z13 - y13*z12) * vx6invdef
527  sh11 = (z12*x13 - z13*x12) * vx6invdef
528  sh12 = (x12*y13 - x13*y12) * vx6invdef
529 
530 ! ASSEMBLE THE INTERNAL FORCE VECTOR
531 !
532 ! local node 1
533 
534  rnet(k1n1) = rnet(k1n1) - voldef* &
535  ( s11e*sh1 + s12e*sh2 + s13e*sh3 )
536  rnet(k2n1) = rnet(k2n1) - voldef* &
537  ( s12e*sh1 + s22e*sh2 + s23e*sh3 )
538  rnet(k3n1) = rnet(k3n1) - voldef* &
539  ( s13e*sh1 + s23e*sh2 + s33e*sh3 )
540 ! local node 2
541  rnet(k1n2) = rnet(k1n2) - voldef* &
542  ( s11e*sh4 + s12e*sh5 + s13e*sh6 )
543  rnet(k2n2) = rnet(k2n2) - voldef* &
544  ( s12e*sh4 + s22e*sh5 + s23e*sh6 )
545  rnet(k3n2) = rnet(k3n2) - voldef* &
546  ( s13e*sh4 + s23e*sh5 + s33e*sh6 )
547 ! local node 3
548  rnet(k1n3) = rnet(k1n3) - voldef* &
549  ( s11e*sh7 + s12e*sh8 + s13e*sh9 )
550  rnet(k2n3) = rnet(k2n3) - voldef* &
551  ( s12e*sh7 + s22e*sh8 + s23e*sh9 )
552  rnet(k3n3) = rnet(k3n3) - voldef* &
553  ( s13e*sh7 + s23e*sh8 + s33e*sh9 )
554 ! local node 4
555  rnet(k1n4) = rnet(k1n4) - voldef* &
556  ( s11e*sh10 + s12e*sh11 + s13e*sh12 )
557  rnet(k2n4) = rnet(k2n4) - voldef* &
558  ( s12e*sh10 + s22e*sh11 + s23e*sh12 )
559  rnet(k3n4) = rnet(k3n4) - voldef* &
560  ( s13e*sh10 + s23e*sh11 + s33e*sh12 )
561 
562  ENDDO
563 
564  RETURN
565 100 FORMAT(' Negative Jacobian for element: ',i10)
566 
567  END SUBROUTINE v3d4n_nl_arruda_boyce
568 
rocstar_rocfraccomm
Definition: ROCSTAR_RocFracComm.f90:53
k
j indices k indices k
Definition: Indexing.h:6
rocstar_rocfraccomm::rcv_buf
Definition: ROCSTAR_RocFracComm.f90:71
arruda_boyce_cauchy
subroutine arruda_boyce_cauchy(F11, F12, F13, F21, F22, F23, F31, F32, F33, Cchy11, Cchy22, Cchy33, Cchy12, Cchy13, Cchy23, ielem, mu, kappa)
Definition: arruda_boyce.cauchy.f90:53
i
blockLoc i
Definition: read.cpp:79
rocstar_rocfraccomm::send_buf
Definition: ROCSTAR_RocFracComm.f90:75
j
j indices j
Definition: Indexing.h:6
v3d4n_nl_arruda_boyce
subroutine v3d4n_nl_arruda_boyce(numnp, numel, coor, disp, nodes, Rnet, S11ev, S22ev, S33ev, S12ev, S23ev, S13ev, NumElNeigh, ElConn, alpha, Ahat, xmu, xkappa, NumMatVol, nprocs, TotNumNdComm, TotNumNeighProcs, NeighProcList, NumNdComm, neigh_lst, MPI_STATUS_SIZE, MPI_COMM_ROCFRAC, MPI_DOUBLE_PRECISION, ReqRcv, ReqSnd, StatRcv, StatSnd)
Definition: v3d4n_nl_arrudaboyce.f90:53
Color::alpha
unsigned char alpha() const
Definition: Color.h:75