31
32
33
35
36
37
38#include "implicit_f.inc"
39
40
41
42#include "mvsiz_p.inc"
43
44
45
46#include "vect01_c.inc"
47#include "param_c.inc"
48#include "com04_c.inc"
49
50
51
52 INTEGER IXS(NIXS,NUMELS)
53 my_real pm(npropm,nummat), v(3,numnod), w(3,numnod), x(3,numnod), flux(6,*), flu1(*)
54 TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECT
55
56
57
58 INTEGER MAT(MVSIZ),
59 . NC1(MVSIZ), NC2(MVSIZ), NC3(), NC4(MVSIZ), ,,II,IAD2
61 . x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz),
62 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz),
63 . z1(mvsiz), z2(mvsiz), z3(mvsiz), z4(mvsiz),
64 . n1x(mvsiz),n1y(mvsiz),n1z(mvsiz),
65 . n2x(mvsiz),n2y(mvsiz),n2z(mvsiz),
66 . n3x(mvsiz),n3y(mvsiz),n3z(mvsiz),
67 . n4x(mvsiz),n4y(mvsiz),n4z(mvsiz),
68 . flux2(mvsiz), flux4(mvsiz), flux5(mvsiz), flux6(mvsiz),
69 . vx1(mvsiz), vx2(mvsiz), vx3(mvsiz),
70 . vx4(mvsiz),
71 . vy1(mvsiz), vy2(mvsiz), vy3(mvsiz),
72 . vy4(mvsiz),
73 . vz1(mvsiz), vz2(mvsiz), vz3(mvsiz),
74 . vz4(mvsiz),
75 . vdx1(mvsiz),vdx2(mvsiz),vdx3(mvsiz),vdx4(mvsiz),
76 . vdy1(mvsiz),vdy2(mvsiz),vdy3(mvsiz),vdy4(mvsiz),
77 . vdz1(mvsiz),vdz2(mvsiz),vdz3(mvsiz),vdz4(mvsiz),
78 . reduc, upwl(6,mvsiz)
79
80
81
82
83
84
85 DO i=lft,llt
86 ii=i+nft
87 mat(i)=ixs(1,ii)
88
89 nc1(i)=ixs(2,ii)
90 nc2(i)=ixs(4,ii)
91 nc3(i)=ixs(7,ii)
92 nc4(i)=ixs(6,ii)
93
94
95 x1(i)=x(1,nc1(i))
96 y1(i)=x(2,nc1(i))
97 z1(i)=x(3,nc1(i))
98
99 x2(i)=x(1,nc2(i))
100 y2(i)=x(2,nc2(i))
101 z2(i)=x(3,nc2(i))
102
103 x3(i)=x(1,nc3(i))
104 y3(i)=x(2,nc3(i))
105 z3(i)=x(3,nc3(i))
106
107 x4(i)=x(1,nc4(i))
108 y4(i)=x(2,nc4(i))
109 z4(i)=x(3,nc4(i))
110
111
112
113
114
115
116
117 vdx1(i)=v(1,nc1(i)) - w(1,nc1(i))
118 vdy1(i)=v(2,nc1(i)) - w(2,nc1(i))
119 vdz1(i)=v(3,nc1(i)) - w(3,nc1(i))
120
121 vdx2(i)=v(1,nc2(i)) - w(1,nc2(i))
122 vdy2(i)=v(2,nc2(i)) - w(2,nc2(i))
123 vdz2(i)=v(3,nc2(i)) - w(3,nc2(i))
124
125 vdx3(i)=v(1,nc3(i)) - w(1,nc3(i))
126 vdy3(i)=v(2,nc3(i)) - w(2,nc3(i))
127 vdz3(i)=v(3,nc3(i)) - w(3,nc3(i))
128
129 vdx4(i)=v(1,nc4(i)) - w(1,nc4(i))
130 vdy4(i)=v(2,nc4(i)) - w(2,nc4(i))
131 vdz4(i)=v(3,nc4(i)) - w(3,nc4(i))
132 ENDDO
133
134
135
136
137
138
139
140 DO i=lft,llt
141
142 vx1(i)=(vdx1(i)+vdx2(i)+vdx3(i))
143 vx2(i)=(vdx4(i)+vdx2(i)+vdx1(i))
144 vx3(i)=(vdx4(i)+vdx3(i)+vdx2(i))
145 vx4(i)=(vdx4(i)+vdx1(i)+vdx3(i))
146
147 vy1(i)=(vdy1(i)+vdy2(i)+vdy3(i))
148 vy2(i)=(vdy4(i)+vdy2(i)+vdy1(i))
149 vy3(i)=(vdy4(i)+vdy3(i)+vdy2(i))
150 vy4(i)=(vdy4(i)+vdy1(i)+vdy3(i))
151
152 vz1(i)=(vdz1(i)+vdz2(i)+vdz3(i))
153 vz2(i)=(vdz4(i)+vdz2(i)+vdz1(i))
154 vz3(i)=(vdz4(i)+vdz3(i)+vdz2(i))
155 vz4(i)=(vdz4(i)+vdz1(i)+vdz3(i))
156 END DO
157
158
159
160
161
162
163
164 DO i=lft,llt
165
166 n2x(i)=-(y1(i)-y4(i))*(z2(i)-z4(i)) +(z1(i)-z4(i))*(y2(i)-y4(i))
167 n2y(i)=-(z1(i)-z4(i))*(x2(i)-x4(i)) +(x1(i)-x4(i))*(z2(i)-z4(i))
168 n2z(i)=-(x1(i)-x4(i))*(y2(i)-y4(i)) +(y1(i)-y4(i))*(x2(i)-x4(i))
169
170 n3x(i)=-(y2(i)-y4(i))*(z3(i)-z4(i)) +(z2(i)-z4(i))*(y3(i)-y4(i))
171 n3y(i)=-(z2(i)-z4(i))*(x3(i)-x4(i)) +(x2(i)-x4(i))*(z3(i)-z4(i))
172 n3z(i)=-(x2(i)-x4(i))*(y3(i)-y4(i)) +(y2(i)-y4(i))*(x3(i)-x4(i))
173
174 n4x(i)=-(y3(i)-y4(i))*(z1(i)-z4(i)) +(z3(i)-z4(i))*(y1(i)-y4(i))
175 n4y(i)=-(z3(i)-z4(i))*(x1(i)-x4(i)) +(x3(i)-x4(i))*(z1(i)-z4(i))
176 n4z(i)=-(x3(i)-x4(i))*(y1(i)-y4(i)) +(y3(i)-y4(i))*(x1(i)-x4(i))
177
178 n1x(i)= -n2x(i) -n3x(i) -n4x(i)
179 n1y(i)= -n2y(i) -n3y(i) -n4y(i)
180 n1z(i)= -n2z(i) -n3z(i) -n4z(i)
181 END DO
182
183
184
185
186
187
188 DO i=lft,llt
189 flux5(i)=(vx1(i)*n1x(i)+vy1(i)*n1y(i)+vz1(i)*n1z(i))*one_over_6
190 flux6(i)=(vx2(i)*n2x(i)+vy2(i)*n2y(i)+vz2(i)*n2z(i))*one_over_6
191 flux2(i)=(vx3(i)*n3x(i)+vy3(i)*n3y(i)+vz3(i)*n3z(i))*one_over_6
192 flux4(i)=(vx4(i)*n4x(i)+vy4(i)*n4y(i)+vz4(i)*n4z(i))*one_over_6
193 END DO
194
195
196
197
198
199 DO j=1,6
200 DO i=lft,llt
201 upwl(j,i)=pm(16,mat(i))
202 END DO
203 END DO
204
205
206
207
208
209
210
211
212
213
214 DO i=lft,llt
215 reduc=pm(92,mat(i))
216
217 iad2 = ale_connect%ee_connect%iad_connect(i + nft)
218 ii = ale_connect%ee_connect%connected(iad2 + 3 - 1)
219 IF(ii == 0)THEN
220 flux2(i)=flux2(i)*reduc
221 ELSEIF(ii > 0)THEN
222 IF(int(pm(19,ixs(1,ii))) == 11)THEN
223 upwl(2,i)=one
224 flux2(i)=flux2(i)*pm(92,ixs(1,ii))
225 ENDIF
226 ENDIF
227
228 ii = ale_connect%ee_connect%connected(iad2 + 4 - 1)
229 IF(ii == 0)THEN
230 flux4(i)=flux4(i)*reduc
231 ELSEIF(ii > 0)THEN
232 IF(int(pm(19,ixs(1,ii))) == 11)THEN
233 upwl(4,i)=one
234 flux4(i)=flux4(i)*pm(92,ixs(1,ii))
235 ENDIF
236 ENDIF
237
238 ii = ale_connect%ee_connect%connected(iad2 + 1 - 1)
239 IF(ii == 0)THEN
240 flux5(i)=flux5(i)*reduc
241 ELSEIF(ii > 0)THEN
242 IF(int(pm(19,ixs(1,ii))) == 11)THEN
243 upwl(5,i)=one
244 flux5(i)=flux5(i)*pm(92,ixs(1,ii))
245 ENDIF
246 ENDIF
247
248 ii = ale_connect%ee_connect%connected(iad2 + 2 - 1)
249 IF(ii == 0)THEN
250 flux6(i)=flux6(i)*reduc
251 ELSEIF(ii > 0)THEN
252 IF(int(pm(19,ixs(1,ii))) == 11)THEN
253 upwl(6,i)=one
254 flux6(i)=flux6(i)*pm(92,ixs(1,ii))
255 ENDIF
256 ENDIF
257 END DO
258
259
260
261
262
263
264
265
266
267
268
269 DO i=lft,llt
270
271 flux(1,i)=zero
272 flux(2,i)=flux2(i)-upwl(2,i)*abs(flux2(i))
273 flux(3,i)=zero
274 flux(4,i)=flux4(i)-upwl(4,i)*abs(flux4(i))
275 flux(5,i)=flux5(i)-upwl(5,i)*abs(flux5(i))
276 flux(6,i)=flux6(i)-upwl(6,i)*abs(flux6(i))
277
278 flu1(i) =flux2(i)+upwl(2,i)*abs(flux2(i))
279 . +flux4(i)+upwl(4,i)*abs(flux4(i))
280 . +flux5(i)+upwl(5,i)*abs(flux5(i))
281 . +flux6(i)+upwl(6,i)*abs(flux6(i))
282 END DO
283
284 RETURN