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