OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
rlink10.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| rlink10 ../engine/source/constraints/general/rlink/rlink10.F
25!||--- called by ------------------------------------------------------
26!|| resol ../engine/source/engine/resol.F
27!||--- calls -----------------------------------------------------
28!|| rlink1 ../engine/source/constraints/general/rlink/rlink1.F
29!|| rlink2 ../engine/source/constraints/general/rlink/rlink2.F
30!|| spmd_exch_fr6 ../engine/source/mpi/kinematic_conditions/spmd_exch_fr6.F
31!||====================================================================
32 SUBROUTINE rlink10(
33 1 MS ,IN ,A ,AR ,V ,
34 2 VR ,NLINK ,LLINK ,SKEW ,FR_RL,
35 3 WEIGHT,FRL6 )
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(*)
52 my_real
53 . MS(*), IN(*), A(3,*), AR(3,*), V(3,*), VR(3,*), SKEW(LSKEW,*)
54 DOUBLE PRECISION FRL6(15,6,NRLINK)
55C-----------------------------------------------
56C L o c a l V a r i a b l e s
57C-----------------------------------------------
58 INTEGER K, K1, N, ISK, I, KIND(NRLINK)
59C-----------------------------------------------
60C
61C Pre calcul index
62C
63 k = 1
64 DO i = 1, nrlink
65 kind(i) = k
66 k = k + nlink(4*i-3)
67 ENDDO
68C
69 k1=1
70
71C boucle parallele sur les threads SMP
72!$omp DO schedule(dynamic,1)
73 DO n=1,nrlink
74 k1=4*n-3
75 k = kind(n)
76 isk=nlink(k1+3)
77 IF(isk==1)THEN
78 CALL rlink1(
79 1 ms ,in ,a ,ar ,nlink(k1) ,
80 2 nlink(k1+1),nlink(k1+2),llink(k),weight,frl6(1,1,n),
81 3 1 )
82 ELSE
83 CALL rlink2(
84 1 ms ,in ,a ,ar ,v ,
85 2 vr ,nlink(k1),nlink(k1+1),nlink(k1+2),llink(k),
86 3 skew(1,isk),weight ,frl6(1,1,n),1 )
87 ENDIF
88 END DO
89!$OMP END DO
90
91 IF(nspmd > 1) THEN
92!$OMP SINGLE
93
94C routine appelee rlink par rlink a optimiser si besoin
95
96 DO n=1,nrlink
97C routine de calcul Parith/ON de A AR MS IN rigid link
98 IF(fr_rl(ispmd+1,n)/=0)
99 + CALL spmd_exch_fr6(fr_rl(1,n),frl6(1,1,n),15*6)
100 END DO
101
102!$OMP END SINGLE
103 END IF
104
105C boucle parallele sur les threads SMP
106!$OMP DO SCHEDULE(DYNAMIC,1)
107 DO n=1,nrlink
108 k1=4*n-3
109 k = kind(n)
110c DO I=1,N-1
111c K=K+NLINK(4*I-3)
112c ENDDO
113 isk=nlink(k1+3)
114 IF(isk==1)THEN
115 CALL rlink1(
116 1 ms ,in ,a ,ar ,nlink(k1) ,
117 2 nlink(k1+1),nlink(k1+2),llink(k),weight,frl6(1,1,n),
118 3 2 )
119 ELSE
120 CALL rlink2(
121 1 ms ,in ,a ,ar ,v ,
122 2 vr ,nlink(k1),nlink(k1+1),nlink(k1+2),llink(k),
123 3 skew(1,isk),weight ,frl6(1,1,n),2 )
124 ENDIF
125 END DO
126!$OMP END DO
127C
128 RETURN
129 END
130!||====================================================================
131!|| rlink11 ../engine/source/constraints/general/rlink/rlink10.F
132!||--- called by ------------------------------------------------------
133!|| resol ../engine/source/engine/resol.F
134!||--- calls -----------------------------------------------------
135!|| rlink1 ../engine/source/constraints/general/rlink/rlink1.F
136!|| rlink2 ../engine/source/constraints/general/rlink/rlink2.F
137!|| rlink3 ../engine/source/constraints/general/rlink/rlink10.F
138!|| spmd_exch_fr6 ../engine/source/mpi/kinematic_conditions/spmd_exch_fr6.F
139!||====================================================================
140 SUBROUTINE rlink11(
141 1 MS ,IN ,A ,AR ,V ,
142 2 VR ,NNLINK,LLLINK,SKEW ,FR_LL ,
143 3 WEIGHT,FRL6 ,X ,XFRAME )
144C-----------------------------------------------
145C I m p l i c i t T y p e s
146C-----------------------------------------------
147#include "implicit_f.inc"
148#include "comlock.inc"
149C-----------------------------------------------
150C C o m m o n B l o c k s
151C-----------------------------------------------
152#include "com01_c.inc"
153#include "com04_c.inc"
154#include "task_c.inc"
155#include "param_c.inc"
156C-----------------------------------------------
157C D u m m y A r g u m e n t s
158C-----------------------------------------------
159 INTEGER NNLINK(10,*), LLLINK(*), FR_LL(NSPMD+2,*),
160 . WEIGHT(*)
161 my_real
162 . MS(*), IN(*), A(3,*), AR(3,*), V(3,*), VR(3,*), SKEW(LSKEW,*),
163 . XFRAME(NXFRAME,*), X(3,*)
164 DOUBLE PRECISION FRL6(15,6,NLINK)
165C-----------------------------------------------
166C L o c a l V a r i a b l e s
167C-----------------------------------------------
168 INTEGER K, K1, N, ISK, I, IPOL, KIND(NLINK)
169C-----------------------------------------------
170C
171C Pre calcul index
172C
173 k = 1
174 DO i = 1, nlink
175 kind(i) = k
176 k = k + nnlink(1,i)
177 ENDDO
178C
179 k1=1
180
181C boucle parallele sur les threads SMP
182!$OMP DO SCHEDULE(DYNAMIC,1)
183 DO n=1,nlink
184c K = 1
185c DO I=1,N-1
186c K=K+NNLINK(1,I)
187c ENDDO
188 isk=nnlink(5,n)
189 ipol=nnlink(6,n)
190 k = kind(n)
191 ipol=nnlink(6,n)
192 IF(ipol==1)THEN
193 CALL rlink3(
194 1 ms ,in ,x ,a ,ar ,
195 2 v ,vr ,nnlink(1,n),nnlink(3,n),nnlink(4,n),
196 3 lllink(k),xframe(1,isk),weight ,frl6(1,1,n),1 )
197 ELSEIF(isk==1)THEN
198 CALL rlink1(
199 1 ms ,in ,a ,ar ,nnlink(1,n),
200 2 nnlink(3,n),nnlink(4,n),lllink(k),weight ,frl6(1,1,n),
201 3 1)
202 ELSE
203 CALL rlink2(
204 1 ms ,in ,a ,ar ,v ,
205 2 vr ,nnlink(1,n),nnlink(3,n),nnlink(4,n),lllink(k),
206 3 skew(1,isk),weight ,frl6(1,1,n),1 )
207 ENDIF
208 END DO
209!$OMP END DO
210
211 IF(nspmd > 1) THEN
212!$OMP SINGLE
213
214C routine appelee rlink par rlink a optimiser si besoin
215
216 DO n=1,nlink
217C routine de calcul Parith/ON de A AR MS IN rigid link
218 IF(fr_ll(ispmd+1,n)/=0)
219 + CALL spmd_exch_fr6(fr_ll(1,n),frl6(1,1,n),15*6)
220 END DO
221
222!$OMP END SINGLE
223 END IF
224
225C boucle parallele sur les threads SMP
226!$OMP DO SCHEDULE(DYNAMIC,1)
227 DO n=1,nlink
228 isk=nnlink(5,n)
229 ipol=nnlink(6,n)
230 k = kind(n)
231 ipol=nnlink(6,n)
232 IF(ipol==1)THEN
233 CALL rlink3(
234 1 ms ,in ,x ,a ,ar ,
235 2 v ,vr ,nnlink(1,n),nnlink(3,n),nnlink(4,n),
236 3 lllink(k),xframe(1,isk),weight ,frl6(1,1,n),2 )
237 ELSEIF(isk==1)THEN
238 CALL rlink1(
239 1 ms ,in ,a ,ar ,nnlink(1,n),
240 2 nnlink(3,n),nnlink(4,n),lllink(k),weight ,frl6(1,1,n),
241 3 2 )
242 ELSE
243 CALL rlink2(
244 1 ms ,in ,a ,ar ,v ,
245 2 vr ,nnlink(1,n),nnlink(3,n),nnlink(4,n),lllink(k),
246 3 skew(1,isk),weight ,frl6(1,1,n),2 )
247 ENDIF
248 END DO
249!$OMP END DO
250
251C
252 RETURN
253 END
254
255!||====================================================================
256!|| rlink3 ../engine/source/constraints/general/rlink/rlink10.F
257!||--- called by ------------------------------------------------------
258!|| rlink11 ../engine/source/constraints/general/rlink/rlink10.F
259!||--- calls -----------------------------------------------------
260!|| sum_6_float ../engine/source/system/parit.F
261!||====================================================================
262 SUBROUTINE rlink3(MS ,IN ,X ,A ,AR ,
263 2 V ,VR ,NSN ,IC ,ICR ,
264 3 NOD,XFRAME,WEIGHT,FRL6 ,IFLAG)
265C-----------------------------------------------
266C I m p l i c i t T y p e s
267C-----------------------------------------------
268#include "implicit_f.inc"
269C-----------------------------------------------
270C C o m m o n B l o c k s
271C-----------------------------------------------
272#include "com01_c.inc"
273#include "com08_c.inc"
274C-----------------------------------------------
275C D u m m y A r g u m e n t s
276C-----------------------------------------------
277 INTEGER NSN, IC, ICR, IFLAG
278 INTEGER NOD(*),WEIGHT(*)
279 my_real
280 . MS(*), IN(*), X(3,*), A(3,*), AR(3,*), V(3,*), VR(3,*),
281 . XFRAME(*)
282 DOUBLE PRECISION FRL6(15,6)
283C-----------------------------------------------
284C L o c a l V a r i a b l e s
285C-----------------------------------------------
286 INTEGER :: IC1, ICC, IC2, IC3, I, N, K
287 my_real :: mass, iner, ax, ay, az, vx, vy, vz,
288 . rx, ry, rz, sx, sy, sz,
289 . tx, ty, tz, ox, oy, oz, aa, mr2,vtmr2,atmr2,rx0, ry0, rz0,r,
290 . f1(nsn), f2(nsn), f3(nsn), f4(nsn), f5(nsn), f6(nsn), f7(nsn)
291C-----------------------------------------------
292C
293 IF(iflag == 1)THEN
294
295C Init Parith/ON
296 DO k = 1, 6
297 frl6(1,k) = zero
298 frl6(2,k) = zero
299 frl6(3,k) = zero
300 frl6(4,k) = zero
301 frl6(5,k) = zero
302 frl6(6,k) = zero
303 frl6(7,k) = zero
304 frl6(8,k) = zero
305 frl6(9,k) = zero
306 frl6(10,k) = zero
307 frl6(11,k) = zero
308 frl6(12,k) = zero
309 frl6(13,k) = zero
310 frl6(14,k) = zero
311 frl6(15,k) = zero
312 END DO
313
314 sx = xframe(1)
315 sy = xframe(2)
316 sz = xframe(3)
317 aa = one / sqrt(sx*sx + sy*sy + sz*sz)
318 sx = sx * aa
319 sy = sy * aa
320 sz = sz * aa
321 ox = xframe(10)
322 oy = xframe(11)
323 oz = xframe(12)
324c MR2 = ZERO
325c VTMR2 = ZERO
326c ATMR2 = ZERO
327C-----------------------------------------------------------------------
328C repere polaire axe S = X(frame) , rayon R , angle T
329C-----------------------------------------------------------------------
330 DO i=1,nsn
331 n = nod(i)
332 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
333 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
334 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
335 tx = ry0 * sz - rz0 * sy
336 ty = rz0 * sx - rx0 * sz
337 tz = rx0 * sy - ry0 * sx
338 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
339 tx = tx * aa
340 ty = ty * aa
341 tz = tz * aa
342 rx = sy * tz - sz * ty
343 ry = sz * tx - sx * tz
344 rz = sx * ty - sy * tx
345 r = rx * rx0 + ry * ry0 + rz * rz0
346 r = max(r,em20)
347C-----------------------------------------------------------------------
348C changement de repere global => polaire
349C-----------------------------------------------------------------------
350 ax = a(1,n)
351 ay = a(2,n)
352 az = a(3,n)
353 a(1,n) = rx * ax + ry * ay + rz * az
354 a(2,n) = sx * ax + sy * ay + sz * az
355 a(3,n) = (tx * ax + ty * ay + tz * az) / r
356 ax = v(1,n)
357 ay = v(2,n)
358 az = v(3,n)
359 v(1,n) = rx * ax + ry * ay + rz * az
360 v(2,n) = sx * ax + sy * ay + sz * az
361 v(3,n) = (tx * ax + ty * ay + tz * az) / r
362 IF(iroddl/=0)THEN
363 ax = ar(1,n)
364 ay = ar(2,n)
365 az = ar(3,n)
366 ar(1,n) = rx * ax + ry * ay + rz * az
367 ar(2,n) = sx * ax + sy * ay + sz * az
368 ar(3,n) = tx * ax + ty * ay + tz * az
369 ax = vr(1,n)
370 ay = vr(2,n)
371 az = vr(3,n)
372 vr(1,n) = rx * ax + ry * ay + rz * az
373 vr(2,n) = sx * ax + sy * ay + sz * az
374 vr(3,n) = tx * ax + ty * ay + tz * az
375 ENDIF
376 IF(weight(n)==1) THEN
377 f1(i) = r*r*ms(n)
378 f2(i) = r*r*ms(n)*v(3,n)
379 f3(i) = r*r*ms(n)*a(3,n)
380c MR2 = MR2 + R*R*MS(N)
381c VTMR2 = VTMR2 + R*R*MS(N)*V(3,N)
382c ATMR2 = ATMR2 + R*R*MS(N)*A(3,N)
383 ELSE
384 f1(i) = zero
385 f2(i) = zero
386 f3(i) = zero
387 ENDIF
388 ENDDO
389C
390 IF(ic/=0)THEN
391C
392C Traitement Parith/ON avant echange
393C
394 CALL sum_6_float(1 ,nsn ,f1, frl6(1,1), 15)
395 CALL sum_6_float(1 ,nsn ,f2, frl6(2,1), 15)
396 CALL sum_6_float(1 ,nsn ,f3, frl6(3,1), 15)
397
398C-----------------------------------------------------------------------
399C masse quantite de mouvement
400C-----------------------------------------------------------------------
401 ax =zero
402 ay =zero
403 vx =zero
404 vy =zero
405 mass=zero
406 DO i=1,nsn
407 n = nod(i)
408 IF(weight(n)==1) THEN
409 f1(i)=ms(n)
410 f2(i)=ms(n)*a(1,n)
411 f3(i)=ms(n)*a(2,n)
412 f4(i)=ms(n)*v(1,n)
413 f5(i)=ms(n)*v(2,n)
414 ELSE
415 f1(i)=zero
416 f2(i)=zero
417 f3(i)=zero
418 f4(i)=zero
419 f5(i)=zero
420 ENDIF
421 ENDDO
422C
423C Traitement Parith/ON avant echange
424C
425 CALL sum_6_float(1 ,nsn ,f1, frl6(4,1), 15)
426 CALL sum_6_float(1 ,nsn ,f2, frl6(5,1), 15)
427 CALL sum_6_float(1 ,nsn ,f3, frl6(6,1), 15)
428 CALL sum_6_float(1 ,nsn ,f4, frl6(7,1), 15)
429 CALL sum_6_float(1 ,nsn ,f5, frl6(8,1), 15)
430 ENDIF
431C-----------------------------------------------------------------------
432C DDL de rotation
433C-----------------------------------------------------------------------
434 IF(iroddl/=0.AND.icr/=0)THEN
435C-----------------------------------------------------------------------
436C inertie quantite de mouvement
437C-----------------------------------------------------------------------
438c AX =ZERO
439c AY =ZERO
440c AZ =ZERO
441c VX =ZERO
442c VY =ZERO
443c VZ =ZERO
444c INER=ZERO
445C
446 DO i=1,nsn
447 n = nod(i)
448 IF(weight(n)==1) THEN
449 f1(i)=in(n)
450 f2(i)=in(n)*ar(1,n)
451 f3(i)=in(n)*ar(2,n)
452 f4(i)=in(n)*ar(3,n)
453 f5(i)=in(n)*vr(1,n)
454 f6(i)=in(n)*vr(2,n)
455 f7(i)=in(n)*vr(3,n)
456 ENDIF
457 ENDDO
458C
459C
460C Traitement Parith/ON avant echange
461C
462 CALL sum_6_float(1 ,nsn ,f1, frl6(9,1), 15)
463 CALL sum_6_float(1 ,nsn ,f2, frl6(10,1), 15)
464 CALL sum_6_float(1 ,nsn ,f3, frl6(11,1), 15)
465 CALL sum_6_float(1 ,nsn ,f4, frl6(12,1), 15)
466 CALL sum_6_float(1 ,nsn ,f5, frl6(13,1), 15)
467 CALL sum_6_float(1 ,nsn ,f6, frl6(14,1), 15)
468 CALL sum_6_float(1 ,nsn ,f7, frl6(15,1), 15)
469 ENDIF
470C
471 ELSEIF(iflag == 2)THEN
472 sx = xframe(1)
473 sy = xframe(2)
474 sz = xframe(3)
475 aa = one / sqrt(sx*sx + sy*sy + sz*sz)
476 sx = sx * aa
477 sy = sy * aa
478 sz = sz * aa
479 ox = xframe(10)
480 oy = xframe(11)
481 oz = xframe(12)
482C
483 IF(ic/=0)THEN
484 ic1=ic/4
485 icc=ic-4*ic1
486 ic2=icc/2
487 ic3=icc-2*ic2
488C
489 mr2 = frl6(1,1)+frl6(1,2)+frl6(1,3)+
490 + frl6(1,4)+frl6(1,5)+frl6(1,6)
491 vtmr2= frl6(2,1)+frl6(2,2)+frl6(2,3)+
492 + frl6(2,4)+frl6(2,5)+frl6(2,6)
493 atmr2= frl6(3,1)+frl6(3,2)+frl6(3,3)+
494 + frl6(3,4)+frl6(3,5)+frl6(3,6)
495 mass = frl6(4,1)+frl6(4,2)+frl6(4,3)+
496 + frl6(4,4)+frl6(4,5)+frl6(4,6)
497 ax = frl6(5,1)+frl6(5,2)+frl6(5,3)+
498 + frl6(5,4)+frl6(5,5)+frl6(5,6)
499 ay = frl6(6,1)+frl6(6,2)+frl6(6,3)+
500 + frl6(6,4)+frl6(6,5)+frl6(6,6)
501 vx = frl6(7,1)+frl6(7,2)+frl6(7,3)+
502 + frl6(7,4)+frl6(7,5)+frl6(7,6)
503 vy = frl6(8,1)+frl6(8,2)+frl6(8,3)+
504 + frl6(8,4)+frl6(8,5)+frl6(8,6)
505
506 ax= ax/mass
507 ay= ay/mass
508 az= atmr2/mr2
509 vx= vx/mass
510 vy= vy/mass
511 vz= vtmr2/mr2
512C-----------------------------------------------------------------------
513C link
514C-----------------------------------------------------------------------
515 DO i=1,nsn
516 n = nod(i)
517 a(1,n) =a(1,n)-ic1*(a(1,n)-ax)
518 a(2,n) =a(2,n)-ic2*(a(2,n)-ay)
519 a(3,n) =a(3,n)-ic3*(a(3,n)-az)
520C
521 v(1,n) =v(1,n)-ic1*(v(1,n)-vx)
522 v(2,n) =v(2,n)-ic2*(v(2,n)-vy)
523 v(3,n) =v(3,n)-ic3*(v(3,n)-vz)
524C
525 ENDDO
526 END IF
527C-----------------------------------------------------------------------
528C DDL de rotation
529C-----------------------------------------------------------------------
530 IF(iroddl/=0.AND.icr/=0)THEN
531 ic1=icr/4
532 icc=icr-4*ic1
533 ic2=icc/2
534 ic3=icc-2*ic2
535C
536 iner= frl6(9,1)+frl6(9,2)+frl6(9,3)+
537 + frl6(9,4)+frl6(9,5)+frl6(9,6)
538 ax = frl6(10,1)+frl6(10,2)+frl6(10,3)+
539 + frl6(10,4)+frl6(10,5)+frl6(10,6)
540 ay = frl6(11,1)+frl6(11,2)+frl6(11,3)+
541 + frl6(11,4)+frl6(11,5)+frl6(11,6)
542 az = frl6(12,1)+frl6(12,2)+frl6(12,3)+
543 + frl6(12,4)+frl6(12,5)+frl6(12,6)
544 vx = frl6(13,1)+frl6(13,2)+frl6(13,3)+
545 + frl6(13,4)+frl6(13,5)+frl6(13,6)
546 vy = frl6(14,1)+frl6(14,2)+frl6(14,3)+
547 + frl6(14,4)+frl6(14,5)+frl6(14,6)
548 vz = frl6(15,1)+frl6(15,2)+frl6(15,3)+
549 + frl6(15,4)+frl6(15,5)+frl6(15,6)
550 IF(iner/=zero)THEN
551 ax=ax/iner
552 ay=ay/iner
553 az=az/iner
554 vx=vx/iner
555 vy=vy/iner
556 vz=vz/iner
557C-----------------------------------------------------------------------
558C link
559C-----------------------------------------------------------------------
560 DO i=1,nsn
561 n = nod(i)
562 ar(1,n) =ar(1,n)-ic1*(ar(1,n)-ax)
563 ar(2,n) =ar(2,n)-ic2*(ar(2,n)-ay)
564 ar(3,n) =ar(3,n)-ic3*(ar(3,n)-az)
565C
566 vr(1,n) =vr(1,n)-ic1*(vr(1,n)-vx)
567 vr(2,n) =vr(2,n)-ic2*(vr(2,n)-vy)
568 vr(3,n) =vr(3,n)-ic3*(vr(3,n)-vz)
569C
570 ENDDO
571 ENDIF
572 ENDIF
573C-----------------------------------------------------------------------
574C repere polaire axe S , rayon R , angle T
575C-----------------------------------------------------------------------
576 DO i=1,nsn
577 n = nod(i)
578 rx0 = x(1,n) + half * dt2 * v(1,n) - ox
579 ry0 = x(2,n) + half * dt2 * v(2,n) - oy
580 rz0 = x(3,n) + half * dt2 * v(3,n) - oz
581 tx = ry0 * sz - rz0 * sy
582 ty = rz0 * sx - rx0 * sz
583 tz = rx0 * sy - ry0 * sx
584 aa = one / sqrt(tx*tx + ty*ty + tz*tz)
585 tx = tx * aa
586 ty = ty * aa
587 tz = tz * aa
588 rx = sy * tz - sz * ty
589 ry = sz * tx - sx * tz
590 rz = sx * ty - sy * tx
591 r = rx * rx0 + ry * ry0 + rz * rz0
592 r = max(r,em20)
593C-----------------------------------------------------------------------
594C changement de repere polaire => global
595C-----------------------------------------------------------------------
596 ax = rx * a(1,n) + sx * a(2,n) + tx * a(3,n) * r
597 ay = ry * a(1,n) + sy * a(2,n) + ty * a(3,n) * r
598 az = rz * a(1,n) + sz * a(2,n) + tz * a(3,n) * r
599 a(1,n) = ax
600 a(2,n) = ay
601 a(3,n) = az
602 ax = rx * v(1,n) + sx * v(2,n) + tx * v(3,n) * r
603 ay = ry * v(1,n) + sy * v(2,n) + ty * v(3,n) * r
604 az = rz * v(1,n) + sz * v(2,n) + tz * v(3,n) * r
605 v(1,n) = ax
606 v(2,n) = ay
607 v(3,n) = az
608 IF(iroddl/=0)THEN
609 ax = rx * ar(1,n) + sx * ar(2,n) + tx * ar(3,n)
610 ay = ry * ar(1,n) + sy * ar(2,n) + ty * ar(3,n)
611 az = rz * ar(1,n) + sz * ar(2,n) + tz * ar(3,n)
612 ar(1,n) = ax
613 ar(2,n) = ay
614 ar(3,n) = az
615 ax = rx * vr(1,n) + sx * vr(2,n) + tx * vr(3,n)
616 ay = ry * vr(1,n) + sy * vr(2,n) + ty * vr(3,n)
617 az = rz * vr(1,n) + sz * vr(2,n) + tz * vr(3,n)
618 vr(1,n) = ax
619 vr(2,n) = ay
620 vr(3,n) = az
621 ENDIF
622 ENDDO
623 END IF
624C
625 RETURN
626 END
max
#define max(a, b)
Definition macros.h:21
sum_6_float
subroutine sum_6_float(jft, jlt, f, f6, n)
Definition parit.F:64
rlink11
subroutine rlink11(ms, in, a, ar, v, vr, nnlink, lllink, skew, fr_ll, weight, frl6, x, xframe)
Definition rlink10.F:144
rlink10
subroutine rlink10(ms, in, a, ar, v, vr, nlink, llink, skew, fr_rl, weight, frl6)
Definition rlink10.F:36
rlink3
subroutine rlink3(ms, in, x, a, ar, v, vr, nsn, ic, icr, nod, xframe, weight, frl6, iflag)
Definition rlink10.F:265
rlink1
subroutine rlink1(ms, in, a, ar, nsn, ic, icr, nod, weight, frl6, iflag)
Definition rlink1.F:34
rlink2
subroutine rlink2(ms, in, a, ar, v, vr, nsn, ic, icr, nod, skew, weight, frl6, iflag)
Definition rlink2.F:34
spmd_exch_fr6
subroutine spmd_exch_fr6(fr, fs6, len)
Definition spmd_exch_fr6.F:51