OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sms_rlink.F File Reference
#include "implicit_f.inc"
#include "comlock.inc"
#include "scr03_c.inc"
#include "com01_c.inc"
#include "task_c.inc"
#include "param_c.inc"
#include "com04_c.inc"
#include "com08_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine sms_rlink10 (ms, a, nlink, llink, skew, fr_rl, weight, frl6, idown, tag_lnk_sms, itab, frl)
subroutine sms_rlink11 (ms, a, nnlink, lllink, skew, fr_ll, weight, frl6, x, xframe, v, idown, tag_lnk_sms, itab, frl)
subroutine sms_rlink1 (ms, a, nsn, ic, nod, weight, frl6, iflag, idown, pmain, frl, tag_lnk, itab)
subroutine sms_rlink2 (ms, a, nsn, ic, nod, skew, weight, frl6, iflag, idown, pmain, frl, tag_lnk, itab)
subroutine sms_rlink3 (ms, x, a, v, nsn, ic, nod, xframe, weight, frl6, iflag, idown, pmain, frl, tag_lnk, itab)

Function/Subroutine Documentation

◆ sms_rlink1()

subroutine sms_rlink1 ( ms,
a,
integer nsn,
integer ic,
integer, dimension(*) nod,
integer, dimension(*) weight,
double precision, dimension(4,6) frl6,
integer iflag,
integer idown,
integer pmain,
frl,
integer tag_lnk,
integer, dimension(*) itab )

Definition at line 286 of file sms_rlink.F.

289C-----------------------------------------------
290C I m p l i c i t T y p e s
291C-----------------------------------------------
292#include "implicit_f.inc"
293C-----------------------------------------------
294C D u m m y A r g u m e n t s
295C-----------------------------------------------
296 INTEGER NSN, IC, IFLAG, IDOWN, PMAIN, TAG_LNK
297 INTEGER NOD(*),WEIGHT(*), ITAB(*)
298C REAL
299 my_real
300 . ms(*), a(3,*), frl(4)
301 DOUBLE PRECISION FRL6(4,6)
302C-----------------------------------------------
303C L o c a l V a r i a b l e s
304C-----------------------------------------------
305 INTEGER I, N, K
306C REAL
307 my_real
308 . mass, ax, ay, az,
309 . f1(nsn), f2(nsn), f3(nsn), f4(nsn)
310C-----------------------------------------------
311 SELECT CASE(idown)
312
313 CASE(0)
314C
315C Remontee
316 IF(iflag == 1)THEN
317
318C Init Parith/ON
319 DO k = 1, 6
320 frl6(1,k) = zero
321 frl6(2,k) = zero
322 frl6(3,k) = zero
323 frl6(4,k) = zero
324 END DO
325C
326 DO i=1,nsn
327 n = nod(i)
328 IF(weight(n)==1) THEN
329 f2(i)=a(1,n)
330 f3(i)=a(2,n)
331 f4(i)=a(3,n)
332 ELSE
333 f2(i)=zero
334 f3(i)=zero
335 f4(i)=zero
336 ENDIF
337 ENDDO
338C
339C Parith treatment/We before exchange
340C
341 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1), 4)
342 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1), 4)
343 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1), 4)
344
345 RETURN
346C
347 ELSEIF(iflag == 2)THEN
348C
349 ax = frl6(2,1)+frl6(2,2)+frl6(2,3)+
350 + frl6(2,4)+frl6(2,5)+frl6(2,6)
351 ay = frl6(3,1)+frl6(3,2)+frl6(3,3)+
352 + frl6(3,4)+frl6(3,5)+frl6(3,6)
353 az = frl6(4,1)+frl6(4,2)+frl6(4,3)+
354 + frl6(4,4)+frl6(4,5)+frl6(4,6)
355 DO i=1,nsn
356 n = nod(i)
357C
358C transmits force in the dion to the nd
359 IF(ic==1.OR.ic==3.OR.ic==5.OR.ic==7)THEN
360 a(3,n) =az
361 ENDIF
362 IF(ic==2.OR.ic==3.OR.ic==6.OR.ic==7)THEN
363 a(2,n) =ay
364 ENDIF
365 IF(ic==4.OR.ic==5.OR.ic==6.OR.ic==7)THEN
366 a(1,n) =ax
367 ENDIF
368 END DO
369C
370 ENDIF
371
372 CASE(1)
373C
374C Redescente
375 IF(iflag == 1)THEN
376C
377C Init Parith/ON
378 DO k = 1, 6
379 frl6(1,k) = zero
380 frl6(2,k) = zero
381 frl6(3,k) = zero
382 frl6(4,k) = zero
383 END DO
384C
385 DO i=1,nsn
386 n = nod(i)
387 IF(weight(n)==1) THEN
388 f1(i)=ms(n)
389 f2(i)=ms(n)*a(1,n)
390 f3(i)=ms(n)*a(2,n)
391 f4(i)=ms(n)*a(3,n)
392 ELSE
393 f1(i)=zero
394 f2(i)=zero
395 f3(i)=zero
396 f4(i)=zero
397 ENDIF
398 ENDDO
399C
400C Parith treatment/We before exchange
401C
402 CALL sum_6_float(1 ,nsn ,f1, frl6(1,1), 4)
403 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1), 4)
404 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1), 4)
405 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1), 4)
406
407 RETURN
408C
409 ELSEIF(iflag == 2)THEN
410C
411 mass = frl6(1,1)+frl6(1,2)+frl6(1,3)+
412 + frl6(1,4)+frl6(1,5)+frl6(1,6)
413 ax = frl6(2,1)+frl6(2,2)+frl6(2,3)+
414 + frl6(2,4)+frl6(2,5)+frl6(2,6)
415 ay = frl6(3,1)+frl6(3,2)+frl6(3,3)+
416 + frl6(3,4)+frl6(3,5)+frl6(3,6)
417 az = frl6(4,1)+frl6(4,2)+frl6(4,3)+
418 + frl6(4,4)+frl6(4,5)+frl6(4,6)
419 ax=ax/max(em30,mass)
420 ay=ay/max(em30,mass)
421 az=az/max(em30,mass)
422 DO i=1,nsn
423 n = nod(i)
424 IF(ic==1.OR.ic==3.OR.ic==5.OR.ic==7)THEN
425 a(3,n) =az
426 ENDIF
427 IF(ic==2.OR.ic==3.OR.ic==6.OR.ic==7)THEN
428 a(2,n) =ay
429 ENDIF
430 IF(ic==4.OR.ic==5.OR.ic==6.OR.ic==7)THEN
431 a(1,n) =ax
432 ENDIF
433 END DO
434C
435 ENDIF
436C
437 END SELECT
438 RETURN
#define my_real
Definition cppsort.cpp:32
#define max(a, b)
Definition macros.h:21
subroutine sum_6_float(jft, jlt, f, f6, n)
Definition parit.F:65

