47
48
49
50 USE sensor_mod
51
52
53
54#include "implicit_f.inc"
55
56
57
58#include "mvsiz_p.inc"
59
60
61
62#include "param_c.inc"
63#include "com_xfem1.inc"
64
65
66
67 INTEGER, INTENT(IN) :: IGRE
68 INTEGER JFT, JLT, IEXPAN,
69 . IPARTTG(*),IXTG(NIXTG,*),GRTH(*),(*),IXFEM,
70 . ILEV,IEL_CRK(*),IADTG_CRK(3,*),NFT1,ITASK
71 INTEGER, INTENT(IN) :: MAT(MVSIZ),ACTIFXFEM
72
74 . pm(npropm,*), v(3,*), thk(*), eint(jlt,2),partsav(npsav,*)
76 . rho(mvsiz), vol00(mvsiz),x(3,*),vr(3,*) ,thk02(*),
area(*),
77 . gresav(*), off(*),eintth(*)
78 my_real,
INTENT(IN) :: gvol(mvsiz)
79 type (sensors_),INTENT(INOUT) :: SENSORS
80 INTEGER, INTENT(IN) :: NEL,G_WPLA
81 my_real,
DIMENSION(NEL*G_WPLA),
INTENT(IN) :: wpla
82
83
84
85 INTEGER I, MX, II, J, IC, JST(MVSIZ+1),FLAG
86
88 . in25,xx,yy,zz,xy,yz,zx, coef, inel,
89 . vx1(mvsiz), vx2(mvsiz), vx3(mvsiz),
90 . vy1(mvsiz), vy2(mvsiz), vy3(mvsiz),
91 . vz1(mvsiz), vz2(mvsiz), vz3(mvsiz),
92 . ei(mvsiz),xmas(mvsiz),rei(mvsiz),rek(mvsiz),
93 . ek(mvsiz), xm(mvsiz), ym(mvsiz), zm(mvsiz),
94 . vxa(mvsiz), vya(mvsiz), vza(mvsiz) , va2, xmas25(mvsiz),
95 . xxm(mvsiz), yym(mvsiz), zzm(mvsiz),
96 . xcg(mvsiz), ycg(mvsiz), zcg(mvsiz),
97 . ixx(mvsiz), iyy(mvsiz), izz(mvsiz),
98 . ixy(mvsiz), iyz(mvsiz), izx(mvsiz),
99 . vrx1(mvsiz), vrx2(mvsiz), vrx3(mvsiz),
100 . vry1(mvsiz), vry2(mvsiz), vry3(mvsiz),
101 . vrz1(mvsiz), vrz2(mvsiz), vrz3(mvsiz),
102 . x1(mvsiz), x2(mvsiz), x3(mvsiz),
103 . y1(mvsiz), y2(mvsiz), y3(mvsiz),
104 . z1(mvsiz), z2(mvsiz), z3(mvsiz),
105 . vl1(mvsiz,3),vl2(mvsiz,3),vl3(mvsiz,3),
106 . vrl1(mvsiz,3),vrl2(mvsiz,3),vrl3(mvsiz,3)
107
108
109 IF(ixfem == 0)THEN
110 DO i=jft,jlt
111 vx1(i)=v(1,ixtg(2,i))
112 vy1(i)=v(2,ixtg(2,i))
113 vz1(i)=v(3,ixtg(2,i))
114 vx2(i)=v(1,ixtg(3,i))
115 vy2(i)=v(2,ixtg(3,i))
116 vz2(i)=v(3,ixtg(3,i))
117 vx3(i)=v(1,ixtg(4,i))
118 vy3(i)=v(2,ixtg(4,i))
119 vz3(i)=v(3,ixtg(4,i))
120 ENDDO
121 ELSEIF(ixfem==1)THEN
123 1 jft ,jlt ,nft1 ,ilev ,iel_crk,
124 2 x ,v ,vr ,vl1 ,vl2 ,
125 3 vl3 ,vrl1 ,vrl2 ,vrl3 ,x1 ,
126 4 x2 ,x3 ,y1 ,y2 ,y3 ,
127 5 z1 ,z2 ,z3 ,iadtg_crk)
128 DO i=jft,jlt
129 vx1(i)=vl1(i,1)
130 vy1(i)=vl1(i,2)
131 vz1(i)=vl1(i,3)
132 vx2(i)=vl2(i,1)
133 vy2(i)=vl2(i,2)
134 vz2(i)=vl2(i,3)
135 vx3(i)=vl3(i,1)
136 vy3(i)=vl3(i,2)
137 vz3(i)=vl3(i,3)
138 vrx1(i)=vrl1(i,1)
139 vry1(i)=vrl1(i,2)
140 vrz1(i)=vrl1(i,3)
141 vrx2(i)=vrl2(i,1)
142 vry2(i)=vrl2(i,2)
143 vrz2(i)=vrl2(i,3)
144 vrx3(i)=vrl3(i,1)
145 vry3(i)=vrl3(i,2)
146 vrz3(i)=vrl3(i,3)
147 ENDDO
148 ENDIF
149
150 mx = mat(jft)
151 DO i=jft,jlt
152 vxa(i)=vx1(i)+vx2(i)+vx3(i)
153 vya(i)=vy1(i)+vy2(i)+vy3(i)
154 vza(i)=vz1(i)+vz2(i)+vz3(i)
155 xmas(i)=pm(1,mx)*gvol(i)
156 va2 =vx1(i)*vx1(i)+vx2(i)*vx2(i)+vx3(i)*vx3(i)
157 1 +vy1(i)*vy1(i)+vy2(i)*vy2(i)+vy3(i)*vy3(i)
158 2 +vz1(i)*vz1(i)+vz2(i)*vz2(i)+vz3(i)*vz3(i)
159 ei(i) = eint(i,1) + eint(i,2)
160 ek(i) = xmas(i)*va2*one_over_6
161 xmas25(i)= xmas(i)*third
162 xm(i) = xmas25(i)*vxa(i)
163 ym(i) = xmas25(i)*vya(i)
164 zm(i) = xmas25(i)*vza(i)
165 ENDDO
166
167
168
169
170 IF(npsav>=21)THEN
171 IF(ixfem == 0)THEN
172 DO i=jft,jlt
173 xx= x(1,ixtg(2,i))+x(1,ixtg(3,i))+x(1,ixtg(4,i))
174 yy= x(2,ixtg(2,i))+x(2,ixtg(3,i))+x(2,ixtg(4,i))
175 zz= x(3,ixtg(2,i))+x(3,ixtg(3,i))+x(3,ixtg(4,i))
176 xcg(i)= xmas25(i)*xx
177 ycg(i)= xmas25(i)*yy
178 zcg(i)= xmas25(i)*zz
179 in25 = xmas25(i)*(thk02(i)+
area(i))*one_over_12
180 inel = three*in25
181 xx = third*xx
182 yy = third*yy
183 zz = third*zz
184 ixy(i) = -xcg(i)*yy
185 iyz(i) = -ycg(i)*zz
186 izx(i) = -zcg(i)*xx
187 xx = xcg(i)*xx
188 yy = ycg(i)*yy
189 zz = zcg(i)*zz
190 ixx(i)= inel + yy + zz
191 iyy(i)= inel + zz + xx
192 izz(i)= inel + xx + yy
193 vxa(i)=third*vxa(i)
194 vya(i)=third*vya(i)
195 vza(i)=third*vza(i)
196 xxm(i)= vza(i)*ycg(i)-vya(i)*zcg(i)
197 . +in25*
198 . (vr(1,ixtg(2,i))+vr(1,ixtg(3,i))+vr(1,ixtg(4,i)))
199 yym(i)= vxa(i)*zcg(i)-vza(i)*xcg(i)
200 . +in25*
201 . (vr(2,ixtg(2,i))+vr(2,ixtg(3,i))+vr(2,ixtg(4,i)))
202 zzm(i)= vya(i)*xcg(i)-vxa(i)*ycg(i)
203 . +in25*
204 . (vr(3,ixtg(2,i))+vr(3,ixtg(3,i))+vr(3,ixtg(4,i)))
205 va2 =
206 . vr(1,ixtg(2,i))*vr(1,ixtg(2,i))+vr(1,ixtg(3,i))*vr(1,ixtg(3,i))
207 . +vr(1,ixtg(4,i))*vr(1,ixtg(4,i))
208 . +vr(2,ixtg(2,i))*vr(2,ixtg(2,i))+vr(2,ixtg(3,i))*vr(2,ixtg(3,i))
209 . +vr(2,ixtg(4,i))*vr(2,ixtg(4,i))
210 . +vr(3,ixtg(2,i))*vr(3,ixtg(2,i))+vr(3,ixtg(3,i))*vr(3,ixtg(3,i))
211 . +vr(3,ixtg(4,i))*vr(3,ixtg(4,i))
212 rei(i)= eint(i,2)
213 rek(i)= in25*va2*half
214 ENDDO
215 ELSEIF(ixfem==1)THEN
216 DO i=jft,jlt
217 xx= x1(i)+x2(i)+x3(i)
218 yy= y1(i)+y2(i)+y3(i)
219 zz= z1(i)+z2(i)+z3(i)
220 xcg(i)= xmas25(i)*xx
221 ycg(i)= xmas25(i)*yy
222 zcg(i)= xmas25(i)*zz
223 in25 = xmas25(i)*(thk02(i)+
area(i))*one_over_12
224 inel = three*in25
225 xx = third*xx
226 yy = third
227 zz = third*zz
228 ixy(i) = -xcg(i)*yy
229 iyz(i) = -ycg(i)*zz
230 izx(i) = -zcg(i)*xx
231 xx = xcg(i)*xx
232 yy = ycg(i)*yy
233 zz = zcg(i)*zz
234 ixx(i)= inel + yy + zz
235 iyy(i)= inel + zz + xx
236 izz(i)= inel + xx + yy
237 vxa(i)=third*vxa(i)
238 vya(i)=third*vya(i)
239 vza(i)=third*vza(i)
240 xxm(i)= vza(i)*ycg(i)-vya(i)*zcg(i)
241 . +in25*(vrx1(i)+vrx2(i)+vrx3(i))
242 yym(i)= vxa(i)*zcg(i)-vza(i)*xcg(i)
243 . +in25*(vry1(i)+vry2(i)+vry3(i))
244 zzm(i)= vya(i)*xcg(i)-vxa(i)*ycg(i)
245 . +in25*(vrz1(i)+vrz2(i)+vrz3(i))
246 va2 = vrx1(i)*vrx1(i)+vrx2(i)*vrx2(i)+vrx3(i)*vrx3(i)
247 . +vry1(i)*vry1(i)+vry2(i)*vry2(i)+vry3(i)*vry3(i)
248 . +vrz1(i)*vrz1(i)+vrz2(i)*vrz2(i)+vrz3(i)*vrz3(i)
249 rei(i)= eint(i,2)
250 rek(i)= in25*va2*half
251 ENDDO
252 ENDIF
253
254 flag = 1
255
256 IF (igre /= 0) THEN
257 CALL grelem_sav(jft ,jlt ,gresav,igrth ,grth ,
258 2 off ,ei ,ek ,xm ,ym ,
259 3 zm ,xmas ,xcg ,ycg ,zcg ,
260 4 xxm ,yym ,zzm ,ixx ,iyy ,
261 5 izz ,ixy ,iyz ,izx ,rei ,
262 6 rek ,flag)
263 ENDIF
264
265 ic=1
266 jst(ic)=jft
267 DO j=jft+1,jlt
268 IF (iparttg(j)/=iparttg(j-1)) THEN
269 ic=ic+1
270 jst(ic)=j
271 ENDIF
272 ENDDO
273 jst(ic+1)=jlt+1
274
275 IF (ic==1) THEN
276 mx = iparttg(jft)
277 DO i=jft,jlt
278 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
279 IF(off(i)/=zero)THEN
280 partsav(1,mx)=partsav(1,mx) + ei(i)
281 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
282 partsav(2,mx)=partsav(2,mx) + ek(i)
283 partsav(3,mx)=partsav(3,mx) + xm(i)
284 partsav(4,mx)=partsav(4,mx) + ym(i)
285 partsav(5,mx)=partsav(5,mx) + zm(i)
286 ENDIF
287 ELSE
288 partsav(1,mx)=partsav(1,mx) + ei(i)
289 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
290 partsav(2,mx)=partsav(2,mx) + ek(i)
291 partsav(3,mx)=partsav(3,mx) + xm(i)
292 partsav(4,mx)=partsav(4,mx) + ym(i)
293 partsav(5,mx)=partsav(5,mx) + zm(i)
294 ENDIF
295 IF(off(i)/=zero)THEN
296 partsav(6,mx)=partsav(6,mx) + xmas(i)
297 ENDIF
298 partsav(9,mx)=partsav(9,mx) + xcg(i)
299 partsav(10,mx)=partsav(10,mx) + ycg(i)
300 partsav(11,mx)=partsav(11,mx) + zcg(i)
301 partsav(12,mx)=partsav(12,mx) + xxm(i)
302 partsav(13,mx)=partsav(13,mx) + yym(i)
303 partsav(14,mx)=partsav(14,mx) + zzm(i)
304 partsav(15,mx)=partsav(15,mx) + ixx(i)
305 partsav(16,mx)=partsav(16,mx) + iyy(i)
306 partsav(17,mx)=partsav(17,mx) + izz(i)
307 partsav(18,mx)=partsav(18,mx) + ixy(i)
308 partsav(19,mx)=partsav(19,mx) + iyz(i)
309 partsav(20,mx)=partsav(20,mx) + izx(i)
310 partsav(21,mx)=partsav(21,mx) + rei(i)
311 partsav(22,mx)=partsav(22,mx) + rek(i)
312 ENDDO
313 ELSE
314 DO ii=1,ic
315 mx=iparttg(jst(ii))
316 IF (jst(ii+1)-jst(ii)>15) THEN
317 DO i=jst(ii),jst(ii+1)-1
318 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
319 IF(off(i)/=zero)THEN
320 partsav(1,mx)=partsav(1,mx) + ei(i)
321 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
322 partsav(2,mx)=partsav(2,mx) + ek(i)
323 partsav(3,mx)=partsav(3,mx) + xm(i)
324 partsav(4,mx)=partsav(4,mx) + ym(i)
325 partsav(5,mx)=partsav(5,mx) + zm(i)
326 ENDIF
327 ELSE
328 partsav(1,mx)=partsav(1,mx) + ei(i)
329 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
330 partsav(2,mx)=partsav(2,mx) + ek(i)
331 partsav(3,mx)=partsav(3,mx) + xm(i)
332 partsav(4,mx)=partsav(4,mx) + ym(i)
333 partsav(5,mx)=partsav(5,mx) + zm(i)
334 ENDIF
335 IF(off(i)/=zero)THEN
336 partsav(6,mx)=partsav(6,mx) + xmas(i)
337 ENDIF
338 partsav(9,mx)=partsav(9,mx) + xcg(i)
339 partsav(10,mx)=partsav(10,mx) + ycg(i)
340 partsav(11,mx)=partsav(11,mx) + zcg(i)
341 partsav(12,mx)=partsav(12,mx) + xxm(i)
342 partsav(13,mx)=partsav(13,mx) + yym(i)
343 partsav(14,mx)=partsav(14,mx) + zzm(i)
344 partsav(15,mx)=partsav(15,mx) + ixx(i)
345 partsav(16,mx)=partsav(16,mx) + iyy(i)
346 partsav(17,mx)=partsav(17,mx) + izz(i)
347 partsav(18,mx)=partsav(18,mx) + ixy(i)
348 partsav(19,mx)=partsav(19,mx) + iyz(i)
349 partsav(20,mx)=partsav(20,mx) + izx(i)
350 partsav(21,mx)=partsav(21,mx) + rei(i)
351 partsav(22,mx)=partsav(22,mx) + rek(i)
352 ENDDO
353 ELSE
354 DO i=jst(ii),jst(ii+1)-1
355 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
356 IF(off(i)/=zero)THEN
357 partsav(1,mx)=partsav(1,mx) + ei(i)
358 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
359 partsav(2,mx)=partsav(2,mx) + ek(i)
360 partsav(3,mx)=partsav(3,mx) + xm(i)
361 partsav(4,mx)=partsav(4,mx) + ym(i)
362 partsav(5,mx)=partsav(5,mx) + zm(i)
363 ENDIF
364 ELSE
365 partsav(1,mx)=partsav(1,mx) + ei(i)
366 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
367 partsav(2,mx)=partsav(2,mx) + ek(i)
368 partsav(3,mx)=partsav(3,mx) + xm(i)
369 partsav(4,mx)=partsav(4,mx) + ym(i)
370 partsav(5,mx)=partsav(5,mx) + zm(i)
371 ENDIF
372 IF(off(i)/=zero)THEN
373 partsav(6,mx)=partsav(6,mx) + xmas(i)
374 ENDIF
375 partsav(9,mx)=partsav(9,mx) + xcg(i)
376 partsav(10,mx)=partsav(10,mx) + ycg(i)
377 partsav(11,mx)=partsav(11,mx) + zcg(i)
378 partsav(12,mx)=partsav(12,mx) + xxm(i)
379 partsav(13,mx)=partsav(13,mx) + yym(i)
380 partsav(14,mx)=partsav(14,mx) + zzm(i)
381 partsav(15,mx)=partsav(15,mx) + ixx(i)
382 partsav(16,mx)=partsav(16,mx) + iyy(i)
383 partsav(17,mx)=partsav(17,mx) + izz(i)
384 partsav(18,mx)=partsav(18,mx) + ixy(i)
385 partsav(19,mx)=partsav(19,mx) + iyz(i)
386 partsav(20,mx)=partsav(20,mx) + izx(i)
387 partsav(21,mx)=partsav(21,mx) + rei(i)
388 partsav(22,mx)=partsav(22,mx) + rek(i)
389 ENDDO
390 ENDIF
391 ENDDO
392 ENDIF
393 ELSE
394
395 ic=1
396 jst(ic)=jft
397 DO j=jft+1,jlt
398 IF (iparttg(j)/=iparttg(j-1)) THEN
399 ic=ic+1
400 jst(ic)=j
401 ENDIF
402 ENDDO
403 jst(ic+1)=jlt+1
404
405 IF (ic==1) THEN
406 mx = iparttg(jft)
407 DO i=jft,jlt
408 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
409 IF(off(i)/=zero)THEN
410 partsav(1,mx)=partsav(1,mx) + ei(i)
411 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
412 partsav(2,mx)=partsav(2,mx) + ek(i)
413 partsav(3,mx)=partsav(3,mx) + xm(i)
414 partsav(4,mx)=partsav(4,mx) + ym(i)
415 partsav(5,mx)=partsav(5,mx) + zm(i)
416 ENDIF
417 ELSE
418 partsav(1,mx)=partsav(1,mx) + ei(i)
419 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
420 partsav(2,mx)=partsav(2,mx) + ek(i)
421 partsav(3,mx)=partsav(3,mx) + xm(i)
422 partsav(4,mx)=partsav(4,mx) + ym(i)
423 partsav(5,mx)=partsav(5,mx) + zm(i)
424 ENDIF
425 IF(off(i)/=zero)THEN
426 partsav(6,mx)=partsav(6,mx) + xmas(i)
427 ENDIF
428 ENDDO
429 ELSE
430 DO ii=1,ic
431 mx=iparttg(jst(ii))
432 IF (jst(ii+1)-jst(ii)>15) THEN
433 DO i=jst(ii),jst(ii+1)-1
434 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
435 IF(off(i)/=zero)THEN
436 partsav(1,mx)=partsav(1,mx) + ei(i)
437 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
438 partsav(2,mx)=partsav(2,mx) + ek(i)
439 partsav(3,mx)=partsav(3,mx) + xm(i)
440 partsav(4,mx)=partsav(4,mx) + ym(i)
441 partsav(5,mx)=partsav(5,mx) + zm(i)
442 ENDIF
443 ELSE
444 partsav(1,mx)=partsav(1,mx) + ei(i)
445 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
446 partsav(2,mx)=partsav(2,mx) + ek(i)
447 partsav(3,mx)=partsav(3,mx) + xm(i)
448 partsav(4,mx)=partsav(4,mx) + ym(i)
449 partsav(5,mx)=partsav(5,mx) + zm(i)
450 ENDIF
451 IF(off(i)/=zero)THEN
452 partsav(6,mx)=partsav(6,mx) + xmas(i)
453 ENDIF
454 ENDDO
455 ELSE
456 DO i=jst(ii),jst(ii+1)-1
457 IF (icrack3d > 0 .AND. ixfem > 0 .AND. actifxfem > 0) THEN
458 IF(off(i)/=zero)THEN
459 partsav(1,mx)=partsav(1,mx) + ei(i)
460 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
461 partsav(2,mx)=partsav(2,mx) + ek(i)
462 partsav(3,mx)=partsav(3,mx) + xm(i)
463 partsav(4,mx)=partsav(4,mx) + ym(i)
464 partsav(5,mx)=partsav(5,mx) + zm(i)
465 ENDIF
466 ELSE
467 partsav(1,mx)=partsav(1,mx) + ei(i)
468 IF (g_wpla > 0) partsav(29,mx)=partsav(29,mx) + wpla(i)
469 partsav(2,mx)=partsav(2,mx) + ek(i)
470 partsav(3,mx)=partsav(3,mx) + xm(i)
471 partsav(4,mx)=partsav(4,mx) + ym(i)
472 partsav(5,mx)=partsav(5,mx) + zm(i)
473 ENDIF
474 IF(off(i)/=zero)THEN
475 partsav(6,mx)=partsav(6,mx) + xmas(i)
476 ENDIF
477 ENDDO
478 ENDIF
479 ENDDO
480 ENDIF
481 ENDIF
482
483 IF(iexpan > 0) THEN
484 DO i=jft,jlt
485 mx = iparttg(i)
486 IF(off(i)/=0.0)THEN
487 partsav(27,mx)=partsav(27,mx) + eintth(i)
488 ENDIF
489 ENDDO
490 ENDIF
491
492 DO i = jft,jlt
493 mx = iparttg(i)
494 IF (off(i)==zero) THEN
495 partsav(25,mx) = partsav(25,mx) + one
496 ENDIF
497 ENDDO
498
500
501 RETURN
subroutine c3coor3_crk2(jft, jlt, nft, ilev, iel_crk, x, v, vr, vl1, vl2, vl3, vrl1, vrl2, vrl3, x1, x2, x3, y1, y2, y3, z1, z2, z3, iadtg_crk)
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 area(d1, x, x2, y, y2, eint, stif0)
subroutine sensor_energy_bilan(jft, jlt, ei, ek, off, ipart, itask, sensors)