35
36
37
39
40
41
42#include "implicit_f.inc"
43
44
45
46#include "mvsiz_p.inc"
47
48
49
50 INTEGER NEL
52 . jac_i(10,*)
53 double precision
54 . x1(*), x2(*), x3(*), x4(*),
55 . y1(*), y2(*), y3(*), y4(*),
56 . z1(*), z2(*), z3(*), z4(*)
57
58
59
60 INTEGER I
62 . b1(mvsiz), b2(mvsiz), b3(mvsiz), b4(mvsiz),
63 . c1(mvsiz), c2(mvsiz), c3(mvsiz), c4(mvsiz),
64 . d1(mvsiz), d2(mvsiz), d3(mvsiz), d4(mvsiz),vol(mvsiz)
65 double precision
66 . x41, y41, z41, x42, y42, z42, x43, y43, z43,b1dp,c1dp,d1dp,voldp
68 . d ,px1,px2,px3,px4,py1,py2,py3,py4,pz1,pz2,pz3,pz4
69
70 DO i=1,nel
71 x43 = x4(i) - x3(i)
72 y43 = y4(i) - y3(i)
73 z43 = z4(i) - z3(i)
74 x41 = x4(i) - x1(i)
75 y41 = y4(i) - y1(i)
76 z41 = z4(i) - z1(i)
77 x42 = x4(i) - x2(i)
78 y42 = y4(i) - y2(i)
79 z42 = z4(i) - z2(i)
80
81 b1dp = y43*z42 - y42*z43
82 b1(i) = b1dp
83 b2(i) = y41*z43 - y43*z41
84 b3(i) = y42*z41 - y41*z42
85 b4(i) = -(b1(i) + b2(i) + b3(i))
86
87 c1dp = z43*x42 - z42*x43
88 c1(i) = c1dp
89 c2(i) = z41*x43 - z43*x41
90 c3(i) = z42*x41 - z41*x42
91 c4(i) = -(c1(i) + c2(i) + c3(i))
92
93 d1dp = x43*y42 - x42*y43
94 d1(i) = d1dp
95 d2(i) = x41*y43 - x43*y41
96 d3(i) = x42*y41 - x41*y42
97 d4(i) = -(d1(i) + d2(i) + d3(i))
98
99 voldp = (x41*b1dp + y41*c1dp + z41*d1dp)*one_over_6
100 vol(i) = voldp
101 ENDDO
102
103 DO i=1,nel
104 d = one/vol(i)/six
105 px1=-b1(i)*d
106 py1=-c1(i)*d
107 pz1=-d1(i)*d
108 px2=-b2(i)*d
109 py2=-c2(i)*d
110 pz2=-d2(i)*d
111 px3=-b3(i)*d
112 py3=-c3(i)*d
113 pz3=-d3(i)*d
114 px4=-b4(i)*d
115 py4=-c4(i)*d
116 pz4=-d4(i)*d
117 jac_i(1,i) = px1
118 jac_i(2,i) = px2
119 jac_i(3,i) = px3
120 jac_i(4,i) = py1
121 jac_i(5,i) = py2
122 jac_i(6,i) = py3
123 jac_i(7,i) = pz1
124 jac_i(8,i) = pz2
125 jac_i(9,i) = pz3
126 jac_i(10,i) = vol(i)
127 ENDDO
128
129 RETURN