◆ sms_rlink10()

subroutine sms_rlink10 ( ms,
a,
integer, dimension(*) nlink,
integer, dimension(*) llink,
skew,
integer, dimension(nspmd+2,*) fr_rl,
integer, dimension(*) weight,
double precision, dimension(4,6,nrlink) frl6,
integer idown,
integer, dimension(*) tag_lnk_sms,
integer, dimension(*) itab,
frl )

Definition at line 32 of file sms_rlink.F.

36C-----------------------------------------------
37C I m p l i c i t T y p e s
38C-----------------------------------------------
39#include "implicit_f.inc"
40#include "comlock.inc"
41C-----------------------------------------------
42C C o m m o n B l o c k s
43C-----------------------------------------------
44#include "scr03_c.inc"
45#include "com01_c.inc"
46#include "task_c.inc"
47#include "param_c.inc"
48C-----------------------------------------------
49C D u m m y A r g u m e n t s
50C-----------------------------------------------
51 INTEGER NLINK(*),LLINK(*),FR_RL(NSPMD+2,*),WEIGHT(*), IDOWN,
52 . TAG_LNK_SMS(*), ITAB(*)
54 . ms(*), a(3,*), skew(lskew,*)
55 DOUBLE PRECISION FRL6(4,6,NRLINK)
57 . frl(4,nrlink)
58C-----------------------------------------------
59C L o c a l V a r i a b l e s
60C-----------------------------------------------
61 INTEGER K, K1, N, ISK, I, IC, KIND(NRLINK)
62C-----------------------------------------------
63C
64C Pre calcul index
65C
66 k = 1
67 DO i = 1, nrlink
68 kind(i) = k
69 k = k + nlink(4*i-3)
70 ENDDO
71C
72C parallel loop on SMP threads
73!$OMP DO SCHEDULE(DYNAMIC,1)
74 DO n=1,nrlink
75 frl(1,n)=zero
76 frl(2,n)=zero
77 frl(3,n)=zero
78 frl(4,n)=zero
79 IF(tag_lnk_sms(n) < 0)cycle
80
81 k1=4*n-3
82 ic=nlink(k1+1)
83 IF(ic==0) cycle
84
85 k = kind(n)
86 isk=nlink(k1+3)
87 IF(isk==1)THEN
88 CALL sms_rlink1(
89 1 ms ,a ,nlink(k1) ,nlink(k1+1),llink(k),
90 2 weight,frl6(1,1,n),1 ,idown ,fr_rl(nspmd+2,n),
91 3 frl(1,n),tag_lnk_sms(n),itab)
92 ELSE
93 CALL sms_rlink2(
94 1 ms ,a ,nlink(k1),nlink(k1+1),llink(k),
95 2 skew(1,isk),weight ,frl6(1,1,n),1 ,idown ,
96 3 fr_rl(nspmd+2,n),frl(1,n),tag_lnk_sms(n),itab)
97 ENDIF
98 END DO
99!$OMP END DO
100
101 IF(nspmd > 1) THEN
102!$OMP SINGLE
103
104C routine called rlink by rlink to optimize if needed
105 DO n=1,nrlink
106 IF(tag_lnk_sms(n)==0)cycle
107C calculation routine Parith/ON of A rigid link
108 IF(fr_rl(ispmd+1,n)/=0)
109 + CALL spmd_exch_fr6(fr_rl(1,n),frl6(1,1,n),4*6)
110 END DO
111
112!$OMP END SINGLE
113 END IF
114
115C parallel loop on SMP threads
116!$OMP DO SCHEDULE(DYNAMIC,1)
117 DO n=1,nrlink
118 IF(tag_lnk_sms(n) < 0)cycle
119
120 k1=4*n-3
121 ic=nlink(k1+1)
122 IF(ic==0) cycle
123
124 k = kind(n)
125 isk=nlink(k1+3)
126 IF(isk==1)THEN
127 CALL sms_rlink1(
128 1 ms ,a ,nlink(k1) ,nlink(k1+1),llink(k),
129 2 weight,frl6(1,1,n),2 ,idown ,fr_rl(nspmd+2,n),
130 3 frl(1,n),tag_lnk_sms(n),itab)
131 ELSE
132 CALL sms_rlink2(
133 1 ms ,a ,nlink(k1),nlink(k1+1),llink(k),
134 2 skew(1,isk),weight ,frl6(1,1,n),2 ,idown ,
135 3 fr_rl(nspmd+2,n),frl(1,n),tag_lnk_sms(n),itab)
136 ENDIF
137 END DO
138!$OMP END DO
139C
140 RETURN
subroutine spmd_exch_fr6(fr, fs6, len)

