56
57
58
60
61
62
63#include "implicit_f.inc"
64
65
66
67#include "mvsiz_p.inc"
68
69
70
71#include "param_c.inc"
72#include "com04_c.inc"
73
74
75
76 INTEGER, INTENT(IN) :: NEL
77 INTEGER, INTENT(IN) :: MTN
79 . pm(npropm,nummat), rho(*),volgp(lveul,*),
80 . f11(*),f21(*),f31(*),f12(*),f22(*),f32(*),
81 . f13(*),f23(*),f33(*),f14(*),f24(*),f34(*),
82 . f15(*),f25(*),f35(*),f16(*),f26(*),f36(*),
83 . f17(*),f27(*),f37(*),f18(*),f28(*),f38(*),
84 . px1(*),px2(*),px3(*),px4(*),
85 . py1(*),py2(*),py3(*),py4(*),
86 . pz1(*),pz2(*),pz3(*),pz4(*),
87 . px5(*),px6(*),px7(*),px8(*),
88 . py5(*),py6(*),py7(*),py8(*),
89 . pz5(*),pz6(*),pz7(*),pz8(*),
90 . dxx(*),dxy(*),dxz(*),dyx(*),dyy(*),dyz(*),dzx(*),dzy(*),dzz(*),
91 . vdx(*),vdy(*),vdz(*)
93 . qmv(12,*),bufmat(*),
94 . vx1(*),vx2(*),vx3(*),vx4(*),vx5(*),vx6(*),vx7(*),vx8(*),
95 . vy1(*),vy2(*),vy3(*),vy4(*),vy5(*),vy6(*),vy7(*),vy8
96 . vz1(*),vz2(*),vz3(*),vz4(*),vz5(*),vz6(*),vz7(*),vz8(*)
97 INTEGER IFL, MX
98 INTEGER MAT(*),IPARG1(*)
99 INTEGER,INTENT(IN)::IPM(NPROPMI,NUMMAT)
100
101
102
104 . f1(mvsiz),f2(mvsiz),f3(mvsiz),gam(mvsiz),
105 . a1(mvsiz),a2(mvsiz),a3(mvsiz),a4(mvsiz),
106 . a5(mvsiz),a6(mvsiz),a7(mvsiz),a8(mvsiz),fac
108 INTEGER IFLG,IADBUF
109 INTEGER I
110
111
112
113
114 IF(
ale%UPWIND%UPWM==0.AND.iparg1(64)==1)
RETURN
115
116
117
118
119
120 IF(mtn==51 .AND.
ale%UPWIND%UPWM/=3 .AND. iparg1(64)==0)
THEN
121 iadbuf = ipm(27,mat(1))
122 iflg = nint(bufmat(31+iadbuf-1))
123
124 IF(iflg>1)RETURN
125
126 IF(
ale%UPWIND%UPWM==0.)
THEN
127 mx = mat(1)
128 DO i=1,nel
129 gam(i)= pm(15,mx)
130 ENDDO
131 ELSE
132 DO i=1,nel
133 gam(i)=
ale%UPWIND%CUPWM
134 ENDDO
135 ENDIF
136
137 DO i=1,nel
138 aaa = qmv(1,i)+qmv(2,i)+qmv(3,i)+qmv(4,i)+qmv(5,i)+qmv(6,i)
139 aaa = one_over_6 * aaa
140
141 qmv(1,i) = one_over_8 * (qmv(1,i) - aaa)
142 qmv(2,i) = one_over_8 * (qmv(2,i) - aaa)
143 qmv(3,i) = one_over_8 * (qmv(3,i) - aaa)
144 qmv(4,i) = one_over_8 * (qmv(4,i) - aaa)
145 qmv(5,i) = one_over_8 * (qmv(5,i) - aaa)
146 qmv(6,i) = one_over_8 * (qmv(6,i) - aaa)
147 dmx = zero
148 dmy = zero
149 dmz = zero
150 dm = qmv(1,i)+qmv(6,i)+qmv(4,i)
151 dmx = dmx + vx1(i)*dm
152 dmy = dmy + vy1(i)*dm
153 dmz = dmz + vz1(i)*dm
154 dm = qmv(1,i)+qmv(4,i)+qmv(5,i)
155 dmx = dmx + vx2(i)*dm
156 dmy = dmy + vy2(i)*dm
157 dmz = dmz + vz2(i)*dm
158 dm = qmv(1,i)+qmv(5,i)+qmv(2,i)
159 dmx = dmx + vx3(i)*dm
160 dmy = dmy + vy3(i)*dm
161 dmz = dmz + vz3(i)*dm
162 dm = qmv(1,i)+qmv(2,i)+qmv(6,i)
163 dmx = dmx + vx4(i)*dm
164 dmy = dmy + vy4(i)*dm
165 dmz = dmz + vz4(i)*dm
166 dm = qmv(3,i)+qmv(6,i)+qmv(4,i)
167 dmx = dmx + vx5(i)*dm
168 dmy = dmy + vy5(i)*dm
169 dmz = dmz + vz5(i)*dm
170 dm = qmv(3,i)+qmv(4,i)+qmv(5,i)
171 dmx = dmx + vx6(i)*dm
172 dmy = dmy + vy6(i)*dm
173 dmz = dmz + vz6(i)*dm
174 dm = qmv(3,i)+qmv(5,i)+qmv(2,i)
175 dmx = dmx + vx7(i)*dm
176 dmy = dmy + vy7(i)*dm
177 dmz = dmz + vz7(i)*dm
178 dm = qmv(3,i)+qmv(2,i)+qmv(6,i)
179 dmx = dmx + vx8(i)*dm
180 dmy = dmy + vy8(i)*dm
181 dmz = dmz + vz8(i)*dm
182
183 dmx = -0.125 * dmx
184 dmy = -0.125 * dmy
185 dmz = -0.125 * dmz
186
187 a1(i) = px1(i)*vdx(i)+py1(i)*vdy(i)+pz1(i)*vdz(i)
188 a2(i) = px2(i)*vdx(i)+py2(i)*vdy(i)+pz2(i)*vdz(i)
189 a3(i) = px3(i)*vdx(i)+py3(i)*vdy(i)+pz3(i)*vdz(i)
190 a4(i) = px4(i)*vdx(i)+py4(i)*vdy(i)+pz4(i)*vdz(i)
191 a1(i) = sign(gam(i),a1(i))
192 a2(i) = sign(gam(i),a2(i))
193 a3(i) = sign(gam(i),a3(i))
194 a4(i) = sign(gam(i),a4(i))
195
196 f11(i) = f11(i) - (one+a1(i))*dmx
197 f12(i) = f12(i) - (one+a2(i))*dmx
198 f13(i) = f13(i) - (one+a3(i))*dmx
199 f14(i) = f14(i) - (one+a4(i))*dmx
200 f15(i) = f15(i) - (one-a3(i))*dmx
201 f16(i) = f16(i) - (one-a4(i))*dmx
202 f17(i) = f17(i) - (one-a1(i))*dmx
203 f18(i) = f18(i) - (one-a2(i))*dmx
204
205 f21(i) = f21(i) - (one+a1(i))*dmy
206 f22(i) = f22(i) - (one+a2(i))*dmy
207 f23(i) = f23(i) - (one+a3(i))*dmy
208 f24(i) = f24(i) - (one+a4(i))*dmy
209 f25(i) = f25(i) - (one-a3(i))*dmy
210 f26(i) = f26(i) - (one-a4(i))*dmy
211 f27(i) = f27(i) - (one-a1(i))*dmy
212 f28(i) = f28(i) - (one-a2(i))*dmy
213
214 f31(i) = f31(i) - (one+a1(i))*dmz
215 f32(i) = f32(i) - (one+a2(i))*dmz
216 f33(i) = f33(i) - (one+a3(i))*dmz
217 f34(i) = f34(i) - (one+a4(i))*dmz
218 f35(i) = f35(i) - (one-a3(i))*dmz
219 f36(i) = f36(i) - (one-a4(i))*dmz
220 f37(i) = f37(i) - (one-a1(i))*dmz
221 f38(i) = f38(i) - (one-a2(i))*dmz
222 ENDDO
223 RETURN
224 ENDIF
225
226
227
228
229 mx = mat(1)
230 DO i=1,nel
231 gam(i)= pm(15,mx)
232 ENDDO
233 DO i=1,nel
234 a1(i) = px1(i)*vdx(i)+py1(i)*vdy(i)+pz1(i)*vdz(i)
235 a2(i) = px2(i)*vdx(i)+py2(i)*vdy(i)+pz2(i)*vdz(i)
236 a3(i) = px3(i)*vdx(i)+py3(i)*vdy(i)+pz3(i)*vdz(i)
237 a4(i) = px4(i)*vdx(i)+py4(i)*vdy(i)+pz4(i)*vdz(i)
238 a5(i) = px5(i)*vdx(i)+py5(i)*vdy(i)+pz5(i)*vdz(i)
239 a6(i) = px6(i)*vdx(i)+py6(i)*vdy(i)+pz6(i)*vdz(i)
240 a7(i) = px7(i)*vdx(i)+py7(i)*vdy(i)+pz7(i)*vdz(i)
241 a8(i) = px8(i)*vdx(i)+py8(i)*vdy(i)+pz8(i)*vdz(i)
242 a1(i) = sign(gam(i),a1(i))
243 a2(i) = sign(gam(i),a2(i))
244 a3(i) = sign(gam(i),a3(i))
245 a4(i) = sign(gam(i),a4(i))
246 a5(i) = sign(gam(i),a5(i))
247 a6(i) = sign(gam(i),a6(i))
248 a7(i) = sign(gam(i),a7(i))
249 a8(i) = sign(gam(i),a8(i))
250 ENDDO
251 DO i=1,nel
252 fac = rho(i)
253 f1(i) = (vdx(i)*dxx(i)+vdy(i)*dxy(i)+vdz(i)*dxz(i))*fac
254 f2(i) = (vdx(i)*dyx(i)+vdy(i)*dyy(i)+vdz(i)*dyz(i))*fac
255 f3(i) = (vdx(i)*dzx(i)+vdy(i)*dzy(i)+vdz(i)*dzz(i))*fac
256 ENDDO
257 DO i=1,nel
258 f11(i) = f11(i) - (one+a1(i))*f1(i)*volgp(1,i)
259 f12(i) = f12(i) - (one+a2(i))*f1(i)*volgp(2,i)
260 f13(i) = f13(i) - (one+a3(i))*f1(i)*volgp(3,i)
261 f14(i) = f14(i) - (one+a4(i))*f1(i)*volgp(4,i)
262 f15(i) = f15(i) - (one+a5(i))*f1(i)*volgp(5,i)
263 f16(i) = f16(i) - (one+a6(i))*f1(i)*volgp(6,i)
264 f17(i) = f17(i) - (one+a7(i))*f1(i)*volgp(7,i)
265 f18(i) = f18(i) - (one+a8(i))*f1(i)*volgp(8,i)
266
267 f21(i) = f21(i) - (one+a1(i))*f2(i)*volgp(1,i)
268 f22(i) = f22(i) - (one+a2(i))*f2(i)*volgp(2,i)
269 f23(i) = f23(i) - (one+a3(i))*f2(i)*volgp(3,i)
270 f24(i) = f24(i) - (one+a4(i))*f2(i)*volgp(4,i)
271 f25(i) = f25(i) - (one+a5(i))*f2(i)*volgp(5,i)
272 f26(i) = f26(i) - (one+a6(i))*f2(i)*volgp(6,i)
273 f27(i) = f27(i) - (one+a7(i))*f2(i)*volgp(7,i)
274 f28(i) = f28(i) - (one+a8(i))*f2(i)*volgp(8,i)
275
276 f31(i) = f31(i) - (one+a1(i))*f3(i)*volgp(1,i)
277 f32(i) = f32(i) - (one+a2(i))*f3(i)*volgp(2,i)
278 f33(i) = f33(i) - (one+a3(i))*f3(i)*volgp(3,i)
279 f34(i) = f34(i) - (one+a4(i))*f3(i)*volgp(4,i)
280 f35(i) = f35(i) - (one+a5(i))*f3(i)*volgp(5,i)
281 f36(i) = f36(i) - (one+a6(i))*f3(i)*volgp(6,i)
282 f37(i) = f37(i) - (one+a7(i))*f3(i)*volgp(7,i)
283 f38(i) = f38(i) - (one+a8(i))*f3(i)*volgp(8,i)
284 ENDDO
285
286 RETURN