OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
rskew33.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!|| rskew33 ../engine/source/elements/joint/rskew33.F
25!||--- called by ------------------------------------------------------
26!|| rgjoint ../engine/source/elements/joint/rgjoint.F
27!||--- calls -----------------------------------------------------
28!|| get_u_skew ../engine/source/user_interface/uaccess.F
29!|| inv3 ../engine/source/elements/joint/rskew33.F
30!|| prod_ab ../engine/source/elements/joint/rskew33.F
31!|| qrot33 ../engine/source/elements/joint/rskew33.F
32!||--- uses -----------------------------------------------------
33!|| element_mod ../common_source/modules/elements/element_mod.F90
34!|| sensor_mod ../common_source/modules/sensor_mod.F90
35!||====================================================================
36 SUBROUTINE rskew33(JFT ,JLT ,IXR ,IOUT ,IPROP,
37 . NUVAR ,UVAR ,RBY ,X ,XL ,
38 . ROT1 ,ROT2 ,DX ,DY ,DZ ,
39 . RX ,RY ,RZ ,VR ,IGTYP,
40 . NSENSOR,SENSOR_TAB,ISENS ,NC1 ,NC2 ,
41 . XDP)
42C-----------------------------------------------
43C M o d u l e s
44C-----------------------------------------------
45 USE sensor_mod
46 use element_mod , only : nixr
47C-----------------------------------------------
48C IXR | 5*NEL | I | R | SPRING CONNECTIVITY
49C | IXR(1,I) IPROP
50C | IXR(2,I) NODE 1 ID
51C | IXR(3,I) NODE 2 ID
52C | IXR(4,I) OPTIONAL NODE 3 ID
53C | IXR(5,I) Material ID
54C | IXR(6,I) SPRING ID
55C-----------------------------------------------
56C I m p l i c i t T y p e s
57C-----------------------------------------------
58#include "implicit_f.inc"
59C-----------------------------------------------
60C G l o b a l P a r a m e t e r s
61C-----------------------------------------------
62#include "mvsiz_p.inc"
63#include "param_c.inc"
64#include "com08_c.inc"
65#include "com04_c.inc"
66#include "com01_c.inc"
67#include "scr05_c.inc"
68C-----------------------------------------------
69C D u m m y A r g u m e n t s
70C-----------------------------------------------
71 INTEGER ,INTENT(IN) :: NSENSOR
72 INTEGER JFT, JLT, IOUT, NUVAR, IPROP, IXR(NIXR,*),IGTYP,
73 . ISENS,NC1(*),NC2(*)
74C REAL
75 my_real UVAR(NUVAR,*),X(3,*),
76 . rot1(3,mvsiz),rot2(3,mvsiz),rby(*),
77 . dx(*), dy(*), dz(*), rx(*), ry(*), rz(*),vr(3,*)
78 DOUBLE PRECISION XDP(3,*),XL(MVSIZ,3)
79 TYPE (SENSOR_STR_) ,DIMENSION(NSENSOR) :: SENSOR_TAB
80C-----------------------------------------------
81C L o c a l V a r i a b l e s
82C-----------------------------------------------
83 INTEGER I, J, K, JTYP, IERR, SKFLG,
84 . idsk1,idsk2,isk1,isk2,
85 . n1,n2,n3,usens,insens,
86 . isens_old,isens_act
87 my_real co,si,nx,
88 . u(lskew),v(lskew),ex(lskew),
89 . r1(3),r2(3),rm(3),
90 . a0(lskew),
91 . a(lskew),exi(lskew),
92 . get_u_geo,nr,dt(3),dex(lskew),exprec(lskew)
93 DOUBLE PRECISION X21(MVSIZ),Y21(MVSIZ),Z21(MVSIZ)
94C-----------------------------------------------
95C E x t e r n a l F u n c t i o n s
96C-----------------------------------------------
97 INTEGER GET_U_SKEW
98 EXTERNAL get_u_skew
99C=======================================================================
100 DO i=jft,jlt
101 nc1(i)= ixr(2,i)
102 nc2(i)= ixr(3,i)
103 ENDDO
104C
105 ierr = 0
106 isens = 0
107 isens_act = 0
108 jtyp = nint(get_u_geo(1,iprop))
109 IF (igtyp==33) THEN
110C---------------- skew initialisation for kjoints
111 idsk1 = nint(get_u_geo(2,iprop))
112 idsk2 = nint(get_u_geo(3,iprop))
113 skflg = nint(get_u_geo(14,iprop))
114 isk1 = get_u_skew(idsk1,n1,n2,n2,u)
115 isk2 = get_u_skew(idsk2,n1,n2,n3,v)
116 ELSE
117 skflg = 0
118C---------------- sensor check for kjoints2
119 usens = nint(get_u_geo(2,iprop))
120 isens_act = 0
121 IF (usens>0) THEN
122 isens_old = nint(uvar(16,jft))
123 DO k=1,nsensor
124 IF(usens==sensor_tab(k)%SENS_ID) insens=k
125 ENDDO
126C
127 IF (tt>sensor_tab(insens)%TSTART) THEN
128 isens = 1
129 uvar(16,jft) = isens
130 ENDIF
131 IF (isens/=isens_old) isens_act = 1
132 ENDIF
133 ENDIF
134C
135C ----------------
136C
137 DO i=jft,jlt
138C
139C---------------- initizialisation at time 0
140C
141 IF ((ncycle==0).AND.(tt==0)) THEN
142 IF (igtyp==33) THEN
143 DO j=1,9
144 a0(j)= uvar(3+j,i)
145 END DO
146 ELSE
147 rx(i)= uvar(7,i)
148 ry(i)= uvar(8,i)
149 rz(i)= uvar(9,i)
150 ENDIF
151 ENDIF
152C
153C---------------- stockage of displacements if sensor is activated
154C
155 IF ((igtyp==45).AND.(isens_act==1)) THEN
156 uvar(4,i) = dx(i)
157 uvar(5,i) = dy(i)
158 uvar(6,i) = dz(i)
159 uvar(7,i) = rx(i)
160 uvar(8,i) = ry(i)
161 uvar(9,i) = rz(i)
162 ENDIF
163C
164C---------------- local frame
165C
166 IF (jtyp==5) THEN
167C---------------------------------------------
168C---- universal joint - orthogonalized ----
169C---------------------------------------------
170 IF (igtyp==45) THEN
171C-------------------- nod 1
172 dt(1)= vr(1,nc1(i))*dt1
173 dt(2)= vr(2,nc1(i))*dt1
174 dt(3)= vr(3,nc1(i))*dt1
175 u(1)=uvar(10,i) - uvar(11,i)*dt(3)+uvar(12,i)*dt(2)
176 u(2)=uvar(11,i) - uvar(12,i)*dt(1)+uvar(10,i)*dt(3)
177 u(3)=uvar(12,i) - uvar(10,i)*dt(2)+uvar(11,i)*dt(1)
178 nr =sqrt(u(1)*u(1)+u(2)*u(2)+u(3)*u(3))
179 IF (nr>0) THEN
180 u(1)=u(1)/nr
181 u(2)=u(2)/nr
182 u(3)=u(3)/nr
183 ENDIF
184C-----
185 uvar(10,i) = u(1)
186 uvar(11,i) = u(2)
187 uvar(12,i) = u(3)
188
189C-------------------- nod 2
190 dt(1)= vr(1,nc2(i))*dt1
191 dt(2)= vr(2,nc2(i))*dt1
192 dt(3)= vr(3,nc2(i))*dt1
193 v(1)=uvar(13,i) - uvar(14,i)*dt(3)+uvar(15,i)*dt(2)
194 v(2)=uvar(14,i) - uvar(15,i)*dt(1)+uvar(13,i)*dt(3)
195 v(3)=uvar(15,i) - uvar(13,i)*dt(2)+uvar(14,i)*dt(1)
196 nr =sqrt(v(1)*v(1)+v(2)*v(2)+v(3)*v(3))
197 IF (nr>0) THEN
198 v(1)=v(1)/nr
199 v(2)=v(2)/nr
200 v(3)=v(3)/nr
201 ENDIF
202C-----
203 uvar(13,i) = v(1)
204 uvar(14,i) = v(2)
205 uvar(15,i) = v(3)
206 ENDIF
207
208 ex(1) = u(2)*v(3) - u(3)*v(2)
209 ex(2) = u(3)*v(1) - u(1)*v(3)
210 ex(3) = u(1)*v(2) - u(2)*v(1)
211 nx = sqrt(ex(1)*ex(1)+ex(2)*ex(2)+ex(3)*ex(3))
212 ex(1) = ex(1) / nx
213 ex(2) = ex(2) / nx
214 ex(3) = ex(3) / nx
215 ex(4) = u(1)
216 ex(5) = u(2)
217 ex(6) = u(3)
218 ex(7) = v(1)
219 ex(8) = v(2)
220 ex(9) = v(3)
221
222 CALL inv3(ex,exi)
223
224C---- Calculation of the rotation vector R12
225 x21(i) = (vr(1,nc2(i))-vr(1,nc1(i)))*dt1
226 y21(i) = (vr(2,nc2(i))-vr(2,nc1(i)))*dt1
227 z21(i) = (vr(3,nc2(i))-vr(3,nc1(i)))*dt1
228C
229 rm(1) = rx(i)+exi(1)*x21(i)+exi(4)*y21(i)+exi(7)*z21(i)
230 rm(2) = ry(i)+exi(2)*x21(i)+exi(5)*y21(i)+exi(8)*z21(i)
231 rm(3) = rz(i)+exi(3)*x21(i)+exi(6)*y21(i)+exi(9)*z21(i)
232C-----
233 r2(1) = 0.5*rm(1)
234 r2(2) = 0.5*rm(2)
235 r2(3) = 0.5*rm(3)
236 r1(1) = -0.5*rm(1)
237 r1(2) = -0.5*rm(2)
238 r1(3) = -0.5*rm(3)
239C-----
240 x21(i) = x(1,nc2(i))-x(1,nc1(i))
241 y21(i) = x(2,nc2(i))-x(2,nc1(i))
242 z21(i) = x(3,nc2(i))-x(3,nc1(i))
243 xl(i,1)=exi(1)*x21(i)+exi(4)*y21(i)+exi(7)*z21(i)
244 xl(i,2)=exi(2)*x21(i)+exi(5)*y21(i)+exi(8)*z21(i)
245 xl(i,3)=exi(3)*x21(i)+exi(6)*y21(i)+exi(9)*z21(i)
246
247C---------------------------
248 ELSE
249C---------------------------
250 IF (skflg==1) THEN
251C---------------------------------------------
252C---- first skew is used as mean skew -----
253C---------------------------------------------
254
255 DO j=1,9
256 ex(j) = u(j)
257 ENDDO
258
259 ELSE
260C-------------------------------------
261C---- mean skew is calculated -----
262C-------------------------------------
263
264C----- Initialization of EXPREC
265 IF ((ncycle==0).AND.(tt==0).AND.(igtyp==33)) THEN
266 CALL prod_ab(u,a0,a)
267 DO j=1,9
268 exprec(j) = a(j)
269 END DO
270 ELSE
271 DO j=1,9
272 exprec(j) = uvar(21+j,i)
273 END DO
274 ENDIF
275
276C----- Calculation of DEX
277 x21(i) = half*(vr(1,nc2(i))+vr(1,nc1(i)))*dt1
278 y21(i) = half*(vr(2,nc2(i))+vr(2,nc1(i)))*dt1
279 z21(i) = half*(vr(3,nc2(i))+vr(3,nc1(i)))*dt1
280 dt(1)=exprec(1)*x21(i)+exprec(2)*y21(i)+exprec(3)*z21(i)
281 dt(2)=exprec(4)*x21(i)+exprec(5)*y21(i)+exprec(6)*z21(i)
282 dt(3)=exprec(7)*x21(i)+exprec(8)*y21(i)+exprec(9)*z21(i)
283C-----
284 nr =sqrt(dt(1)*dt(1)+dt(2)*dt(2)+dt(3)*dt(3))
285 IF (nr>0) THEN
286 dt(1)=dt(1)/nr
287 dt(2)=dt(2)/nr
288 dt(3)=dt(3)/nr
289 ENDIF
290 co = cos(nr)
291 si = sin(nr)
292 CALL qrot33(dex, dt, co, si)
293
294C----- Calculation of new EX
295 CALL prod_ab(exprec,dex,ex)
296
297C-------------------------------
298 ENDIF
299C-------------------------------
300C---- Calculation of the rotation vector R12
301 x21(i) = (vr(1,nc2(i))-vr(1,nc1(i)))*dt1
302 y21(i) = (vr(2,nc2(i))-vr(2,nc1(i)))*dt1
303 z21(i) = (vr(3,nc2(i))-vr(3,nc1(i)))*dt1
304C
305 rm(1) = rx(i)+ex(1)*x21(i)+ex(2)*y21(i)+ex(3)*z21(i)
306 rm(2) = ry(i)+ex(4)*x21(i)+ex(5)*y21(i)+ex(6)*z21(i)
307 rm(3) = rz(i)+ex(7)*x21(i)+ex(8)*y21(i)+ex(9)*z21(i)
308C-----
309 r2(1) = 0.5*rm(1)
310 r2(2) = 0.5*rm(2)
311 r2(3) = 0.5*rm(3)
312 r1(1) = -0.5*rm(1)
313 r1(2) = -0.5*rm(2)
314 r1(3) = -0.5*rm(3)
315C-----
316 IF (iresp == 1) THEN
317C- simple precision - extended sple precsion only for translational dof
318 x21(i) = xdp(1,nc2(i))-xdp(1,nc1(i))
319 y21(i) = xdp(2,nc2(i))-xdp(2,nc1(i))
320 z21(i) = xdp(3,nc2(i))-xdp(3,nc1(i))
321 ELSE
322C- Double precision
323 x21(i) = x(1,nc2(i))-x(1,nc1(i))
324 y21(i) = x(2,nc2(i))-x(2,nc1(i))
325 z21(i) = x(3,nc2(i))-x(3,nc1(i))
326 ENDIF
327C
328 xl(i,1)=ex(1)*x21(i)+ex(2)*y21(i)+ex(3)*z21(i)
329 xl(i,2)=ex(4)*x21(i)+ex(5)*y21(i)+ex(6)*z21(i)
330 xl(i,3)=ex(7)*x21(i)+ex(8)*y21(i)+ex(9)*z21(i)
331
332C-------------------------------
333 ENDIF
334C-------------------------------
335C
336 rot1(1,i) = r1(1)
337 rot1(2,i) = r1(2)
338 rot1(3,i) = r1(3)
339 rot2(1,i) = r2(1)
340 rot2(2,i) = r2(2)
341 rot2(3,i) = r2(3)
342C
343 DO j=1,9
344 uvar(21+j,i) = ex(j)
345 END DO
346 END DO ! DO I=JFT,JLT
347
348C-----------
349 RETURN
350 END
351!||====================================================================
352!|| inv3 ../engine/source/elements/joint/rskew33.F
353!||--- called by ------------------------------------------------------
354!|| i2curvf ../engine/source/interfaces/interf/i2curvf.F
355!|| i2curvfp ../engine/source/interfaces/interf/i2curvfp.F
356!|| i2curvv ../engine/source/interfaces/interf/i2curvv.F
357!|| i2loceq ../common_source/interf/i2loceq.F
358!|| i2loceq_27 ../common_source/interf/i2loceq.F
359!|| rskew33 ../engine/source/elements/joint/rskew33.F
360!||====================================================================
361 SUBROUTINE inv3(A,B)
362C----------------------------------------------------------
363C I m p l i c i t T y p e s
364C-----------------------------------------------
365#include "implicit_f.inc"
366#include "param_c.inc"
367C----------------------------------------------------------
368C D u m m y A r g u m e n t s a n d F u n c t i o n
369C----------------------------------------------------------
370 my_real a(lskew),b(lskew)
371C-----------------------------------------------
372C L o c a l V a r i a b l e s
373C-----------------------------------------------
374 my_real det
375C=======================================================================
376
377 det = a(1)*a(5)*a(9)+a(4)*a(8)*a(3)+a(7)*a(2)*a(6)
378 . - a(4)*a(2)*a(9)-a(1)*a(8)*a(6)-a(7)*a(5)*a(3)
379
380 b(1) = (a(5)*a(9)-a(6)*a(8)) / det
381 b(4) = -(a(4)*a(9)-a(6)*a(7)) / det
382 b(7) = (a(4)*a(8)-a(5)*a(7)) / det
383
384 b(2) = -(a(2)*a(9)-a(3)*a(8)) / det
385 b(5) = (a(1)*a(9)-a(3)*a(7)) / det
386 b(8) = -(a(1)*a(8)-a(2)*a(7)) / det
387
388 b(3) = (a(2)*a(6)-a(3)*a(5)) / det
389 b(6) = -(a(1)*a(6)-a(3)*a(4)) / det
390 b(9) = (a(1)*a(5)-a(2)*a(4)) / det
391C-----------------------------------------------
392 RETURN
393 END
394
395!||====================================================================
396!|| rotq33 ../engine/source/elements/joint/rskew33.F
397!||====================================================================
398 SUBROUTINE rotq33(SKEW, T, ROT, C, S)
399C-------------------------------------------------------------------------
400C I m p l i c i t T y p e s
401C-----------------------------------------------
402#include "implicit_f.inc"
403#include "param_c.inc"
404C----------------------------------------------------------
405C D u m m y A r g u m e n t s a n d F u n c t i o n
406C----------------------------------------------------------
407 my_real rot(3), t(3), skew(lskew), s, c
408C-----------------------------------------------
409C L o c a l V a r i a b l e s
410C-----------------------------------------------
411 my_real e11,e22,e33,e12,e21,e13,e31,e23,e32,nr,ksi,
412 . pi1
413C=======================================================================
414CCsm45a2 PI= 2.*ATAN2(ONE,ZERO)
415 pi1 = 2.*atan2(one,zero)
416C--- first skew rotation
417 e11 = skew(1)
418 e12 = skew(4)
419 e13 = skew(7)
420 e21 = skew(2)
421 e22 = skew(5)
422 e23 = skew(8)
423 e31 = skew(3)
424 e32 = skew(6)
425 e33 = skew(9)
426 c = half * (e11+e22+e33 - one)
427 c = min(c,one)
428 c = max(c,-one)
429 ksi = acos(c)
430 s = sin(ksi)
431 IF(s/=zero) s = half / s
432 t(1) = (e32 - e23) * s
433 t(2) = (e13 - e31) * s
434 t(3) = (e21 - e12) * s
435 nr = sqrt(t(1)*t(1)+t(2)*t(2)+t(3)*t(3))
436 IF (nr/=zero) nr = one/nr
437 t(1) = t(1)*nr
438 t(2) = t(2)*nr
439 t(3) = t(3)*nr
440 rot(1) = t(1)*ksi
441 rot(2) = t(2)*ksi
442 rot(3) = t(3)*ksi
443C-----------------------------------------------
444 RETURN
445 END
446!||====================================================================
447!|| rot12 ../engine/source/elements/joint/rskew33.F
448!||====================================================================
449 SUBROUTINE rot12(SKEW, RVEC, C, S)
450C-------------------------------------------------------------------------
451C I m p l i c i t T y p e s
452C-----------------------------------------------
453#include "implicit_f.inc"
454#include "param_c.inc"
455C----------------------------------------------------------
456C D u m m y A r g u m e n t s a n d F u n c t i o n
457C----------------------------------------------------------
458 my_real skew(lskew), rvec(3), s, c
459C-----------------------------------------------
460C L o c a l V a r i a b l e s
461C-----------------------------------------------
462 my_real e11,e22,e33,e12,e21,e13,e31,e23,e32,nr,ksi,
463 . pi1
464C=======================================================================
465C
466 pi1 = 2.*atan2(one,zero)
467C--- first skew rotation
468 e11 = skew(1)
469 e12 = skew(4)
470 e13 = skew(7)
471 e21 = skew(2)
472 e22 = skew(5)
473 e23 = skew(8)
474 e31 = skew(3)
475 e32 = skew(6)
476 e33 = skew(9)
477 c = half * (e11+e22+e33 - one)
478 c = min(c,one)
479 c = max(c,-one)
480 ksi = acos(c)
481 s = sin(ksi)
482 IF(s/=zero) s = half / s
483 rvec(1) = (e32 - e23) * s
484 rvec(2) = (e13 - e31) * s
485 rvec(3) = (e21 - e12) * s
486 nr = sqrt(rvec(1)*rvec(1)+rvec(2)*rvec(2)+rvec(3)*rvec(3))
487 IF (nr/=zero) nr = one/nr
488 rvec(1) = rvec(1)*nr
489 rvec(2) = rvec(2)*nr
490 rvec(3) = rvec(3)*nr
491C
492 c = half*(c+ one)
493 s = sqrt(one-c)
494 c = sqrt(c)
495C---------------------------
496 RETURN
497 END
498!||====================================================================
499!|| qrot33 ../engine/source/elements/joint/rskew33.F
500!||--- called by ------------------------------------------------------
501!|| rskew33 ../engine/source/elements/joint/rskew33.F
502!||====================================================================
503 SUBROUTINE qrot33(SKEW, T, C, S)
504C-------------------------------------------------------------------------
505C I m p l i c i t T y p e s
506C-----------------------------------------------
507#include "implicit_f.inc"
508#include "param_c.inc"
509C----------------------------------------------------------
510C D u m m y A r g u m e n t s a n d F u n c t i o n
511C----------------------------------------------------------
512 my_real t(3), skew(lskew)
513C-----------------------------------------------
514C L o c a l V a r i a b l e s
515C-----------------------------------------------
516 my_real e11,e22,e33,e12,e21,e13,e31,e23,e32,
517 . u1,u2,u3, u1s, u2s, u3s, s,c,ci
518C=======================================================================
519 ci = one - c
520 u1 = t(1)
521 u2 = t(2)
522 u3 = t(3)
523 u1s = u1*s
524 u2s = u2*s
525 u3s = u3*s
526C
527 e11 = u1 * u1 *ci + c
528 e22 = u2 * u2 *ci + c
529 e33 = u3 * u3 *ci + c
530
531 e12 = u1 * u2 * ci
532 e21 = e12 + u3s
533 e12 = e12 - u3s
534 e13 = u1 * u3 * ci
535 e31 = e13 - u2s
536 e13 = e13 + u2s
537 e23 = u2 * u3 * ci
538 e32 = e23 + u1s
539 e23 = e23 - u1s
540C
541 skew(1) = e11
542 skew(4) = e12
543 skew(7) = e13
544 skew(2) = e21
545 skew(5) = e22
546 skew(8) = e23
547 skew(3) = e31
548 skew(6) = e32
549 skew(9) = e33
550C-----------------------------------------------
551 RETURN
552 END
553C====================================================
554!||====================================================================
555!|| prod_abt ../engine/source/elements/joint/rskew33.F
556!||====================================================================
557 SUBROUTINE prod_abt(A,B,X)
558#include "implicit_f.inc"
559#include "param_c.inc"
560 my_real a(lskew),b(lskew),x(lskew)
561C
562 x(1)=a(1)*b(1)+a(4)*b(4)+a(7)*b(7)
563 x(2)=a(2)*b(1)+a(5)*b(4)+a(8)*b(7)
564 x(3)=a(3)*b(1)+a(6)*b(4)+a(9)*b(7)
565 x(4)=a(1)*b(2)+a(4)*b(5)+a(7)*b(8)
566 x(5)=a(2)*b(2)+a(5)*b(5)+a(8)*b(8)
567 x(6)=a(3)*b(2)+a(6)*b(5)+a(9)*b(8)
568 x(7)=a(1)*b(3)+a(4)*b(6)+a(7)*b(9)
569 x(8)=a(2)*b(3)+a(5)*b(6)+a(8)*b(9)
570 x(9)=a(3)*b(3)+a(6)*b(6)+a(9)*b(9)
571 RETURN
572 END
573C====================================================
574!||====================================================================
575!|| prod_atb ../engine/source/elements/joint/rskew33.F
576!||====================================================================
577 SUBROUTINE prod_atb(A,B,X)
578#include "implicit_f.inc"
579#include "param_c.inc"
580 my_real a(lskew),b(lskew),x(lskew)
581C
582 x(1)=a(1)*b(1)+a(2)*b(2)+a(3)*b(3)
583 x(2)=a(4)*b(1)+a(5)*b(2)+a(6)*b(3)
584 x(3)=a(7)*b(1)+a(8)*b(2)+a(9)*b(3)
585 x(4)=a(1)*b(4)+a(2)*b(5)+a(3)*b(6)
586 x(5)=a(4)*b(4)+a(5)*b(5)+a(6)*b(6)
587 x(6)=a(7)*b(4)+a(8)*b(5)+a(9)*b(6)
588 x(7)=a(1)*b(7)+a(2)*b(8)+a(3)*b(9)
589 x(8)=a(4)*b(7)+a(5)*b(8)+a(6)*b(9)
590 x(9)=a(7)*b(7)+a(8)*b(8)+a(9)*b(9)
591 RETURN
592 END
593C====================================================
594!||====================================================================
595!|| prod_ab ../engine/source/elements/joint/rskew33.F
596!||--- called by ------------------------------------------------------
597!|| rskew33 ../engine/source/elements/joint/rskew33.F
598!||====================================================================
599 SUBROUTINE prod_ab(A,B,X)
600#include "implicit_f.inc"
601#include "param_c.inc"
602 my_real a(lskew),b(lskew),x(lskew)
603C
604 x(1)=a(1)*b(1)+a(4)*b(2)+a(7)*b(3)
605 x(2)=a(2)*b(1)+a(5)*b(2)+a(8)*b(3)
606 x(3)=a(3)*b(1)+a(6)*b(2)+a(9)*b(3)
607 x(4)=a(1)*b(4)+a(4)*b(5)+a(7)*b(6)
608 x(5)=a(2)*b(4)+a(5)*b(5)+a(8)*b(6)
609 x(6)=a(3)*b(4)+a(6)*b(5)+a(9)*b(6)
610 x(7)=a(1)*b(7)+a(4)*b(8)+a(7)*b(9)
611 x(8)=a(2)*b(7)+a(5)*b(8)+a(8)*b(9)
612 x(9)=a(3)*b(7)+a(6)*b(8)+a(9)*b(9)
613 RETURN
614 END
615C====================================================
#define my_real
Definition cppsort.cpp:32
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine qrot33(skew, t, c, s)
Definition rskew33.F:504
subroutine prod_atb(a, b, x)
Definition rskew33.F:578
subroutine rot12(skew, rvec, c, s)
Definition rskew33.F:450
subroutine prod_ab(a, b, x)
Definition rskew33.F:600
subroutine prod_abt(a, b, x)
Definition rskew33.F:558
subroutine rskew33(jft, jlt, ixr, iout, iprop, nuvar, uvar, rby, x, xl, rot1, rot2, dx, dy, dz, rx, ry, rz, vr, igtyp, nsensor, sensor_tab, isens, nc1, nc2, xdp)
Definition rskew33.F:42
subroutine rotq33(skew, t, rot, c, s)
Definition rskew33.F:399
subroutine inv3(a, b)
Definition rskew33.F:362