◆ sms_rlink11()

subroutine sms_rlink11 ( ms,
a,
integer, dimension(10,*) nnlink,
integer, dimension(*) lllink,
skew,
integer, dimension(nspmd+2,*) fr_ll,
integer, dimension(*) weight,
double precision, dimension(4,6,nlink) frl6,
x,
xframe,
v,
integer idown,
integer, dimension(*) tag_lnk_sms,
integer, dimension(*) itab,
frl )

Definition at line 152 of file sms_rlink.F.

156C-----------------------------------------------
157C I m p l i c i t T y p e s
158C-----------------------------------------------
159#include "implicit_f.inc"
160#include "comlock.inc"
161C-----------------------------------------------
162C C o m m o n B l o c k s
163C-----------------------------------------------
164#include "scr03_c.inc"
165#include "com01_c.inc"
166#include "com04_c.inc"
167#include "task_c.inc"
168#include "param_c.inc"
169C-----------------------------------------------
170C D u m m y A r g u m e n t s
171C-----------------------------------------------
172 INTEGER NNLINK(10,*), LLLINK(*), FR_LL(NSPMD+2,*),
173 . WEIGHT(*), IDOWN, TAG_LNK_SMS(*),ITAB(*)
174 my_real
175 . ms(*), a(3,*), skew(lskew,*),
176 . xframe(nxframe,*), x(3,*), v(3,*)
177 DOUBLE PRECISION FRL6(4,6,NLINK)
178 my_real
179 . frl(4,nlink)
180C-----------------------------------------------
181C L o c a l V a r i a b l e s
182C-----------------------------------------------
183 INTEGER K, N, ISK, I, IC, IPOL, KIND(NLINK)
184C-----------------------------------------------
185C
186C Pre calcul index
187C
188 k = 1
189 DO i = 1, nlink
190 kind(i) = k
191 k = k + nnlink(1,i)
192 ENDDO
193C parallel loop on SMP threads
194!$OMP DO SCHEDULE(DYNAMIC,1)
195 DO n=1,nlink
196 frl(1,n)=zero
197 frl(2,n)=zero
198 frl(3,n)=zero
199 frl(4,n)=zero
200 IF(tag_lnk_sms(nrlink+n) < 0)cycle
201
202 ic=nnlink(3,n)
203 IF(ic==0)cycle
204
205 isk=nnlink(5,n)
206 ipol=nnlink(6,n)
207 k = kind(n)
208 ipol=nnlink(6,n)
209 IF(ipol==1)THEN
210 CALL sms_rlink3(
211 1 ms ,x ,a ,v ,nnlink(1,n),
212 2 nnlink(3,n),lllink(k),xframe(1,isk),weight ,frl6(1,1,n),
213 3 1 ,idown ,fr_ll(nspmd+2,n),frl(1,n),tag_lnk_sms(nrlink+n),
214 4 itab)
215 ELSEIF(isk==1)THEN
216 CALL sms_rlink1(
217 1 ms ,a ,nnlink(1,n),nnlink(3,n),lllink(k),
218 2 weight ,frl6(1,1,n),1 ,idown ,fr_ll(nspmd+2,n),
219 3 frl(1,n) ,tag_lnk_sms(nrlink+n),itab)
220 ELSE
221 CALL sms_rlink2(
222 1 ms ,a ,nnlink(1,n),nnlink(3,n),lllink(k),
223 2 skew(1,isk), weight ,frl6(1,1,n),1 ,idown ,
224 3 fr_ll(nspmd+2,n),frl(1,n),tag_lnk_sms(nrlink+n),itab)
225 ENDIF
226 END DO
227!$OMP END DO
228 IF(nspmd > 1) THEN
229!$OMP SINGLE
230
231C routine called rlink by rlink to optimize if needed
232
233 DO n=1,nlink
234 IF(tag_lnk_sms(nrlink+n) < 0)cycle
235C calculation routine Parith/ON of A rigid link
236 IF(fr_ll(ispmd+1,n)/=0)
237 + CALL spmd_exch_fr6(fr_ll(1,n),frl6(1,1,n),4*6)
238 END DO
239
240!$OMP END SINGLE
241 END IF
242
243C parallel loop on SMP threads
244!$OMP DO SCHEDULE(DYNAMIC,1)
245 DO n=1,nlink
246 IF(tag_lnk_sms(nrlink+n) < 0)cycle
247
248 ic=nnlink(3,n)
249 IF(ic==0)cycle
250
251 isk=nnlink(5,n)
252 ipol=nnlink(6,n)
253 k = kind(n)
254 ipol=nnlink(6,n)
255 IF(ipol==1)THEN
256 CALL sms_rlink3(
257 1 ms ,x ,a ,v ,nnlink(1,n),
258 2 nnlink(3,n),lllink(k),xframe(1,isk),weight ,frl6(1,1,n),
259 3 2 ,idown ,fr_ll(nspmd+2,n),frl(1,n),tag_lnk_sms(nrlink+n),
260 4 itab)
261 ELSEIF(isk==1)THEN
262 CALL sms_rlink1(
263 1 ms ,a ,nnlink(1,n),nnlink(3,n),lllink(k),
264 2 weight ,frl6(1,1,n),2 ,idown ,fr_ll(nspmd+2,n),
265 3 frl(1,n) ,tag_lnk_sms(nrlink+n),itab)
266 ELSE
267 CALL sms_rlink2(
268 1 ms ,a ,nnlink(1,n),nnlink(3,n),lllink(k),
269 2 skew(1,isk),weight ,frl6(1,1,n),2 ,idown ,
270 3 fr_ll(nspmd+2,n),frl(1,n),tag_lnk_sms(nrlink+n),itab)
271 ENDIF
272 END DO
273!$OMP END DO
274C
275 RETURN

