32
33
34
36 use element_mod , only : nixs
37
38
39
40#include "implicit_f.inc"
41
42
43
44#include "mvsiz_p.inc"
45
46
47
48#include "param_c.inc"
49#include "vect01_c.inc"
50#include "tabsiz_c.inc"
51
52
53
54
55
56
57
58
59
60
61
62 INTEGER,INTENT(IN) :: IXS(NIXS,SIXS/NIXS)
63 my_real,
INTENT(IN) :: x(3,sx/3), veul(lveul,*)
64 my_real,
INTENT(INOUT) :: grad(6,*)
65 TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECT
66
67
68
69 INTEGER I, II, IE, IV1, IV2, IV3, IV4, IV5, IV6
70 my_real x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz), x5(mvsiz), x6(mvsiz), x7(mvsiz), x8(mvsiz),
71 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz), y5(mvsiz), y6(mvsiz), y7(mvsiz), y8(mvsiz), z1(mvsiz),
72 . z2(mvsiz), z3(mvsiz), z4(mvsiz), z5(mvsiz), z6(mvsiz), z7(mvsiz), z8(mvsiz), xc(mvsiz), yc(mvsiz),
73 . zc(mvsiz), n1x(mvsiz), n2x(mvsiz), n3x(mvsiz), n4x(mvsiz), n5x(mvsiz), n6x(mvsiz), n1y(mvsiz),
74 . n2y(mvsiz), n3y(mvsiz), n4y(mvsiz), n5y(mvsiz), n6y(mvsiz), n1z(mvsiz), n2z(mvsiz), n3z(mvsiz),
75 . n4z(mvsiz), n5z(mvsiz), n6z(mvsiz), dd1(mvsiz), dd2(mvsiz), dd3(mvsiz), dd4(mvsiz), dd5(mvsiz),
76 . dd6(mvsiz), d1x(mvsiz), d2x(mvsiz), d3x(mvsiz), d4x(mvsiz), d5x(mvsiz), d6x(mvsiz), d1y(mvsiz),
77 . d2y(mvsiz), d3y(mvsiz), d4y(mvsiz), d5y(mvsiz), d6y(mvsiz), d1z(mvsiz), d2z(mvsiz), d3z(mvsiz),
78 . d4z(mvsiz), d5z(mvsiz), d6z(mvsiz)
79 INTEGER :: IAD2
80
81
82
83 DO i=lft,llt
84 ii=i+nft
85 x1(i)=x(1,ixs(2,ii))
86 y1(i)=x(2,ixs(2,ii))
87 z1(i)=x(3,ixs(2,ii))
88
89 x2(i)=x(1,ixs(3,ii))
90 y2(i)=x(2,ixs(3,ii))
91 z2(i)=x(3,ixs(3,ii))
92
93 x3(i)=x(1,ixs(4,ii))
94 y3(i)=x(2,ixs(4,ii))
95 z3(i)=x(3,ixs(4,ii))
96
97 x4(i)=x(1,ixs(5,ii))
98 y4(i)=x(2,ixs(5,ii))
99 z4(i)=x(3,ixs(5,ii))
100
101 x5(i)=x(1,ixs(6,ii))
102 y5(i)=x(2,ixs(6,ii))
103 z5(i)=x(3,ixs(6,ii))
104
105 x6(i)=x(1,ixs(7,ii))
106 y6(i)=x(2,ixs(7,ii))
107 z6(i)=x(3,ixs(7,ii))
108
109 x7(i)=x(1,ixs(8,ii))
110 y7(i)=x(2,ixs(8,ii))
111 z7(i)=x(3,ixs(8,ii))
112
113 x8(i)=x(1,ixs(9,ii))
114 y8(i)=x(2,ixs(9,ii))
115 z8(i)=x(3,ixs(9,ii))
116 ENDDO
117
118
119
120 DO i=lft,llt
121 ii=i+nft
122
123 n1x(i)=veul(14,ii)
124 n2x(i)=veul(15,ii)
125 n3x(i)=veul(16,ii)
126 n4x(i)=veul(17,ii)
127 n5x(i)=veul(18,ii)
128 n6x(i)=veul(19,ii)
129
130 n1y(i)=veul(20,ii)
131 n2y(i)=veul(21,ii)
132 n3y(i)=veul(22,ii)
133 n4y(i)=veul(23,ii)
134 n5y(i)=veul(24,ii)
135 n6y(i)=veul(25,ii)
136
137 n1z(i)=veul(26,ii)
138 n2z(i)=veul(27,ii)
139 n3z(i)=veul(28,ii)
140 n4z(i)=veul(29,ii)
141 n5z(i)=veul(30,ii)
142 n6z(i)=veul(31,ii)
143
144 xc(i) = (x1(i)+x2(i)+x3(i)+x4(i)+x5(i)+x6(i)+x7(i)+x8(i))
145 yc(i) = (y1(i)+y2(i)+y3(i)+y4(i)+y5(i)+y6(i)+y7(i)+y8(i))
146 zc(i) = (z1(i)+z2(i)+z3(i)+z4(i)+z5(i)+z6(i)+z7(i)+z8(i))
147 ENDDO
148
149
150
151 DO i=lft,llt
152 ie =nft+i
153 iad2 = ale_connect%ee_connect%iad_connect(ie)
154 iv1 = ale_connect%ee_connect%connected(iad2 + 1 - 1)
155 IF(iv1 <= 0) iv1=ie
156 d1x(i) = - xc(i)
157 . +x(1,ixs(2,iv1))+x(1,ixs(3,iv1))+x(1,ixs(4,iv1))+x(1,ixs(5,iv1))
158 . +x(1,ixs(6,iv1))+x(1,ixs(7,iv1))+x(1,ixs(8,iv1))+x(1,ixs(9,iv1))
159 d1y(i) = - yc(i)
160 . +x(2,ixs(2,iv1))+x(2,ixs(3,iv1))+x(2,ixs(4,iv1))+x(2,ixs(5,iv1))
161 . +x(2,ixs(6,iv1))+x(2,ixs(7,iv1))+x(2,ixs(8,iv1))+x(2,ixs(9,iv1))
162 d1z(i) = - zc(i)
163 . +x(3,ixs(2,iv1))+x(3,ixs(3,iv1))+x(3,ixs(4,iv1))+x(3,ixs(5,iv1))
164 . +x(3,ixs(6,iv1))+x(3,ixs(7,iv1))+x(3,ixs(8,iv1))+x(3,ixs(9,iv1))
165 ENDDO
166 DO i=lft,llt
167 ie =nft+i
168 iad2 = ale_connect%ee_connect%iad_connect(ie)
169 iv2 = ale_connect%ee_connect%connected(iad2 + 2 - 1)
170 IF(iv2 <= 0) iv2=ie
171 d2x(i) = - xc(i)
172 . +x(1,ixs(2,iv2))+x(1,ixs(3,iv2))+x(1,ixs(4,iv2))+x(1,ixs(5,iv2))
173 . +x(1,ixs(6,iv2))+x(1,ixs(7,iv2))+x(1,ixs(8,iv2))+x(1,ixs(9,iv2))
174 d2y(i) = - yc(i)
175 . +x(2,ixs(2,iv2))+x(2,ixs(3,iv2))+x(2,ixs(4,iv2))+x(2,ixs(5,iv2))
176 . +x(2,ixs(6,iv2))+x(2,ixs(7,iv2))+x(2,ixs(8,iv2))+x(2,ixs(9,iv2))
177 d2z(i) = - zc(i)
178 . +x(3,ixs(2,iv2))+x(3,ixs(3,iv2))+x(3,ixs(4,iv2))+x(3,ixs(5,iv2))
179 . +x(3,ixs(6,iv2))+x(3,ixs(7,iv2))+x(3,ixs(8,iv2))+x(3,ixs(9,iv2))
180 ENDDO
181 DO i=lft,llt
182 ie =nft+i
183 iad2 = ale_connect%ee_connect%iad_connect(ie)
184 iv3 = ale_connect%ee_connect%connected(iad2 + 3 - 1)
185 IF(iv3 <= 0) iv3=ie
186 d3x(i) = - xc(i)
187 . +x(1,ixs(2,iv3))+x(1,ixs(3,iv3))+x(1,ixs(4,iv3))+x(1,ixs(5,iv3))
188 . +x(1,ixs(6,iv3))+x(1,ixs(7,iv3))+x(1,ixs(8,iv3))+x(1,ixs(9,iv3))
189 d3y(i) = - yc(i)
190 . +x(2,ixs(2,iv3))+x(2,ixs(3,iv3))+x(2,ixs(4,iv3))+x(2,ixs(5,iv3))
191 . +x(2,ixs(6,iv3))+x(2,ixs(7,iv3))+x(2,ixs(8,iv3))+x(2,ixs(9,iv3))
192 d3z(i) = - zc(i)
193 . +x(3,ixs(2,iv3))+x(3,ixs(3,iv3))+x(3,ixs(4,iv3))+x(3,ixs(5,iv3))
194 . +x(3,ixs(6,iv3))+x(3,ixs(7,iv3))+x(3,ixs(8,iv3))+x(3,ixs(9,iv3))
195 ENDDO
196 DO i=lft,llt
197 ie =nft+i
198 iad2 = ale_connect%ee_connect%iad_connect(ie)
199 iv4 = ale_connect%ee_connect%connected(iad2 + 4 - 1)
200 IF(iv4 <= 0) iv4=ie
201 d4x(i) = - xc(i)
202 . +x(1,ixs(2,iv4))+x(1,ixs(3,iv4))+x(1,ixs(4,iv4))+x(1,ixs(5,iv4))
203 . +x(1,ixs(6,iv4))+x(1,ixs(7,iv4))+x(1,ixs(8,iv4))+x(1,ixs(9,iv4))
204 d4y(i) = - yc(i)
205 . +x(2,ixs(2,iv4))+x(2,ixs(3,iv4))+x(2,ixs(4,iv4))+x(2,ixs(5,iv4))
206 . +x(2,ixs(6,iv4))+x(2,ixs(7,iv4))+x(2,ixs(8,iv4))+x(2,ixs(9,iv4))
207 d4z(i) = - zc(i)
208 . +x(3,ixs(2,iv4))+x(3,ixs(3,iv4))+x(3,ixs(4,iv4))+x(3,ixs(5,iv4))
209 . +x(3,ixs(6,iv4))+x(3,ixs(7,iv4))+x(3,ixs(8,iv4))+x(3,ixs(9,iv4))
210 ENDDO
211 DO i=lft,llt
212 ie =nft+i
213 iad2 = ale_connect%ee_connect%iad_connect(ie)
214 iv5 = ale_connect%ee_connect%connected(iad2 + 5 - 1)
215 IF(iv5 <= 0) iv5=ie
216 d5x(i) = - xc(i)
217 . +x(1,ixs(2,iv5))+x(1,ixs(3,iv5))+x(1,ixs(4,iv5))+x(1,ixs(5,iv5))
218 . +x(1,ixs(6,iv5))+x(1,ixs(7,iv5))+x(1,ixs(8,iv5))+x(1,ixs(9,iv5))
219 d5y(i) = - yc(i)
220 . +x(2,ixs(2,iv5))+x(2,ixs(3,iv5))+x(2,ixs(4,iv5))+x(2,ixs(5,iv5))
221 . +x(2,ixs(6,iv5))+x(2,ixs(7,iv5))+x(2,ixs(8,iv5))+x(2,ixs(9,iv5))
222 d5z(i) = - zc(i)
223 . +x(3,ixs(2,iv5))+x(3,ixs(3,iv5))+x(3,ixs(4,iv5))+x(3,ixs(5,iv5))
224 . +x(3,ixs(6,iv5))+x(3,ixs(7,iv5))+x(3,ixs(8,iv5))+x(3,ixs(9,iv5))
225 ENDDO
226 DO i=lft,llt
227 ie =nft+i
228 iad2 = ale_connect%ee_connect%iad_connect(ie)
229 iv6 = ale_connect%ee_connect%connected(iad2 + 6 - 1)
230 IF(iv6 <= 0) iv6=ie
231 d6x(i) = - xc(i)
232 . +x(1,ixs(2,iv6))+x(1,ixs(3,iv6))+x(1,ixs(4,iv6))+x(1,ixs(5,iv6))
233 . +x(1,ixs(6,iv6))+x(1,ixs(7,iv6))+x(1,ixs(8,iv6))+x(1,ixs(9,iv6))
234 d6y(i) = - yc(i)
235 . +x(2,ixs(2,iv6))+x(2,ixs(3,iv6))+x(2,ixs(4,iv6))+x(2,ixs(5,iv6))
236 . +x(2,ixs(6,iv6))+x(2,ixs(7,iv6))+x(2,ixs(8,iv6))+x(2,ixs(9,iv6))
237 d6z(i) = - zc(i)
238 . +x(3,ixs(2,iv6))+x(3,ixs(3,iv6))+x(3,ixs(4,iv6))+x(3,ixs(5,iv6))
239 . +x(3,ixs(6,iv6))+x(3,ixs(7,iv6))+x(3,ixs(8,iv6))+x(3,ixs(9,iv6))
240 ENDDO
241
242
243
244
245
246 DO i=lft,llt
247 dd1(i)=d1x(i)**2+d1y(i)**2+d1z(i)**2
248 dd2(i)=d2x(i)**2+d2y(i)**2+d2z(i)**2
249 dd3(i)=d3x(i)**2+d3y(i)**2+d3z(i)**2
250 dd4(i)=d4x(i)**2+d4y(i)**2+d4z(i)**2
251 dd5(i)=d5x(i)**2+d5y(i)**2+d5z(i)**2
252 dd6(i)=d6x(i)**2+d6y(i)**2+d6z(i)**2
253 ENDDO
254
255
256
257
258
259
260
261
262
263 DO i=lft,llt
264 grad(1,i)= four*(d1x(i)*n1x(i)+d1y(i)*n1y(i)+d1z(i)*n1z(i)) /
max(em15,dd1(i))
265 grad(2,i)= four*(d2x(i)*n2x(i)+d2y(i)*n2y(i)+d2z(i)*n2z(i)) /
max(em15,dd2(i))
266 grad(3,i)= four*(d3x(i)*n3x(i)+d3y(i)*n3y(i)+d3z(i)*n3z(i)) /
max(em15,dd3(i))
267 grad(4,i)= four*(d4x(i)*n4x(i)+d4y(i)*n4y(i)+d4z(i)*n4z(i)) /
max(em15,dd4(i))
268 grad(5,i)= four*(d5x(i)*n5x(i)+d5y(i)*n5y(i)+d5z(i)*n5z(i)) /
max(em15,dd5(i))
269 grad(6,i)= four*(d6x(i)*n6x(i)+d6y(i)*n6y(i)+d6z(i)*n6z(i)) /
max(em15,dd6(i))
270 ENDDO
271
272 RETURN