43
44
45
46
47
48
49
50
51
52
53
56
57
58
59#include "implicit_f.inc"
60
61
62
63#include "param_c.inc"
64
65
66
67#include "units_c.inc"
68#include "com04_c.inc"
69#include "r2r_c.inc"
70#include "com01_c.inc"
71
72
73
74 INTEGER NPBY(NNPBY,*), M, ISPH, NRB
75 INTEGER ITAB(*), RBYID,ITAGND(*)
76 INTEGER, TARGET :: LPBY(*)
78 . rby(nrby,*), ms(*), in(*), x(3,*), skew(lskew,*),
79 . b1, b2, b3, b5, b6, b9,totmas ,xgt ,ygt ,
80 . zgt, stifn(*), stifr(*), v(3,*), vr(3,*),
81 . ii1,ii2,ii3,ii4,ii5,ii6,ii7,ii8,ii9,
82 . rby_r2r(9),rby_iniaxis(7,*)
83 INTEGER ID
84 CHARACTER(LEN=NCHARTITLE) :: TITR
85
86
87
88 INTEGER NSL, J, NOSKEW, I, N, NONOD,ICDG, ,NSL_XTRA
90 . xg(3), xm0(3), xmg, xx, xy, xz, yy, yz, zz, xiin, inmin,
91 . masrb,dd,tol,
92 . v1x2, v2x1, v2x3, v3x2, v3x1, v1x3,
93 . inmax,x_msn0(3),dist
94 INTEGER M_I, OPT_MERGE, RBLEVEL
95 INTEGER, DIMENSION(:), POINTER :: LSN, LSN_XTRA
96 my_real,
DIMENSION(:),
ALLOCATABLE :: ms_loc,in_loc
97
98 tol=one+em04
99
100 nsl = npby(2,nrb)
101 lsn => lpby(npby(11,nrb)+1:npby(11,nrb)+nsl)
102 rblevel = npby(12,nrb)
103 ALLOCATE (ms_loc(nsl))
104 ALLOCATE (in_loc(nsl))
105
106
107 flag = 0
108 DO i=1,nsl
109 n=lsn(i)
110 ms_loc(i)=ms(n)
111 in_loc(i)=in(n)
112 ENDDO
113
114 IF ((nsubdom>0).AND.(ipid==0)) THEN
115 IF (
tagno(m+n_part)==3)
THEN
116 flag = 1
117 DO j=1,9
118 rby(j,nrb) = zero
119 END DO
120
121 DO i=1,nsl
122 n=lsn(i)
123 IF (
tagno(n+n_part)>3) ms_loc(i)= 0
124 IF (
tagno(n+n_part)>3) in_loc(i)= 0
125 IF (
tagno(n+n_part)==0) ms_loc(i)= 0
126 IF (
tagno(n+n_part)==0) in_loc(i)= 0
127 ENDDO
128 ENDIF
129 ENDIF
130
131 IF (ns10e>0) THEN
132 DO i=1,nsl
133 n=lsn(i)
134 IF (itagnd(n)/=0) THEN
135 ms_loc(i)= zero
136 in_loc(i)= zero
137 END IF
138 ENDDO
139 END IF
140
141 nsl_xtra = npby(14,nrb)+npby(15,nrb)+npby(16,nrb)
142 nsl = npby(2,nrb) - nsl_xtra
143 lsn => lpby(npby(11,nrb)+1:npby(11,nrb)+nsl)
144
145
146 ms(m)=ms(m)+rby(1,nrb)
147 IF(ms(m)==zero) ms(m)=em20
148
149
150 nonod=itab(m)
151
152 icdg = npby(3,nrb)
153 noskew=npby(9,nrb)
154 rby(1,nrb)=rby(2,nrb)
155 rby(9,nrb)=rby(4,nrb)
156 rby(8,nrb)=rby(6,nrb)
157 rby(2,nrb)=rby(5,nrb)
158 rby(4,nrb)=rby(5,nrb)
159 rby(5,nrb)=rby(3,nrb)
160 rby(3,nrb)=rby(7,nrb)
161
162 IF(noskew/=0) THEN
163 CALL chbas(skew(1,noskew),rby(1,nrb))
164 ENDIF
165
166 rby(1,nrb)=rby(1,nrb)+in(m)
167 rby(5,nrb)=rby(5,nrb)+in(m)
168 rby(9,nrb)=rby(9,nrb)+in(m)
169 in(m) = (rby(1,nrb) + rby(5,nrb) + rby(9,nrb)) * third
170
171
172
173
174
175
176 xmg=ms(m)
177
178
179 IF (rby_iniaxis(1,nrb) > 0) THEN
180 DO j=1,3
181 x_msn0(j)=x(j,m)
182 ENDDO
183 ENDIF
184
185
186 IF(icdg==1)THEN
187 masrb=ms(m)
188 DO j=1,3
189 xg(j)=x(j,m)
190 x(j,m)=x(j,m)*ms(m)
191 ENDDO
192 DO i=1,nsl
193 n=lsn(i)
194 DO j=1,3
195 x(j,m) = x(j,m)+x(j,n)*ms_loc(i)
196 ENDDO
197 masrb = masrb+ms_loc(i)
198 ENDDO
199
200 IF(masrb<=1.e-30) THEN
202 . msgtype=msgerror,
203 . anmode=aninfo_blind_1,
205 . c1=titr)
206 RETURN
207 ENDIF
208
209 DO j=1,3
210 x(j,m)=x(j,m)/masrb
211 ENDDO
212
213
214 ELSEIF(icdg==2)THEN
215 masrb=zero
216 DO j=1,3
217 x(j,m)=zero
218 ENDDO
219 DO i=1,nsl
220 n=lsn(i)
221 DO j=1,3
222 x(j,m) = x(j,m)+x(j,n)*ms_loc(i)
223 ENDDO
224 masrb = masrb+ms_loc(i)
225 ENDDO
226 IF (flag==1) masrb =
max(masrb,em20)
227
228 IF(masrb<=em30) THEN
230 . msgtype=msgerror,
231 . anmode=aninfo_blind_1,
233 . c1=titr,
234 . c2'ON SECONDARY NODES')
235 RETURN
236 ENDIF
237
238 DO j=1,3
239 x(j,m)=x(j,m)/masrb
240 xg(j)=x(j,m)
241 ENDDO
242
243 masrb=masrb+ms(m)
244
245
246 ELSEIF(icdg==3)THEN
247 DO j=1,3
248 xg(j)=x(j,m)
249 ENDDO
250 masrb=ms(m)
251 DO i=1,nsl
252 n=lsn(i)
253 masrb = masrb+ms_loc(i)
254 ENDDO
255
256 IF(masrb<=em30) THEN
258 . msgtype=msgerror,
259 . anmode=aninfo_blind_1,
261 . c1=titr)
262 RETURN
263 ENDIF
264
265
266 ELSEIF(icdg==4)THEN
267 DO j=1,3
268 xg(j)=x(j,m)
269 ENDDO
270 masrb=ms(m)
271
272 IF(masrb<=em30) THEN
274 . msgtype=msgerror,
275 . anmode=aninfo_blind_1,
277 . c1=titr,
278 . c2='ON MAIN NODE')
279 RETURN
280 ENDIF
281
282 ENDIF
283
284
285
286
287 IF(npby(15,nrb) > 0) THEN
288 lsn_xtra => lpby(npby(11,nrb)+nsl+npby(14,nrb)+1:
289 . npby(11,nrb)+nsl+npby(14,nrb)+npby(15,nrb))
290 DO j=1,3
291 xg(j)=x(j,m)
292 x(j,m)=x(j,m)*masrb
293 ENDDO
294 xmg=masrb
295
296 DO i=1,npby(15,nrb)
297 n=lsn_xtra(i)
298 DO j=1,3
299 x(j,m) = x(j,m)+x(j,n)*ms_loc(nsl+i)
300 ENDDO
301 masrb = masrb+ms_loc(nsl+i)
302 ENDDO
303
304 DO j=1,3
305 x(j,m)=x(j,m)/masrb
306 xg(j)=x(j,m)
307 ENDDO
308 ENDIF
309
310 IF(npby(16,nrb) > 0) THEN
311 lsn_xtra => lpby(npby(11,nrb)+nsl+npby(14,nrb)+npby(15,nrb)+1:
312 . npby(11,nrb)+nsl+npby(14,nrb)+npby(15,nrb)+npby(16,nrb))
313 DO i=1,npby(16,nrb)
314 n=lsn_xtra(i)
315 masrb = masrb+ms_loc(nsl+i)
316 ENDDO
317 ENDIF
318
319
320
321
322 IF(icdg<=3)THEN
323 IF(n2d==0)THEN
324
325 xx=(xg(1)-x(1,m))*(xg(1)-x(1,m))
326 xy=(xg(1)-x(1,m))*(xg(2)-x(2,m))
327 xz=(xg(1)-x(1,m))*(xg(3)-x(3,m))
328 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
329 yz=(xg(2)-x(2,m))*(xg(3)-x(3,m))
330 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
331 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
332 rby(2,nrb)=rby(2,nrb)-xy*xmg
333 rby(3,nrb)=rby(3,nrb)-xz*xmg
334 rby(4,nrb)=rby(4,nrb)-xy*xmg
335 rby(5,nrb)=rby(5,nrb)+(zz+xx)*xmg
336 rby(6,nrb)=rby(6,nrb)-yz*xmg
337 rby(7,nrb)=rby(7,nrb)-xz*xmg
338 rby(8,nrb)=rby(8,nrb)-yz*xmg
339 rby(9,nrb)=rby(9,nrb)+(xx+yy)*xmg
340
341 IF (nsl==1) THEN
342 rby(1,nrb)=rby(1,nrb)+em20
343 rby(5,nrb)=rby(5,nrb)+em20
344 rby(9,nrb)=rby(9,nrb)+em20
345 ENDIF
346
347 DO i=1,nsl
348 n=lsn(i)
349 xx=(x(1,n)-x(1,m))*(x(1,n)-x(1,m))
350 xy=(x(1,n)-x(1,m))*(x(2,n)-x(2,m))
351 xz=(x(1,n)-x(1,m))*(x(3,n)-x(3,m))
352 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
353 yz=(x(2,n)-x(2,m))*(x(3,n)-x(3,m))
354 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
355 rby(1,nrb)=rby(1,nrb)+in_loc(i)+(yy+zz)*ms_loc(i)
356 rby(2,nrb)=rby(2,nrb)-xy*ms_loc(i)
357 rby(3,nrb)=rby(3,nrb)-xz*ms_loc(i)
358 rby(4,nrb)=rby(4,nrb)-xy*ms_loc(i)
359 rby(5,nrb)=rby(5,nrb)+in_loc(i)+(zz+xx)*ms_loc(i)
360 rby(6,nrb)=rby(6,nrb)-yz*ms_loc(i)
361 rby(7,nrb)=rby(7,nrb)-xz*ms_loc(i)
362 rby(8,nrb)=rby(8,nrb)-yz*ms_loc(i)
363 rby(9,nrb)=rby(9,nrb)+in_loc(i)+(xx+yy)*ms_loc(i)
364 ENDDO
365
366 ELSEIF(n2d==1) THEN
367
368
369
370
371 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
372 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
373 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
374 rby(2,nrb)=zero
375 rby(3,nrb)=zero
376
377 rby(4,nrb)=zero
378 rby(5,nrb)=rby(5,nrb)+zz*xmg
379 rby(6,nrb)=zero
380
381 rby(7,nrb)=zero
382 rby(8,nrb)=zero
383 rby(9,nrb)=rby(9,nrb)+yy*xmg
384
385 IF (nsl==1) THEN
386 rby(1,nrb)=rby(1,nrb)+em20
387 rby(5,nrb)=rby(5,nrb)+em20
388 rby(9,nrb)=rby(9,nrb)+em20
389 ENDIF
390
391 DO i=1,nsl
392 n=lsn(i)
393 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
394 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
395 rby(1,nrb)=rby(1,nrb)+in_loc(i)+(yy+zz)*ms_loc(i)
396 rby(5,nrb)=rby(5,nrb)+in_loc(i)+zz*ms_loc(i)
397 rby(9,nrb)=rby(9,nrb)+in_loc(i)+yy*ms_loc(i)
398 ENDDO
399 ELSEIF(n2d==2) THEN
400
401
402
403
404 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
405 yz=(xg(2)-x(2,m))*(xg(3)-x(3,m))
406 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
407 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
408 rby(2,nrb)=zero
409 rby(3,nrb)=zero
410
411 rby(4,nrb)=zero
412 rby(5,nrb)=rby(5,nrb)+zz*xmg
413 rby(6,nrb)=rby(6,nrb)-yz*xmg
414
415 rby(7,nrb)=zero
416 rby(8,nrb)=rby(8,nrb)-yz*xmg
417 rby(9,nrb)=rby(9,nrb)+yy*xmg
418
419 IF (nsl==1) THEN
420 rby(1,nrb)=rby(1,nrb)+em20
421 rby(5,nrb)=rby(5,nrb)+em20
422 rby(9,nrb)=rby(9,nrb)+em20
423 ENDIF
424
425 DO i=1,nsl
426 n=lsn(i)
427 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
428 yz=(x(2,n)-x(2,m))*(x(3,n)-x(3,m))
429 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
430 rby(1,nrb)=rby(1,nrb)+in_loc(i)+(yy+zz)*ms_loc(i)
431
432 rby(5,nrb)=rby(5,nrb)+in_loc(i)+zz*ms_loc(i)
433 rby(6,nrb)=rby(6,nrb)-yz*ms_loc(i)
434
435 rby(8,nrb)=rby(8,nrb)-yz*ms_loc(i)
436 rby(9,nrb)=rby(9,nrb)+in_loc(i)+yy*ms_loc(i)
437 ENDDO
438 ENDIF
439 ENDIF
440
441
442
443
444 IF((npby(15,nrb)+npby(16,nrb)) > 0) THEN
445
446 lsn_xtra => lpby(npby(11,nrb)+nsl+npby(14,nrb)+1:
447 . npby(11,nrb)+nsl+npby(14,nrb)+npby(15,nrb)+npby(16,nrb))
448
449 IF(n2d==0)THEN
450 DO i=1,npby(15,nrb)+npby(16,nrb)
451 n=lsn_xtra(i)
452 xx=(x(1,n)-x(1,m))*(x(1,n)-x(1,m))
453 xy=(x(1,n)-x(1,m))*(x(2,n)-x(2,m))
454 xz=(x(1,n)-x(1,m))*(x(3,n)-x(3,m))
455 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
456 yz=(x(2,n)-x(2,m))*(x(3,n)-x(3,m))
457 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
458 rby(1,nrb)=rby(1,nrb)+in_loc(nsl+i)+(yy+zz)*ms_loc(nsl+i)
459 rby(2,nrb)=rby(2,nrb)-xy*ms_loc(nsl+i)
460 rby(3,nrb)=rby(3,nrb)-xz*ms_loc(nsl+i)
461 rby(4,nrb)=rby(4,nrb)-xy*ms_loc(nsl+i)
462 rby(5,nrb)=rby(5,nrb)+in_loc(nsl+i)+(zz+xx)*ms_loc(nsl+i)
463 rby(6,nrb)=rby(6,nrb)-yz*ms_loc(nsl+i)
464 rby(7,nrb)=rby(7,nrb)-xz*ms_loc(nsl+i)
465 rby(8,nrb)=rby(8,nrb)-yz*ms_loc(nsl+i)
466 rby(9,nrb)=rby(9,nrb)+in_loc(nsl+i)+(xx+yy)*ms_loc(nsl+i)
467 ENDDO
468 ELSEIF(n2d==1) THEN
469 DO i=1,npby(15,nrb)+npby(16,nrb)
470 n=lsn_xtra(i)
471 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
472 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
473 rby(1,nrb)=rby(1,nrb)+in_loc(nsl+i)+(yy+zz)*ms_loc(nsl+i)
474 rby(5,nrb)=rby(5,nrb)+in_loc(nsl+i)+zz*ms_loc(nsl+i)
475 rby(9,nrb)=rby(9,nrb)+in_loc(nsl+i)+yy*ms_loc(nsl+i)
476 ENDDO
477 ELSEIF(n2d==1) THEN
478 DO i=1,npby(15,nrb)+npby(16,nrb)
479 n=lsn_xtra(i)
480 yy=(x(2,n)-x(2,m))*(x(2,n)-x(2,m))
481 yz=(x(2,n)-x(2,m))*(x(3,n)-x(3,m))
482 zz=(x(3,n)-x(3,m))*(x(3,n)-x(3,m))
483 rby(1,nrb)=rby(1,nrb)+in_loc(nsl+i)+(yy+zz)*ms_loc(nsl+i)
484 rby(5,nrb)=rby(5,nrb)+in_loc(nsl+i)+zz*ms_loc(nsl+i)
485 rby(6,nrb)=rby(6,nrb)-yz*ms_loc(nsl+i)
486 rby(8,nrb)=rby(8,nrb)-yz*ms_loc(nsl+i)
487 rby(9,nrb)=rby(9,nrb)+in_loc(nsl+i)+yy*ms_loc(nsl+i)
488 ENDDO
489 ENDIF
490 ENDIF
491
492
493
494
495 DO i=nrb-1,1,-1
496 IF((npby(12,i)) == rblevel-1) THEN
497 DO j=1,3
498 xg(j)=x(j,m)
499 ENDDO
500 xmg=masrb
501
502 opt_merge = npby(13,i)
503 npby(11,nrb) = npby(11,i)
504 npby(2,nrb) = npby(2,nrb)+npby(2,i)
505 m_i=npby(1,i)
506
507 IF(opt_merge == 2) THEN
508 DO j=1,3
509 x(j,m)=x(j,m)*masrb+x(j,m_i)*rby(14,i)
510 ENDDO
511 masrb = masrb + rby(14,i)
512 DO j=1,3
513 x(j,m)=x(j,m)/masrb
514 ENDDO
515 ELSEIF(opt_merge == 3) THEN
516 masrb = masrb + rby(14,i)
517 ENDIF
518
519
520
521
522 IF(opt_merge == 2) THEN
523 IF(n2d==0)THEN
524 xx=(xg(1)-x(1,m))*(xg(1)-x(1,m))
525 xy=(xg(1)-x(1,m))*(xg(2)-x(2,m))
526 xz=(xg(1)-x(1,m))*(xg(3)-x(3,m))
527 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
528 yz=(xg(2)-x(2,m))*(xg(3)-x(3,m))
529 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
530 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
531 rby(2,nrb)=rby(2,nrb)-xy*xmg
532 rby(3,nrb)=rby(3,nrb)-xz*xmg
533 rby(4,nrb)=rby(4,nrb)-xy*xmg
534 rby(5,nrb)=rby(5,nrb)+(zz+xx)*xmg
535 rby(6,nrb)=rby(6,nrb)-yz*xmg
536 rby(7,nrb)=rby(7,nrb)-xz*xmg
537 rby(8,nrb)=rby(8,nrb)-yz*xmg
538 rby(9,nrb)=rby(9,nrb)+(xx+yy)*xmg
539 ELSEIF(n2d==1) THEN
540 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
541 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
542 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
543 rby(2,nrb)=zero
544 rby(3,nrb)=zero
545 rby(4,nrb)=zero
546 rby(5,nrb)=rby(5,nrb)+zz*xmg
547 rby(6,nrb)=zero
548 rby(7,nrb)=zero
549 rby(8,nrb)=zero
550 rby(9,nrb)=rby(9,nrb)+yy*xmg
551 ELSEIF(n2d==2) THEN
552 yy=(xg(2)-x(2,m))*(xg(2)-x(2,m))
553 yz=(xg(2)-x(2,m))*(xg(3)-x(3,m))
554 zz=(xg(3)-x(3,m))*(xg(3)-x(3,m))
555 rby(1,nrb)=rby(1,nrb)+(yy+zz)*xmg
556 rby(2,nrb)=zero
557 rby(3,nrb)=zero
558 rby(4,nrb)=zero
559 rby(5,nrb)=rby(5,nrb)+zz*xmg
560 rby(6,nrb)=rby(6,nrb)-yz*xmg
561 rby(7,nrb)=zero
562 rby(8,nrb)=rby(8,nrb)-yz*xmg
563 rby(9,nrb)=rby(9,nrb)+yy*xmg
564 ENDIF
565 ENDIF
566 IF(opt_merge > 1) THEN
567 IF(n2d==0)THEN
568 xx=(x(1,m_i)-x(1,m))*(x(1,m_i)-x(1,m))
569 xy=(x(1,m_i)-x(1,m))*(x(2,m_i)-x(2,m))
570 xz=(x(1,m_i)-x(1,m))*(x(3,m_i)-x(3,m))
571 yy=(x(2,m_i)-x(2,m))*(x(2,m_i)-x(2,m))
572 yz=(x(2,m_i)-x(2,m))*(x(3,m_i)-x(3,m))
573 zz=(x(3,m_i)-x(3,m))*(x(3,m_i)-x(3,m))
574 rby(1,nrb)=rby(1,nrb)+rby(1,i)+(yy+zz)*rby(14,i)
575 rby(2,nrb)=rby(2,nrb)+rby(2,i)-xy*rby(14,i)
576 rby(3,nrb)=rby(3,nrb)+rby(3,i)-xz*rby(14,i)
577 rby(4,nrb)=rby(4,nrb)+rby(4,i)-xy*rby(14,i)
578 rby(5,nrb)=rby(5,nrb)+rby(5,i)+(zz+xx)*rby(14,i)
579 rby(6,nrb)=rby(6,nrb)+rby(6,i)-yz*rby(14,i)
580 rby(7,nrb)=rby(7,nrb)+rby(7,i)-xz*rby(14,i)
581 rby(8,nrb)=rby(8,nrb)+rby(8,i)-yz*rby(14,i)
582 rby(9,nrb)=rby(9,nrb)+rby(9,i)+(xx+yy)*rby(14,i)
583 ELSEIF(n2d==1) THEN
584 yy=(x(2,m_i)-x(2,m))*(x(2,m_i)-x(2,m))
585 zz=(x(3,m_i)-x(3,m))*(x(3,m_i)-x(3,m))
586 rby(1,nrb)=rby(1,nrb)+rby(1,i)+(yy+zz)*rby(14,i)
587 rby(5,nrb)=rby(5,nrb)+rby(5,i)+zz*rby(14,i)
588 rby(9,nrb)=rby(9,nrb)+rby(9,i)+yy*rby(14,i)
589 ELSEIF(n2d==2) THEN
590 yy=(x(2,m_i)-x(2,m))*(x(2,m_i)-x(2,m))
591 yz=(x(2,m_i)-x(2,m))*(x(3,m_i)-x(3,m))
592 zz=(x(3,m_i)-x(3,m))*(x(3,m_i)-x(3,m))
593 rby(1,nrb)=rby(1,nrb)+rby(1,i)+(yy+zz)*rby(14,i)
594 rby(5,nrb)=rby(5,nrb)+rby(5,i)+zz*rby(14,i)
595 rby(6,nrb)=rby(6,nrb)-yz*rby(14,i)
596 rby(8,nrb)=rby(8,nrb)-yz*rby(14,i)
597 rby(9,nrb)=rby(9,nrb)+rby(9,i)+yy*rby(14,i)
598 ENDIF
599 ENDIF
600
601 rby(1:15,i)=zero
602 in(m_i)=zero
603 ms(m_i)=zero
604 ELSEIF(npby(12,i) == rblevel) THEN
605 EXIT
606 ENDIF
607 ENDDO
608
609
610
611
612
613 IF(rblevel == 0) THEN
614 IF(nrbmerge > 0) THEN
615 DEALLOCATE(ms_loc)
616 DEALLOCATE(in_loc)
617 nsl = npby(2,nrb)
618 lsn => lpby(npby(11,nrb)+1:npby(11,nrb)+nsl)
619 ALLOCATE (ms_loc(nsl))
620 ALLOCATE (in_loc(nsl))
621
622
623 IF ((nsubdom>0).AND.(ipid==0)) THEN
624 IF (
tagno(m+n_part)==3)
THEN
625
626 DO i=1,nsl
627 n=lsn(i)
628 IF (
tagno(n+n_part)>3) ms_loc(i)= 0
629 IF (
tagno(n+n_part)>3) in_loc(i)= 0
630 IF (
tagno(n+n_part)==0) ms_loc(i)= 0
631 IF (
tagno(n+n_part)==0) in_loc(i)= 0
632 ENDDO
633 ENDIF
634 ENDIF
635
636 IF (ns10e>0) THEN
637 DO i=1,nsl
638 n=lsn(i)
639 IF (itagnd(n)/=0) THEN
640 ms_loc(i)= zero
641 in_loc(i)= zero
642 ELSE
643 ms_loc(i)=ms(n)
644 in_loc(i)=in(n)
645 ENDIF
646 ENDDO
647 ELSE
648 DO i=1,nsl
649 n=lsn(i)
650 ms_loc(i)=ms(n)
651 in_loc(i)=in(n)
652 ENDDO
653 ENDIF
654 ENDIF
655
656
657 IF(nsubdom>0)THEN
658 IF(
tagno(m+n_part)==3)
THEN
659 DO j=1,9
660 rby_r2r(j)=rby(j,nrb)
661 END DO
662 END IF
663 ENDIF
664
665
666
667
668
669
670
671
672
673
674
675
676 WRITE(iout,1000)
677 DO i=1,nsl
678 n=lsn(i)
679 xx=(x(1,n))*(x(1,n))
680 xy=(x(1,n))*(x(2,n))
681 xz=(x(1,n))*(x(3,n))
682 yy=(x(2,n))*(x(2,n))
683 yz=(x(2,n))*(x(3,n))
684 zz=(x(3,n))*(x(3,n))
685 b1 = b1 - in_loc(i)-(yy+zz)*ms_loc(i)
686 b2 = b2 + xy*ms_loc(i)
687 b3 = b3 + xz*ms_loc(i)
688 b5 = b5 - in_loc(i)-(zz+xx)*ms_loc(i)
689 b6 = b6 + yz*ms_loc(i)
690 b9 = b9 - in_loc(i)-(xx+yy)*ms_loc(i)
691 xgt = xgt - ms_loc(i)*x(1,n)
692 ygt = ygt - ms_loc(i)*x(2,n)
693 zgt = zgt - ms_loc(i)*x(3,n)
694 totmas = totmas - ms_loc(i)
695 ENDDO
696 IF(isph/=1)THEN
697 b1 = b1 + rby(1,nrb)
698 b2 = b2 + rby(2,nrb)
699 b3 = b3 + rby(3,nrb)
700 b5 = b5 + rby(5,nrb)
701 b6 = b6 + rby(6,nrb)
702 b9 = b9 + rby(9,nrb)
703 IF(nsubdom>0)THEN
704 IF(
tagno(m+n_part)==3)
THEN
705 WRITE(iout,1300) rbyid
706 WRITE(iout,1400)
707 ELSE
708 WRITE(iout,1100) rbyid,nonod,x(1,m),x(2,m),x(3,m),
709 . masrb,rby(1,nrb),rby(5,nrb),rby(9,nrb),rby(2,nrb),rby(6,nrb),rby(3,nrb)
710 END IF
711 ELSE
712 WRITE(iout,1100) rbyid,nonod,x(1,m),x(2,m),x(3,m),
713 . masrb,rby(1,nrb),rby(5,nrb),rby(9,nrb),rby(2,nrb),rby(6,nrb),rby(3,nrb)
714 END IF
715 ENDIF
716
717
718
719
720 IF(n2d == 1) THEN
721 rby(10,nrb) = rby(1,nrb)
722 rby(11,nrb) = rby(5,nrb)
723 rby(12,nrb) = rby(9,nrb)
724 rby(1,nrb) = one
725 rby(5,nrb) = one
726 rby(9,nrb) = one
727 ELSE
728 CALL inepri(rby(10,nrb),rby(1,nrb))
729 ENDIF
730
731 IF(isph==1)THEN
732 xiin = (rby(10,nrb) + rby(11,nrb) + rby(12,nrb)) * third
733 rby(10,nrb) = xiin
734 rby(11,nrb) = xiin
735 rby(12,nrb) = xiin
736 inmin = xiin
737 b1 = b1 + xiin
738 b5 = b5 + xiin
739 b9 = b9 + xiin
740 IF(nsubdom>0)THEN
741 IF(
tagno(m+n_part)==3)
THEN
742 WRITE(iout,1300) rbyid
743 WRITE(iout,1400)
744 ELSE
745 WRITE(iout,1100) rbyid,nonod,x(1,m),x(2,m),x(3,m),
746 . masrb,xiin,xiin,xiin,zero,zero,zero
747 END IF
748 ELSE
749 WRITE(iout,1100) rbyid,nonod,x(1,m),x(2,m),x(3,m),
750 . masrb,xiin,xiin,xiin,zero,zero,zero
751 ENDIF
752 ELSEIF (isph==2) THEN
753 inmin =
min(rby(10,nrb),rby(11,nrb),rby(12,nrb))
754 inmax =
max(rby(10,nrb),rby(11,nrb),rby(12,nrb))
755 IF(inmin<=1.e-3*inmax)THEN
756 IF(rby(10,nrb)/inmax<em03) rby(10,nrb)=rby(10,nrb)+em01*inmax
757 IF(rby(11,nrb)/inmax<em03) rby(11,nrb)=rby(11,nrb)+em01*inmax
758 IF(rby(12,nrb)/inmax<em03) rby(12,nrb)=rby(12,nrb)+em01*inmax
759 IF(nsubdom>0)THEN
760 IF(
tagno(m+n_part) /= 3)
THEN
762 . msgtype=msgwarning,
763 . anmode=aninfo_blind_1,
765 . c1=titr)
766 ENDIF
767 ELSE
769 . msgtype=msgwarning,
770 . anmode=aninfo_blind_1,
772 . c1=titr)
773 ENDIF
774 ENDIF
775 ENDIF
776
777
778 IF(nsubdom>0)THEN
779 IF(
tagno(m+n_part)==3)
GOTO 350
780 END IF
781
782
783 inmin =
min(rby(10,nrb),rby(11,nrb),rby(12,nrb))
784 WRITE(iout,1200) rby(10,nrb),rby(11,nrb),rby(12,nrb)
785 WRITE(iout,1101)
786 WRITE(iout,1102) (itab(lpby(i+npby(11,nrb))),i=1,nsl)
787
788 IF(rby(10,nrb)>=rby(11,nrb).AND.rby(10,nrb)>=rby(12,nrb))THEN
789 IF(rby(10,nrb)>(rby(11,nrb)+rby(12,nrb))*tol)THEN
791 . msgtype=msgwarning,
792 . anmode=aninfo_blind_1,
794 . c1=titr,
795 . r1=rby(10,nrb),
796 . r2=rby(11,nrb),
797 . r3=rby(12,nrb))
798 ENDIF
799 ELSEIF(rby(11,nrb)>=rby(10,nrb).AND.rby(11,nrb)>=rby(12,nrb))THEN
800 IF(rby(11,nrb)>(rby(10,nrb)+rby(12,nrb))*tol)THEN
802 . msgtype=msgwarning,
803 . anmode=aninfo_blind_1,
805 . c1=titr,
806 . r1=rby(11,nrb),
807 . r2=rby(10,nrb),
808 . r3=rby(12,nrb))
809 ENDIF
810 ELSEIF(rby(12,nrb)>=rby(10,nrb).AND.rby(12,nrb)>=rby(11,nrb))THEN
811 IF(rby(12,nrb)>(rby(10,nrb)+rby(11,nrb))*tol)THEN
813 . msgtype=msgwarning,
814 . anmode=aninfo_blind_1,
816 . c1=titr,
817 . r1=rby(12,nrb),
818 . r2=rby(10,nrb),
819 . r3=rby(11,nrb))
820 ENDIF
821 ENDIF
822
823 IF(inmin<=0.0)THEN
825 . msgtype=msgerror,
826 . anmode=aninfo_blind_1,
828 . c1=titr)
829 ENDIF
830 ENDIF
831
832350 CONTINUE
833
834 rby(13,nrb)=in(m)
835 rby(14,nrb)=masrb
836 rby(15,nrb)=ms(m)
837 ms(m) = masrb
838 in(m) =
min(rby(10,nrb),rby(11,nrb),rby(12,nrb))
839
840 DEALLOCATE(ms_loc)
841 DEALLOCATE(in_loc)
842
843
844
845
846
847 IF(rblevel == 0) THEN
848
849 IF(n2d == 0) THEN
850 xx=(x(1,m))*(x(1,m))
851 xy=(x(1,m))*(x(2,m))
852 xz=(x(1,m))*(x(3,m))
853 yy=(x(2,m))*(x(2,m))
854 yz=(x(2,m))*(x(3,m))
855 zz=(x(3,m))*(x(3,m))
856 b1 = b1 - in(m)
857 b5 = b5 - in(m)
858 b9 = b9 - in(m)
859 totmas = totmas - ms(m) + masrb
860 xgt = xgt - ms(m)*x(1,m) + masrb*x(1,m)
861 ygt = ygt - ms(m)*x(2,m) + masrb*x(2,m)
862 zgt = zgt - ms(m)*x(3,m) + masrb*x(3,m)
863
864
865 IF (ns10e>0) THEN
866 DO i=1,nsl
867 n = lsn(i)
868 IF (itagnd(n)/=0) cycle
869 stifn(m)= stifn(m)+stifn(n)
870 dd = (x(1,n)-x(1,m))**2+(x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
871 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
872 stifr(n)= em20
873 stifn(n)= em20
874 END DO
875 ELSE
876 DO i=1,nsl
877 n = lsn(i)
878 stifn(m)= stifn(m)+stifn(n)
879 dd = (x(1,n)-x(1,m))**2+(x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
880 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
881 stifr(n)= em20
882 stifn(n)= em20
883 END DO
884 END IF
885
886
887 ii1=rby(10,nrb)*rby(1,nrb)
888 ii2=rby(10,nrb)*rby(2,nrb)
889 ii3=rby(10,nrb)*rby(3,nrb)
890 ii4=rby(11,nrb)*rby(4,nrb)
891 ii5=rby(11,nrb)*rby(5,nrb)
892 ii6=rby(11,nrb)*rby(6,nrb)
893 ii7=rby(12,nrb)*rby(7,nrb)
894 ii8=rby(12,nrb)*rby(8,nrb)
895 ii9=rby(12,nrb)*rby(9,nrb)
896
897 rby(17,nrb)=rby(1,nrb)*ii1 + rby(4,nrb)*ii4 + rby(7,nrb)*ii7
898 rby(18,nrb)=rby(1,nrb)*ii2 + rby(4,nrb)*ii5 + rby(7,nrb)*ii8
899 rby(19,nrb)=rby(1,nrb)*ii3 + rby(4,nrb)*ii6 + rby(7,nrb)*ii9
900 rby(20,nrb)=rby(2,nrb)*ii1 + rby(5,nrb)*ii4 + rby
901 rby(21,nrb)=rby(2,nrb)*ii2 + rby(5,nrb)*ii5 + rby(8,nrb)*ii8
902 rby(22,nrb)=rby(2,nrb)*ii3 + rby(5,nrb)*ii6 + rby(8,nrb)*ii9
903 rby(23,nrb)=rby(3,nrb)*ii1 + rby(6,nrb)*ii4 + rby(9,nrb)*ii7
904 rby(24,nrb)=rby(3,nrb)*ii2 + rby(6,nrb)*ii5 + rby(9,nrb)*ii8
905 rby(25,nrb)=rby(3,nrb)*ii3 + rby(6,nrb)*ii6 + rby
906
907 ELSEIF(n2d == 1) THEN
908
909 b1 = b1 - in(m)
910 b5 = b5 - in(m)
911 b9 = b9 - in(m)
912 totmas = totmas - ms(m) + masrb
913 xgt = zero
914 ygt = ygt - ms(m)*x(2,m) + masrb*x(2,m)
915 zgt = zgt - ms(m)*x(3,m) + masrb*x(3,m)
916
917
918
919 IF (ns10e>0) THEN
920 DO i=1,nsl
921 n = lsn(i)
922 IF (itagnd(n)/=0) cycle
923 stifn(m)= stifn(m)+stifn(n)
924 dd =(x(1,n)-x(1,m))**2+(x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
925 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
926 stifr(n)= em20
927 stifn(n)= em20
928 END DO
929 ELSE
930 DO i=1,nsl
931 n = lsn(i)
932 stifn(m)= stifn(m)+stifn(n)
933 dd = (x(1,n)-x(1,m))**2+(x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
934 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
935 stifr(n)= em20
936 stifn(n)= em20
937 END DO
938 END IF
939
940
941 rby(17,nrb)=rby(10,nrb)
942 rby(18,nrb)=zero
943 rby(19,nrb)=zero
944 rby(20,nrb)=zero
945 rby(21,nrb)=rby(11,nrb)
946 rby(22,nrb)=zero
947 rby(23,nrb)=zero
948 rby(24,nrb)=zero
949 rby(25,nrb)=rby(12,nrb)
950
951 ELSEIF(n2d == 2) THEN
952 b1 = b1 - in(m)
953 b5 = b5 - in(m)
954 b9 = b9 - in(m)
955 totmas = totmas - ms(m) + masrb
956 xgt = zero
957 ygt = ygt - ms(m)*x(2,m) + masrb*x(2,m)
958 zgt = zgt - ms(m)*x(3,m) + masrb*x(3,m)
959
960
961 IF (ns10e>0) THEN
962 DO i=1,nsl
963 n = lsn(i)
964 IF (itagnd(n)/=0) cycle
965 stifn(m)= stifn(m)+stifn(n)
966 dd = (x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
967 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
968 stifr(n)= em20
969 stifn(n)= em20
970 END DO
971 ELSE
972 DO i=1,nsl
973 n = lsn(i)
974 stifn(m)= stifn(m)+stifn(n)
975 dd = (x(2,n)-x(2,m))**2+(x(3,n)-x(3,m))**2
976 stifr(m)= stifr(m)+(stifr(n)+dd*stifn(n))
977 stifr(n)= em20
978 stifn(n)= em20
979 END DO
980 END IF
981
982
983 ii1=rby(10,nrb)*rby(1,nrb)
984 ii2=zero
985 ii3=zero
986 ii4=rby(11,nrb)*rby(4,nrb)
987 ii5=rby(11,nrb)*rby(5,nrb)
988 ii6=rby(11,nrb)*rby(6,nrb)
989 ii7=zero
990 ii8=rby(12,nrb)*rby(8,nrb)
991 ii9=rby(12,nrb)*rby(9,nrb)
992
993 rby(17,nrb)=rby(1,nrb)*ii1 + rby(4,nrb)*ii4
994 rby(18,nrb)=rby(4,nrb)*ii5 + rby(7,nrb)*ii8
995 rby(19,nrb)=rby(4,nrb)*ii6 + rby(7,nrb)*ii9
996 rby(20,nrb)=rby(2,nrb)*ii1 + rby(5,nrb)*ii4
997 rby(21,nrb)=rby(5,nrb)*ii5 + rby(8,nrb)*ii8
998 rby(22,nrb)=rby(5,nrb)*ii6 + rby(8,nrb)*ii9
999 rby(23,nrb)=rby(3,nrb)*ii1 + rby(6,nrb)*ii4
1000 rby(24,nrb)=rby(6,nrb)*ii5 + rby(9,nrb)*ii8
1001 rby(25,nrb)=rby(6,nrb)*ii6 + rby(9,nrb)*ii9
1002 ENDIF
1003
1004
1005 IF(nsubdom>0)THEN
1006 IF(
tagno(m+n_part)==3)
THEN
1007 DO j=1,9
1008 rby(16+j,nrb)=rby_r2r(j)
1009 END DO
1010 ENDIF
1011 END IF
1012
1013
1014 IF (rby_iniaxis(1,nrb) > 0) THEN
1015 dist = sqrt((x(1,m)-x_msn0(1))**2+(x(2,m)-x_msn0(2))**2+(x(3,m)-x_msn0(3))**2)
1016 IF ((dist > zero).AND.
1017 . v(1,m)==rby_iniaxis(2,nrb).AND.
1018 . v(2,m)==rby_iniaxis(3,nrb).AND.
1019 . v(3,m)==rby_iniaxis(4,nrb)) THEN
1020
1021 v1x2=rby_iniaxis(5,nrb)*(x(2,m)-x_msn0(2))
1022 v2x1=rby_iniaxis(6,nrb)*(x(1,m)-x_msn0(1))
1023 v2x3=rby_iniaxis(6,nrb)*(x(3,m)-x_msn0(3))
1024 v3x2=rby_iniaxis(7,nrb)*(x(2,m)-x_msn0(2))
1025 v3x1=rby_iniaxis(7,nrb)*(x(1,m)-x_msn0(1))
1026 v1x3=rby_iniaxis(5,nrb)*(x(3,m)-x_msn0(3))
1027 v(1,m)= v(1,m)+v2x3-v3x2
1028 v(2,m)= v(2,m)+v3x1-v1x3
1029 v(3,m)= v(3,m)+v1x2-v2x1
1030 ENDIF
1031 ENDIF
1032
1033 IF(n2d == 0) THEN
1034 DO i=1,nsl
1035 n=lsn(i)
1036 v1x2=vr(1,m)*(x(2,n)-x(2,m))
1037 v2x1=vr(2,m)*(x(1,n)-x(1,m))
1038 v2x3=vr(2,m)*(x(3,n)-x(3,m))
1039 v3x2=vr(3,m)*(x(2,n)-x(2,m))
1040 v3x1=vr(3,m)*(x(1,n)-x(1,m))
1041 v1x3=vr(1,m)*(x(3,n)-x(3,m))
1042 v(1,n)= v(1,m)+v2x3-v3x2
1043 v(2,n)= v(2,m)+v3x1-v1x3
1044 v(3,n)= v(3,m)+v1x2-v2x1
1045 vr(1,n)= vr(1,m)
1046 vr(2,n)= vr(2,m)
1047 vr(3,n)= vr(3,m)
1048 ENDDO
1049 ELSEIF(n2d == 1) THEN
1050 v(1,m)= zero
1051 vr(1,m)= zero
1052 vr(2,m)= zero
1053 DO i=1,nsl
1054 n=lsn(i)
1055 v3x2=vr(3,m)*(x(2,n)-x(2,m))
1056 v3x1=vr(3,m)*(x(1,n)-x(1,m))
1057 v(1,n)= v(1,m)-v3x2
1058 v(2,n)= v(2,m)+v3x1
1059 v(3,n)= v(3,m)
1060 vr(1,n)= vr(1,m)
1061 vr(2,n)= vr(2,m)
1062 vr(3,n)= vr(3,m)
1063 ENDDO
1064 ELSEIF(n2d == 2) THEN
1065 v(1,m)= zero
1066 vr(2,m)= zero
1067 vr(3,m)= zero
1068 DO i=1,nsl
1069 n=lsn(i)
1070 v1x2=vr(1,m)*(x(2,n)-x(2,m))
1071 v1x3=vr(1,m)*(x(3,n)-x(3,m))
1072 v(1,n)= zero
1073 v(2,n)= v(2,m)-v1x3
1074 v(3,n)= v(3,m)+v1x2
1075 vr(1,n)= vr(1,m)
1076 vr(2,n)= zero
1077 vr(3,n)= zero
1078 ENDDO
1079 ENDIF
1080
1081 ENDIF
1082
1083 RETURN
1084
10851000 FORMAT(//
1086 . ' RIGID BODY INITIALIZATION '/
1087 . ' ------------------------- ')
1088
10891100 FORMAT(/5x,'RIGID BODY ID',i10
1090 . /10x,'PRIMARY NODE ',i10
1091 . /10x,'NEW X,Y,Z ',3g14.7
1092 . /10x,'NEW MASS ',1g14.7
1093 . /10x,'NEW INERTIA xx yy zz ',3g14.7
1094 . /10x,'NEW INERTIA xy yz zx ',3g14.7)
10951101 FORMAT(10x,'SECONDARY NODES ')
10961102 FORMAT( 10x,10i10)
1097
10981200 FORMAT(10x,'PRINCIPAL INERTIA',1p3g20.13)
10991300 FORMAT(/5x,'RIGID BODY ID',i10)
11001400 FORMAT(
1101 & 5x,40hrigid body on multidomains INTERFACE ,/,
1102 & 5x,55h --> mass and inertia matrix are computed in
the engine,/)
end diagonal values have been computed in the(sparse) matrix id.SOL
integer, parameter nchartitle
integer, dimension(:), allocatable tagno
subroutine inepri(xi, bm)
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)