◆ sms_rlink2()

subroutine sms_rlink2 ( ms,
a,
integer nsn,
integer ic,
integer, dimension(*) nod,
skew,
integer, dimension(*) weight,
double precision, dimension(4,6) frl6,
integer iflag,
integer idown,
integer pmain,
frl,
integer tag_lnk,
integer, dimension(*) itab )

Definition at line 449 of file sms_rlink.F.

452C-----------------------------------------------
453C I m p l i c i t T y p e s
454C-----------------------------------------------
455#include "implicit_f.inc"
456C-----------------------------------------------
457C D u m m y A r g u m e n t s
458C-----------------------------------------------
459 INTEGER NSN, IC, IFLAG, IDOWN, PMAIN, TAG_LNK
460 INTEGER NOD(*),WEIGHT(*), ITAB(*)
461C REAL
462 my_real
463 . ms(*), a(3,*),skew(*), frl(4)
464 DOUBLE PRECISION FRL6(4,6)
465C-----------------------------------------------
466C L o c a l V a r i a b l e s
467C-----------------------------------------------
468 INTEGER IC1, ICC, IC2, IC3, I, N, K
469C REAL
470 my_real
471 . mass, ax, ay, az, dax, day, daz, aax, aay,
472 . aaz,
473 . f1(nsn), f2(nsn), f3(nsn), f4(nsn)
474C-----------------------------------------------
475 SELECT CASE(idown)
476
477 CASE(0)
478C
479C Remontee
480 IF(iflag == 1)THEN
481
482C Init Parith/ON
483 DO k = 1, 6
484 frl6(1,k) = zero
485 frl6(2,k) = zero
486 frl6(3,k) = zero
487 frl6(4,k) = zero
488 END DO
489
490C
491c AX =ZERO
492c AY =ZERO
493c AZ =ZERO
494C
495 DO i=1,nsn
496 n = nod(i)
497 IF(weight(n)==1) THEN
498 f2(i)=a(1,n)
499 f3(i)=a(2,n)
500 f4(i)=a(3,n)
501 ELSE
502 f2(i)=zero
503 f3(i)=zero
504 f4(i)=zero
505 ENDIF
506 ENDDO
507C
508C Parith treatment/We before exchange
509C
510 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1), 4)
511 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1), 4)
512 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1), 4)
513C
514 ELSEIF(iflag == 2)THEN
515 ic1=ic/4
516 icc=ic-4*ic1
517 ic2=icc/2
518 ic3=icc-2*ic2
519C
520 ax = frl6(2,1)+frl6(2,2)+frl6(2,3)+
521 + frl6(2,4)+frl6(2,5)+frl6(2,6)
522 ay = frl6(3,1)+frl6(3,2)+frl6(3,3)+
523 + frl6(3,4)+frl6(3,5)+frl6(3,6)
524 az = frl6(4,1)+frl6(4,2)+frl6(4,3)+
525 + frl6(4,4)+frl6(4,5)+frl6(4,6)
526 DO i=1,nsn
527 n = nod(i)
528C
529C transmits force in the dion to the nd
530 dax =a(1,n)-ax
531 day =a(2,n)-ay
532 daz =a(3,n)-az
533 aax =ic1*(skew(1)*dax+skew(2)*day+skew(3)*daz)
534 aay =ic2*(skew(4)*dax+skew(5)*day+skew(6)*daz)
535 aaz =ic3*(skew(7)*dax+skew(8)*day+skew(9)*daz)
536 a(1,n) =a(1,n)-aax*skew(1)-aay*skew(4)-aaz*skew(7)
537 a(2,n) =a(2,n)-aax*skew(2)-aay*skew(5)-aaz*skew(8)
538 a(3,n) =a(3,n)-aax*skew(3)-aay*skew(6)-aaz*skew(9)
539 END DO
540 END IF
541C
542
543 CASE(1)
544C
545C Redescente
546 IF(iflag == 1)THEN
547
548C Init Parith/ON
549 DO k = 1, 6
550 frl6(1,k) = zero
551 frl6(2,k) = zero
552 frl6(3,k) = zero
553 frl6(4,k) = zero
554 END DO
555
556C
557c MASS =ZERO
558c AX =ZERO
559c AY =ZERO
560c AZ =ZERO
561C
562 DO i=1,nsn
563 n = nod(i)
564 IF(weight(n)==1) THEN
565 f1(i)=ms(n)
566 f2(i)=ms(n)*a(1,n)
567 f3(i)=ms(n)*a(2,n)
568 f4(i)=ms(n)*a(3,n)
569 ELSE
570 f1(i)=zero
571 f2(i)=zero
572 f3(i)=zero
573 f4(i)=zero
574 ENDIF
575 ENDDO
576C
577C Parith treatment/We before exchange
578C
579 CALL sum_6_float(1 ,nsn ,f1, frl6(1,1), 4)
580 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1), 4)
581 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1), 4)
582 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1), 4)
583
584C
585 ELSEIF(iflag == 2)THEN
586C
587 ic1=ic/4
588 icc=ic-4*ic1
589 ic2=icc/2
590 ic3=icc-2*ic2
591C
592 mass = frl6(1,1)+frl6(1,2)+frl6(1,3)+
593 + frl6(1,4)+frl6(1,5)+frl6(1,6)
594 ax = frl6(2,1)+frl6(2,2)+frl6(2,3)+
595 + frl6(2,4)+frl6(2,5)+frl6(2,6)
596 ay = frl6(3,1)+frl6(3,2)+frl6(3,3)+
597 + frl6(3,4)+frl6(3,5)+frl6(3,6)
598 az = frl6(4,1)+frl6(4,2)+frl6(4,3)+
599 + frl6(4,4)+frl6(4,5)+frl6(4,6)
600 ax=ax/max(em30,mass)
601 ay=ay/max(em30,mass)
602 az=az/max(em30,mass)
603C
604 DO i=1,nsn
605 n = nod(i)
606 dax =a(1,n)-ax
607 day =a(2,n)-ay
608 daz =a(3,n)-az
609 aax =ic1*(skew(1)*dax+skew(2)*day+skew(3)*daz)
610 aay =ic2*(skew(4)*dax+skew(5)*day+skew(6)*daz)
611 aaz =ic3*(skew(7)*dax+skew(8)*day+skew(9)*daz)
612 a(1,n) =a(1,n)-aax*skew(1)-aay*skew(4)-aaz*skew(7)
613 a(2,n) =a(2,n)-aax*skew(2)-aay*skew(5)-aaz*skew(8)
614 a(3,n) =a(3,n)-aax*skew(3)-aay*skew(6)-aaz*skew(9)
615 END DO
616 END IF
617
618 END SELECT
619C
620 RETURN

