Télécharger la plaquette Cast3M


fr en

Télécharger qzvecp.eso

Retour à la liste

Numérotation des lignes : 

qzvecp

  1. C QZVECP    SOURCE    BP208322  16/01/21    21:15:21     8791
  2.       SUBROUTINE QZVECP (MATA, MATB, ALFR, ALFI, BETA, MATZ)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8 (A-H,O-Z)
  5. *
  6. *************************************************************************
  7. *                                                                       *
  8. *                CALCUL DES VECTEURS PROPRES COMPLEXES                  *
  9. *   ( QUATRIEME PAS DE LA METHODE "QZ" POUR LE CALCUL MODAL COMPLEXE )  *
  10. * _____________________________________________________________________ *
  11. *                                                                       *
  12. *                                                                       *
  13. *     DATE : le 7 Avril 1995                                            *
  14. *     AUTEURS : J.ANTUNES , Nicolas BENECH                              *
  15. * _____________________________________________________________________ *
  16. *                                                                       *
  17. *    MODULE(S) APPELANT(S) : VIBRAC                                     *
  18. *                                                                       *
  19. *    MODULE(S) APPELE(S) :                                              *
  20. * _____________________________________________________________________ *
  21. *                                                                       *
  22. *    EN ENTREE :                                                        *
  23. *       - MATA : matrice reelle quasi triangulaire superieure (XMATRI)  *
  24. *       - MATB : matrice reelle triangulaire superieure (XMATRI)        *
  25. *       - ALFR, ALFI, BETA : liste des valeurs propres                  *
  26. *            (parties reelle et imaginaire du numerateur + denominateur)*
  27. *       - MATZ : produit des transformations effectuees auparavant      *
  28. *                                                                       *
  29. *    EN SORTIE :                                                        *
  30. *       - MATA : matrice identique a la matrice d'entree                *
  31. *       - MATB : matrice detruite                                       *
  32. *       - ALFR, ALFI, BETA : liste des valeurs propres                  *
  33. *            (parties reelle et imaginaire du numerateur + denominateur)*
  34. *       - MATZ : matrice des vecteurs propres, normalises de facon      *
  35. *                telle que le module de la composante la plus large     *
  36. *                soit unitaire.                                         *
  37. *                si ALFI = 0 : le vecteur est dans la colonne J de MATZ *
  38. *                si ALFI <> 0 : la partie reelle est dans la colonne J, *
  39. *                               la partie imag. dans la colonne J+1.    *
  40. *                    (ce vecteur propre correspond a la valeur propre   *
  41. *                     de partie imagianire positive. Le vecteur propre  *
  42. *                     conjugue correspond a la valeur propre conjuguee) *
  43. *                                                                       *
  44. *************************************************************************
  45. *
  46.  
  47. -INC PPARAM
  48. -INC CCOPTIO
  49. -INC SMRIGID
  50. -INC SMLREEL
  51. C
  52.       REAL*8 SWAP
  53.       INTEGER I , J , K , M , EN , II , JJ , NA , NN , ISW ,
  54.      &        ENM2
  55. C-----------------------------------------------------------------------
  56.       POINTEUR MATA.XMATRI, MATB.XMATRI, MATZ.XMATRI
  57. cbp      POINTEUR ALFR.XMATRI, ALFI.XMATRI, BETA.XMATRI
  58.       POINTEUR ALFR.MLREEL, ALFI.MLREEL, BETA.MLREEL
  59. C
  60. c       SEGACT , MATA, ALFR,ALFI,BETA
  61. c       SEGACT , MATB*mod
  62. c       SEGACT , MATZ*mod
  63.       NM=MATA.RE(/1)
  64.       N=MATA.RE(/2)
  65. C
  66. C***********************************************************************
  67. C***********************************************************************
  68. C***********************************************************************
  69. C
  70. C=======================================================================
  71. C========================= DEBUT DES CALCULS ===========================
  72. C=======================================================================
  73. C
  74.       EPSB = MATB.RE(N,1,1)
  75.       ISW = 1
  76. C
  77. C-----------------------------------------------------------------------
  78. C-------------------------- BOUCLE PRINCIPALE --------------------------
  79. C-----------------------------------------------------------------------
  80. C
  81.       DO 800 NN = 1, N
  82.          EN = N + 1 - NN
  83.          NA = EN - 1
  84.          IF ( ISW . EQ . 2 ) GOTO 795
  85.          SWAP = ALFI.PROG(EN)
  86.          IF ( SWAP . NE . 0.0D0 ) GOTO 710
  87. C
  88. C-----------------------------------------------------------------------
  89. C------------------------- VECTEUR PROPRE REEL -------------------------
  90. C-----------------------------------------------------------------------
  91. C
  92.          M = EN
  93.          MATB.RE(EN,EN,1) = 1.0D0
  94.          IF ( NA . EQ . 0 ) GOTO 800
  95.          ALFM = ALFR.PROG(M)
  96.          BETM = BETA.PROG(M)
  97. C
  98.          DO 700 II = 1, NA
  99.             I = EN - II
  100.             W = BETM * MATA.RE(I,I,1) - ALFM * MATB.RE(I,I,1)
  101.             R = 0.0D0
  102. C
  103.             DO 610 J = M, EN
  104.   610       R = R + (BETM * MATA.RE(I,J,1)
  105.      &        - ALFM * MATB.RE(I,J,1)) * MATB.RE(J,EN,1)
  106. C
  107.             IF ( I . EQ . 1 . OR . ISW . EQ . 2 ) GOTO 630
  108.             IF ( BETM * MATA.RE(I,I-1,1) . EQ . 0.0D0 ) GOTO 630
  109.             ZZ = W
  110.             S = R
  111.             GOTO 690
  112.   630       M = I
  113.             IF ( ISW . EQ . 2 ) GOTO 640
  114. C
  115.             T = W
  116.             IF ( W . EQ . 0.0D0 ) T = EPSB
  117.             MATB.RE(I,EN,1) = -R / T
  118.             GOTO 700
  119. C
  120.   640       X = BETM * MATA.RE(I,I+1,1) - ALFM * MATB.RE(I,I+1,1)
  121.             Y = BETM * MATA.RE(I+1,I,1)
  122.             Q = W * ZZ - X * Y
  123.             T = (X * S - ZZ * R) / Q
  124.             MATB.RE(I,EN,1) = T
  125.             IF ( ABS(X) . LE . ABS(ZZ) ) GOTO 650
  126.             MATB.RE(I+1,EN,1) = (-R - W * T) / X
  127.             GOTO 690
  128.   650       MATB.RE(I+1,EN,1) = (-S - Y * T) / ZZ
  129.   690       ISW = 3 - ISW
  130.   700    CONTINUE
  131. C
  132.          GOTO 800
  133. C
  134. C-----------------------------------------------------------------------
  135. C----------------------- VECTEUR PROPRE COMPLEXE -----------------------
  136. C-----------------------------------------------------------------------
  137. C
  138.   710    M = NA
  139.          ALMR = ALFR.PROG(M)
  140.          ALMI = ALFI.PROG(M)
  141.          BETM = BETA.PROG(M)
  142. C
  143.          Y = BETM * MATA.RE(EN,NA,1)
  144.          MATB.RE(NA,EN,1) = (ALMR * MATB.RE(EN,EN,1)
  145.      &                  - BETM * MATA.RE(EN,EN,1)) / Y
  146.          MATB.RE(NA,NA,1) = -ALMI * MATB.RE(EN,EN,1) / Y
  147.          MATB.RE(EN,NA,1) = 0.0D0
  148.          MATB.RE(EN,EN,1) = 1.0D0
  149.          ENM2 = NA - 1
  150.          IF ( ENM2 . EQ . 0 ) GOTO 795
  151. C
  152.          DO 790 II = 1, ENM2
  153.             I = NA - II
  154.             RA = 0.0D0
  155.             SA = 0.0D0
  156. C
  157.             W = BETM * MATA.RE(I,I,1) - ALMR * MATB.RE(I,I,1)
  158.             W1 = -ALMI * MATB.RE(I,I,1)
  159. C
  160.             DO 760 J = M, EN
  161.                X = BETM * MATA.RE(I,J,1) - ALMR * MATB.RE(I,J,1)
  162.                X1 = -ALMI * MATB.RE(I,J,1)
  163.                RA = RA + X * MATB.RE(J,NA,1) - X1 * MATB.RE(J,EN,1)
  164.                SA = SA + X * MATB.RE(J,EN,1) + X1 * MATB.RE(J,NA,1)
  165.   760       CONTINUE
  166. C
  167.             IF ( I . EQ . 1 . OR . ISW . EQ . 2 ) GOTO 770
  168.             SWAP = MATA.RE(I,I-1)
  169.             IF ( BETM * SWAP . EQ . 0.0D0 ) GOTO 770
  170.             ZZ = W
  171.             Z1 = W1
  172.             R = RA
  173.             S = SA
  174.             ISW = 2
  175.             GOTO 790
  176.   770       M = I
  177.             IF ( ISW . EQ . 2 ) GOTO 780
  178. C
  179.             TR = -RA
  180.             TI = -SA
  181.   773       DR = W
  182.             DI = W1
  183. C
  184.   775       IF ( ABS(DI) . GT . ABS(DR) ) GOTO 777
  185.             RR = DI / DR
  186.             D = DR + DI * RR
  187.             T1 = (TR + TI * RR) / D
  188.             T2 = (TI - TR * RR) / D
  189.             GOTO (787,782), ISW
  190.   777       RR = DR / DI
  191.             D = DR * RR + DI
  192.             T1 = (TR * RR + TI) / D
  193.             T2 = (TI * RR - TR) / D
  194.             GOTO (787,782), ISW
  195. C
  196.   780       continue
  197.             X = BETM * MATA.RE(I,I+1,1) - ALMR * MATB.RE(I,I+1,1)
  198.             Y = BETM * MATA.RE(I+1,I,1)
  199.             X1 = -ALMI * MATB.RE(I,I+1,1)
  200.             TR = Y * RA - W * R + W1 * S
  201.             TI = Y * SA - W * S - W1 * R
  202.             DR = W * ZZ - W1 * Z1 - X * Y
  203.             DI = W * Z1 + W1 * ZZ - X1 * Y
  204.             IF ( DR . EQ . 0.0D0 . AND . DI . EQ . 0.0D0 ) DR = EPSB
  205.             GOTO 775
  206.   782       continue
  207.             MATB.RE(I+1,NA,1) = T1
  208.             MATB.RE(I+1,EN,1) = T2
  209.             ISW = 1
  210.             IF ( ABS(Y) . GT . ABS(W) + ABS(W1) ) GOTO 785
  211.             TR = -RA - X * MATB.RE(I+1,NA,1) + X1 * MATB.RE(I+1,EN,1)
  212.             TI = -SA - X * MATB.RE(I+1,EN,1) - X1 * MATB.RE(I+1,NA,1)
  213.             GOTO 773
  214.   785       continue
  215.             T1 = (-R - ZZ * MATB.RE(I+1,NA,1)
  216.      &        + Z1 * MATB.RE(I+1,EN,1)) / Y
  217.             T2 = (-S - ZZ * MATB.RE(I+1,EN,1)
  218.      &        - Z1 * MATB.RE(I+1,NA,1)) / Y
  219.   787       continue
  220.             MATB.RE(I,NA,1) = T1
  221.             MATB.RE(I,EN,1) = T2
  222.   790    CONTINUE
  223. C
  224.   795    ISW = 3 - ISW
  225. C
  226. C-----------------------------------------------------------------------
  227. C------------------- FIN DE LA SUBSTITUTION INVERSE --------------------
  228. C-----------------------------------------------------------------------
  229. C
  230.   800 CONTINUE
  231. C
  232. C-----------------------------------------------------------------------
  233. C------- ON TRANSFORME DANS LE SYSTEME DE COORDONNEES D ORIGINE --------
  234. C-----------------------------------------------------------------------
  235. C
  236.       DO 880 JJ = 1, N
  237.          J = N + 1 - JJ
  238. C
  239.          DO 880 I = 1, N
  240.             ZZ = 0.0D0
  241. C
  242.             DO 860 K = 1, J
  243.   860       ZZ = ZZ + MATZ.RE(I,K,1) * MATB.RE(K,J,1)
  244. C
  245.             MATZ.RE(I,J,1) = ZZ
  246.   880 CONTINUE
  247. C
  248. C-----------------------------------------------------------------------
  249. C----------------- NORMALISATION DES VECTEURS PROPRES ------------------
  250. C-----------------------------------------------------------------------
  251. C
  252.       DO 950 J = 1, N
  253.          D = 0.0D0
  254.          IF ( ISW . EQ . 2 ) GOTO 920
  255.          SWAP = ALFI.PROG(J)
  256.          IF ( SWAP . NE . 0.0D0 ) GOTO 945
  257. C
  258.          DO 890 I = 1, N
  259.             IF ( ABS(MATZ.RE(I,J,1)) . GT . D )
  260.      &       D = ABS(MATZ.RE(I,J,1))
  261.   890    CONTINUE
  262. C
  263.          DO 900 I = 1, N
  264.   900      MATZ.RE(I,J,1) = MATZ.RE(I,J,1) / D
  265. C
  266.          GOTO 950
  267. C
  268.   920       continue
  269.          DO 930 I = 1, N
  270.             R = ABS(MATZ.RE(I,J-1,1)) + ABS(MATZ.RE(I,J,1))
  271.             IF ( R . NE . 0.0D0 )
  272.      &        R =     SQRT((MATZ.RE(I,J-1,1)  )**2
  273.      &          + (MATZ.RE(I,J,1)  )**2)
  274.             IF ( R . GT . D ) D = R
  275.   930    CONTINUE
  276. C
  277.          DO 940 I = 1, N
  278.             MATZ.RE(I,J-1,1) = MATZ.RE(I,J-1,1) / D
  279.             MATZ.RE(I,J,1) = MATZ.RE(I,J,1) / D
  280.   940    CONTINUE
  281. C
  282.   945    ISW = 3 - ISW
  283.   950 CONTINUE
  284. C
  285. c       SEGDES , MATA,MATB,MATZ
  286. c       SEGDES , ALFR,ALFI,BETA
  287.       SEGSUP , MATA,MATB
  288. c
  289. C=======================================================================
  290. C========================== FIN DES CALCULS ============================
  291. C=======================================================================
  292. C
  293.       RETURN
  294.       END
  295.  
  296.  
  297.  
  298.  
  299.  
  300.  
  301.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales