38
39
40
41 USE sensor_mod
42
43
44
45#include "implicit_f.inc"
46
47
48
49#include "mvsiz_p.inc"
50
51
52
53#include "param_c.inc"
54#include "vect01_c.inc"
55
56
57
58 INTEGER IEXPAN,
59 . IPARTS(*),GRTH(*),IGRTH(*),IPARG(*),ITASK
60
62 . partsav(npsav,*), eintg(*), rhog(*),
63 . vx(mvsiz,*),vy(mvsiz,*),vz(mvsiz,*),
64 . volg(*),vol0g(*),gresav(*),eintth(*),
65 . fill(*),
66 . x(mvsiz,*),y(mvsiz,*),z(mvsiz,*)
67 my_real,
INTENT(IN) :: offg(mvsiz)
68 INTEGER, INTENT(IN) :: NEL,G_WPLA
69 my_real,
DIMENSION(NEL*G_WPLA),
INTENT(IN) :: wpla
70
71
72
73 INTEGER I, M,N, FLAG
74
76 . x1, x2, x3, v1, v2, v3, fac,
77 . vxa(mvsiz), vya(mvsiz) , vza(mvsiz),
78 . va2(mvsiz), xmas(mvsiz),
79 . off(mvsiz), ei(mvsiz) , ek(mvsiz),
80 . xm(mvsiz) , ym(mvsiz) , zm(mvsiz),
81 . xx(mvsiz) , yy(mvsiz) , zz(mvsiz),
82 . xy(mvsiz) , yz(mvsiz) , zx(mvsiz),
83 . xx2(mvsiz), yy2(mvsiz), zz2(mvsiz),
84 . xcg(mvsiz), ycg(mvsiz), zcg(mvsiz),
85 . xxm(mvsiz), yym(mvsiz), zzm(mvsiz),
86 . ixx(mvsiz), iyy(mvsiz), izz(mvsiz),
87 . ixy(mvsiz), iyz(mvsiz), izx(mvsiz),
88 . rei(mvsiz), rek(mvsiz), off_l(mvsiz)
89 type (sensors_),INTENT(INOUT) :: SENSORS
90
91 flag = iparg(80)
92
93
94
95 n=1
96 DO i=lft,llt
97 v1=vx(i,n)
98 v2=vy(i,n)
99 v3=vz(i,n)
100 vxa(i)=v1
101 vya(i)=v2
102 vza(i)=v3
103 va2(i)=(v1**2+v2**2+v3**2)
104 ENDDO
105
106 DO n=2,8
107 DO i=lft,llt
108 v1=vx(i,n)
109 v2=vy(i,n)
110 v3=vz(i,n)
111 vxa(i)=vxa(i)+v1
112 vya(i)=vya(i)+v2
113 vza(i)=vza(i)+v3
114 va2(i)=va2(i)+(v1**2+v2**2+v3**2)
115 ENDDO
116 ENDDO
117
118 fac = three/fourteen
119 DO i=lft,llt
120 vxa(i)=vxa(i)*fac
121 vya(i)=vya(i)*fac
122 vza(i)=vza(i)*fac
123 va2(i)=va2(i)*fac
124 ENDDO
125
126 DO n=9,20
127 DO i=lft,llt
128 v1=vx(i,n)
129 v2=vy(i,n)
130 v3=vz(i,n)
131 vxa(i)=vxa(i)+v1
132 vya(i)=vya(i)+v2
133 vza(i)=vza(i)+v3
134 va2(i)=va2(i)+(v1**2+v2**2+v3**2)
135 ENDDO
136 ENDDO
137
138 fac = seven/eighty16
139 DO i=lft,llt
140 vxa(i)=vxa(i)*fac
141 vya(i)=vya(i)*fac
142 vza(i)=vza(i)*fac
143 va2(i)=va2(i)*fac
144 ENDDO
145
146 DO i=lft,llt
147 xmas(i)= fill(i)*rhog(i)*volg(i)
148 ei(i) = fill(i)*eintg(i)*vol0g(i)
149 ek(i) = xmas(i)*va2(i)*half
150 xm(i) = xmas(i)*vxa(i)
151 ym(i) = xmas(i)*vya(i)
152 zm(i) = xmas(i)*vza(i)
153 ENDDO
154
155 DO i=lft,llt
156 m=iparts(i)
157 partsav(1,m)=partsav(1,m) + fill(i)*eintg(i)*vol0g(i)
158 partsav(2,m)=partsav(2,m) + xmas(i)*va2(i)*half
159 partsav(3,m)=partsav(3,m) + xmas(i)*vxa(i)
160 partsav(4,m)=partsav(4,m) + xmas(i)*vya(i)
161 partsav(5,m)=partsav(5,m) + xmas(i)*vza(i)
162 IF (offg(i) >= one) partsav(6,m)=partsav(6,m) + xmas(i)
163 IF (g_wpla > 0) partsav(29,m)=partsav(29,m) + wpla(i)
164 ENDDO
165
166 IF(iexpan/=0)THEN
167 DO i=lft,llt
168 m=iparts(i)
169 partsav(27,m)=partsav(27,m) + eintth(i)*vol0g(i)
170 ENDDO
171 ENDIF
172
173
174
175
176 IF(flag==1) THEN
177 n=1
178 DO i=lft,llt
179 x1=x(i,n)
180 x2=y(i,n)
181 x3=z(i,n)
182 xx(i)=x1
183 yy(i)=x2
184 zz(i)=x3
185 xx2(i)=(x1**2)
186 yy2(i)=(x2**2)
187 zz2(i)=(x3**2)
188 xy(i)=(x1*x2)
189 yz(i)=(x2*x3)
190 zx(i)=(x3*x1)
191 ENDDO
192 DO n=2,8
193 DO i=lft,llt
194 x1=x(i,n)
195 x2=y(i,n)
196 x3=z(i,n)
197 xx(i)=xx(i)+x1
198 yy(i)=yy(i)+x2
199 zz(i)=zz(i)+x3
200 xx2(i)=xx2(i)+(x1**2)
201 yy2(i)=yy2(i)+(x2**2)
202 zz2(i)=zz2(i)+(x3**2)
203 xy(i)=xy(i)+(x1*x2)
204 yz(i)=yz(i)+(x2*x3)
205 zx(i)=zx(i)+(x3*x1)
206 ENDDO
207 ENDDO
208
209 fac = three/fourteen
210 DO i=lft,llt
211 xx(i)=xx(i)*fac
212 yy(i)=yy(i)*fac
213 zz(i)=zz(i)*fac
214 xx2(i)=xx2(i)*fac
215 yy2(i)=yy2(i)*fac
216 zz2(i)=zz2(i)*fac
217 xy(i)=xy(i)*fac
218 yz(i)=yz(i)*fac
219 zx(i)=zx(i)*fac
220 ENDDO
221
222 DO n=9,20
223 DO i=lft,llt
224 x1=x(i,n)
225 x2=y(i,n)
226 x3=z(i,n)
227 xx(i)=xx(i)+x1
228 yy(i)=yy(i)+x2
229 zz(i)=zz(i)+x3
230 xx2(i)=xx2(i)+(x1**2)
231 yy2(i)=yy2(i)+(x2**2)
232 zz2(i)=zz2(i)+(x3**2)
233 xy(i)=xy(i)+(x1*x2)
234 yz(i)=yz(i)+(x2*x3)
235 zx(i)=zx(i)+(x3*x1)
236 ENDDO
237 ENDDO
238
239 fac = seven/eighty16
240 DO i=lft,llt
241 xx(i)=xx(i)*fac
242 yy(i)=yy(i)*fac
243 zz(i)=zz(i)*fac
244 xx2(i)=xx2(i)*fac
245 yy2(i)=yy2(i)*fac
246 zz2(i)=zz2(i)*fac
247 xy(i)=xy(i)*fac
248 yz(i)=yz(i)*fac
249 zx(i)=zx(i)*fac
250 ENDDO
251
252 DO i=lft,llt
253 xcg(i)= xmas(i)*xx(i)
254 ycg(i)= xmas(i)*yy(i)
255 zcg(i)= xmas(i)*zz(i)
256 ixy(i)= -xmas(i)*xy(i)
257 iyz(i)= -xmas(i)*yz(i)
258 izx(i)= -xmas(i)*zx(i)
259 ixx(i)= xmas(i)*(yy2(i) + zz2(i))
260 iyy(i)= xmas(i)*(zz2(i) + xx2(i))
261 izz(i)= xmas(i)*(xx2(i) + yy2(i))
262 xxm(i)= vza(i)*ycg(i)-vya(i)*zcg(i)
263 yym(i)= vxa(i)*zcg(i)-vza(i)*xcg(i)
264 zzm(i)= vya(i)*xcg(i)-vxa(i)*ycg(i)
265 rei(i)= zero
266 rek(i)= zero
267 ENDDO
268
269 DO i=lft,llt
270 m=iparts(i)
271 partsav(9,m) =partsav(9,m) + xcg(i)
272 partsav(10,m)=partsav(10,m) + ycg(i)
273 partsav(11,m)=partsav(11,m) + zcg(i)
274 partsav(12,m)=partsav(12,m) + xxm(i)
275 partsav(13,m)=partsav(13,m) + yym(i)
276 partsav(14,m)=partsav(14,m) + zzm(i)
277 partsav(15,m)=partsav(15,m) + ixx(i)
278 partsav(16,m)=partsav(16,m) + iyy(i)
279 partsav(17,m)=partsav(17,m) + izz(i)
280 partsav(18,m)=partsav(18,m) + ixy(i)
281 partsav(19,m)=partsav(19,m) + iyz(i)
282 partsav(20,m)=partsav(20,m) + izx(i)
283 partsav(21,m)=partsav(21,m) + rei(i)
284 partsav(22,m)=partsav(22,m) + rek(i)
285
286 ENDDO
287 ENDIF
288
289
290 IF (igre /= 0) THEN
291 DO i=lft,llt
292 off(i) = one
293 ENDDO
294 CALL grelem_sav(lft ,llt ,gresav,igrth ,grth ,
295 2 off ,ei ,ek ,xm ,ym ,
296 3 zm ,xmas ,xcg ,ycg ,zcg ,
297 4 xxm ,yym ,zzm ,ixx ,iyy ,
298 5 izz ,ixy ,iyz ,izx ,rei ,
299 6 rek ,flag )
300 ENDIF
301
302
303 off_l(lft:llt) = one
304
305 DO i = lft,llt
306 m=iparts(i)
307 IF (offg(i) < one) THEN
308 partsav(25,m) = partsav(25,m) + one
309 ENDIF
310 ENDDO
311
313
314 RETURN
subroutine grelem_sav(jft, jlt, gresav, igrth, grth, off, ei, ek, xm, ym, zm, xmas, xcg, ycg, zcg, xxm, yym, zzm, ixx, iyy, izz, ixy, iyz, izx, rei, rek, flag)
subroutine sensor_energy_bilan(jft, jlt, ei, ek, off, ipart, itask, sensors)