◆ sms_rlink3()

subroutine sms_rlink3 ( ms,
x,
a,
v,
integer nsn,
integer ic,
integer, dimension(*) nod,
xframe,
integer, dimension(*) weight,
double precision, dimension(5,6) frl6,
integer iflag,
integer idown,
integer pmain,
frl,
integer tag_lnk,
integer, dimension(*) itab )

Definition at line 630 of file sms_rlink.F.

634C-----------------------------------------------
635C I m p l i c i t T y p e s
636C-----------------------------------------------
637#include "implicit_f.inc"
638C-----------------------------------------------
639C C o m m o n B l o c k s
640C-----------------------------------------------
641#include "com08_c.inc"
642C-----------------------------------------------
643C D u m m y A r g u m e n t s
644C-----------------------------------------------
645 INTEGER NSN, IC, IFLAG, IDOWN, PMAIN, TAG_LNK
646 INTEGER NOD(*),WEIGHT(*), ITAB(*)
647C REAL
648 my_real
649 . ms(*), x(3,*), a(3,*), v(3,*),
650 . xframe(*), frl(4)
651 DOUBLE PRECISION FRL6(5,6)
652C-----------------------------------------------
653C L o c a l V a r i a b l e s
654C-----------------------------------------------
655 INTEGER IC1, ICC, IC2, IC3, I, N, K
656C REAL
657 my_real
658 . mass, ax, ay, az,
659 . rx, ry, rz, sx, sy, sz,
660 . tx, ty, tz, ox, oy, oz, aa, r2, atr2,rx0, ry0, rz0, r,
661 . mr2, atmr2,
662 . f1(nsn), f2(nsn), f3(nsn), f4(nsn), f5(nsn)
663C-----------------------------------------------
664 SELECT CASE(idown)
665
666 CASE(0)
667C
668C Remontee
669 sx = xframe(1)
670 sy = xframe(2)
671 sz = xframe(3)
672 aa = one / sqrt(sx*sx + sy*sy + sz*sz)
673 sx = sx * aa
674 sy = sy * aa
675 sz = sz * aa
676 ox = xframe(10)
677 oy = xframe(11)
678 oz = xframe(12)
679
680 IF(iflag == 1)THEN
681
682C Init Parith/ON
683 DO k = 1, 6
684 frl6(1,k) = zero
685 frl6(2,k) = zero
686 frl6(3,k) = zero
687 frl6(4,k) = zero
688 END DO
689
690c R2 = ZERO
691c ATR2 = ZERO
692C-----------------------------------------------------------------------
693C Polar reference axis s' = x (frame), radius r, angle t
694C-----------------------------------------------------------------------
695 DO i=1,nsn
696 n = nod(i)
697 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
698 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
699 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
700 tx = ry0 * sz - rz0 * sy
701 ty = rz0 * sx - rx0 * sz
702 tz = rx0 * sy - ry0 * sx
703 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
704 tx = tx * aa
705 ty = ty * aa
706 tz = tz * aa
707 rx = sy * tz - sz * ty
708 ry = sz * tx - sx * tz
709 rz = sx * ty - sy * tx
710 r = rx * rx0 + ry * ry0 + rz * rz0
711 r = max(r,em20)
712C-----------------------------------------------------------------------
713C Change of overall mark => Polar
714C-----------------------------------------------------------------------
715 ax = a(1,n)
716 ay = a(2,n)
717 az = a(3,n)
718 a(1,n) = rx * ax + ry * ay + rz * az
719 a(2,n) = sx * ax + sy * ay + sz * az
720 a(3,n) = (tx * ax + ty * ay + tz * az) / r
721 IF(weight(n)==1) THEN
722 f1(i) = r*r
723 f2(i) = r*r*a(3,n)
724c R2 = R2 + R*R
725c ATR2 = ATR2 + R*R*A(3,N)
726 ELSE
727 f1(i) = zero
728 f2(i) = zero
729 ENDIF
730 ENDDO
731C
732C Parith treatment/We before exchange
733C
734 CALL sum_6_float(1 ,nsn ,f1, frl6(1,1),5)
735 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1),5)
736
737C-----------------------------------------------------------------------
738C mass quantity of motion
739C-----------------------------------------------------------------------
740 ax =zero
741 ay =zero
742 mass=zero
743 DO i=1,nsn
744 n = nod(i)
745 IF(weight(n)==1) THEN
746 f3(i)=a(1,n)
747 f4(i)=a(2,n)
748 ELSE
749 f3(i)=zero
750 f4(i)=zero
751 ENDIF
752 ENDDO
753C
754C Parith treatment/We before exchange
755C
756 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1),5)
757 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1),5)
758C
759 ELSEIF(iflag == 2)THEN
760C
761 ic1=ic/4
762 icc=ic-4*ic1
763 ic2=icc/2
764 ic3=icc-2*ic2
765C
766 r2 = frl6(1,1)+frl6(1,2)+frl6(1,3)+
767 + frl6(1,4)+frl6(1,5)+frl6(1,6)
768 atr2 = frl6(2,1)+frl6(2,2)+frl6(2,3)+
769 + frl6(2,4)+frl6(2,5)+frl6(2,6)
770 ax = frl6(3,1)+frl6(3,2)+frl6(3,3)+
771 + frl6(3,4)+frl6(3,5)+frl6(3,6)
772 ay = frl6(4,1)+frl6(4,2)+frl6(4,3)+
773 + frl6(4,4)+frl6(4,5)+frl6(4,6)
774
775 az= atr2/r2
776 DO i=1,nsn
777 n = nod(i)
778C
779C transmits force in the dion to nd 1
780 a(1,n) =a(1,n)-ic1*(a(1,n)-ax)
781 a(2,n) =a(2,n)-ic2*(a(2,n)-ay)
782 a(3,n) =a(3,n)-ic3*(a(3,n)-az)
783 END DO
784C-----------------------------------------------------------------------
785C Polar reference axis s, radius r, angle t
786C-----------------------------------------------------------------------
787 DO i=1,nsn
788 n = nod(i)
789 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
790 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
791 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
792 tx = ry0 * sz - rz0 * sy
793 ty = rz0 * sx - rx0 * sz
794 tz = rx0 * sy - ry0 * sx
795 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
796 tx = tx * aa
797 ty = ty * aa
798 tz = tz * aa
799 rx = sy * tz - sz * ty
800 ry = sz * tx - sx * tz
801 rz = sx * ty - sy * tx
802 r = rx * rx0 + ry * ry0 + rz * rz0
803 r = max(r,em20)
804C-----------------------------------------------------------------------
805C Change of polar frame => global
806C-----------------------------------------------------------------------
807 ax = rx * a(1,n) + sx * a(2,n) + tx * a(3,n) * r
808 ay = ry * a(1,n) + sy * a(2,n) + ty * a(3,n) * r
809 az = rz * a(1,n) + sz * a(2,n) + tz * a(3,n) * r
810 a(1,n) = ax
811 a(2,n) = ay
812 a(3,n) = az
813 ENDDO
814
815 END IF
816 CASE(1)
817C
818C Redescente
819
820 sx = xframe(1)
821 sy = xframe(2)
822 sz = xframe(3)
823 aa = one / sqrt(sx*sx + sy*sy + sz*sz)
824 sx = sx * aa
825 sy = sy * aa
826 sz = sz * aa
827 ox = xframe(10)
828 oy = xframe(11)
829 oz = xframe(12)
830C-----------------------------------------------------------------------
831C Polar reference axis s' = x (frame), radius r, angle t
832C-----------------------------------------------------------------------
833 IF(iflag == 1)THEN
834
835C Init Parith/ON
836 DO k = 1, 6
837 frl6(1,k) = zero
838 frl6(2,k) = zero
839 frl6(3,k) = zero
840 frl6(4,k) = zero
841 frl6(5,k) = zero
842 END DO
843
844c R2 = ZERO
845c ATR2 = ZERO
846C-----------------------------------------------------------------------
847C Polar reference axis s' = x (frame), radius r, angle t
848C-----------------------------------------------------------------------
849 DO i=1,nsn
850 n = nod(i)
851 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
852 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
853 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
854 tx = ry0 * sz - rz0 * sy
855 ty = rz0 * sx - rx0 * sz
856 tz = rx0 * sy - ry0 * sx
857 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
858 tx = tx * aa
859 ty = ty * aa
860 tz = tz * aa
861 rx = sy * tz - sz * ty
862 ry = sz * tx - sx * tz
863 rz = sx * ty - sy * tx
864 r = rx * rx0 + ry * ry0 + rz * rz0
865 r = max(r,em20)
866C-----------------------------------------------------------------------
867C Change of overall mark => Polar
868C-----------------------------------------------------------------------
869 ax = a(1,n)
870 ay = a(2,n)
871 az = a(3,n)
872 a(1,n) = rx * ax + ry * ay + rz * az
873 a(2,n) = sx * ax + sy * ay + sz * az
874 a(3,n) = (tx * ax + ty * ay + tz * az) / r
875 IF(weight(n)==1) THEN
876 f1(i) = r*r*ms(n)
877 f2(i) = r*r*ms(n)*a(3,n)
878c MR2 = MR2 + R*R*M
879c ATMR2 = ATMR2 + R*R*MS(N)*A(3,N)
880 ELSE
881 f1(i) = zero
882 f2(i) = zero
883 ENDIF
884 ENDDO
885C
886C Parith treatment/We before exchange
887C
888 CALL sum_6_float(1 ,nsn ,f1, frl6(1,1),5)
889 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1),5)
890
891C-----------------------------------------------------------------------
892C mass quantity of motion
893C-----------------------------------------------------------------------
894 ax =zero
895 ay =zero
896 mass=zero
897 DO i=1,nsn
898 n = nod(i)
899 IF(weight(n)==1) THEN
900 f3(i)=a(1,n)
901 f4(i)=a(2,n)
902 f5(i)=ms(n)
903 ELSE
904 f3(i)=zero
905 f4(i)=zero
906 f5(i)=zero
907 ENDIF
908 ENDDO
909C
910C Parith treatment/We before exchange
911C
912 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1),5)
913 CALL sum_6_float(1 ,nsn ,f4, frl6(4,1),5)
914 CALL sum_6_float(1 ,nsn ,f5, frl6(5,1),5)
915
916C
917C ELSEIF(IFLAG == 2)THEN
918C
919 ic1=ic/4
920 icc=ic-4*ic1
921 ic2=icc/2
922 ic3=icc-2*ic2
923C
924 mr2 = frl6(1,1)+frl6(1,2)+frl6(1,3)+
925 + frl6(1,4)+frl6(1,5)+frl6(1,6)
926 atmr2= frl6(2,1)+frl6(2,2)+frl6(2,3)+
927 + frl6(2,4)+frl6(2,5)+frl6(2,6)
928 ax = frl6(3,1)+frl6(3,2)+frl6(3,3)+
929 + frl6(3,4)+frl6(3,5)+frl6(3,6)
930 ay = frl6(4,1)+frl6(4,2)+frl6(4,3)+
931 + frl6(4,4)+frl6(4,5)+frl6(4,6)
932 mass = frl6(5,1)+frl6(5,2)+frl6(5,3)+
933 + frl6(5,4)+frl6(5,5)+frl6(5,6)
934
935 ax= ax/mass
936 ay= ay/mass
937 az= atmr2/mr2
938C-----------------------------------------------------------------------
939C link
940C-----------------------------------------------------------------------
941 DO i=1,nsn
942 n = nod(i)
943 a(1,n) =a(1,n)-ic1*(a(1,n)-ax)
944 a(2,n) =a(2,n)-ic2*(a(2,n)-ay)
945 a(3,n) =a(3,n)-ic3*(a(3,n)-az)
946 ENDDO
947C-----------------------------------------------------------------------
948C Polar reference axis s, radius r, angle t
949C-----------------------------------------------------------------------
950 DO i=1,nsn
951 n = nod(i)
952 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
953 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
954 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
955 tx = ry0 * sz - rz0 * sy
956 ty = rz0 * sx - rx0 * sz
957 tz = rx0 * sy - ry0 * sx
958 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
959 tx = tx * aa
960 ty = ty * aa
961 tz = tz * aa
962 rx = sy * tz - sz * ty
963 ry = sz * tx - sx * tz
964 rz = sx * ty - sy * tx
965 r = rx * rx0 + ry * ry0 + rz * rz0
966 r = max(r,em20)
967C-----------------------------------------------------------------------
968C Change of polar frame => global
969C-----------------------------------------------------------------------
970 ax = rx * a(1,n) + sx * a(2,n) + tx * a(3,n) * r
971 ay = ry * a(1,n) + sy * a(2,n) + ty * a(3,n) * r
972 az = rz * a(1,n) + sz * a(2,n) + tz * a(3,n) * r
973 a(1,n) = ax
974 a(2,n) = ay
975 a(3,n) = az
976 END DO
977 END IF
978
979 END SELECT
980C
981 RETURN