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