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    PV090527  26/04/30    21:16:05     12529          
  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     J = M, EN
  104.             R = R + (BETM * MATA.RE(I,J,1) 
  105.      &        - ALFM * MATB.RE(I,J,1)) * MATB.RE(J,EN,1)
  106.             enddo
  107. C
  108.             IF ( I . EQ . 1 . OR . ISW . EQ . 2 ) GOTO 630
  109.             IF ( BETM * MATA.RE(I,I-1,1) . EQ . 0.0D0 ) GOTO 630
  110.             ZZ = W
  111.             S = R
  112.             GOTO 690
  113.   630       M = I
  114.             IF ( ISW . EQ . 2 ) GOTO 640
  115. C
  116.             T = W
  117.             IF ( W . EQ . 0.0D0 ) T = EPSB
  118.             MATB.RE(I,EN,1) = -R / T
  119.             GOTO 700
  120. C
  121.   640       X = BETM * MATA.RE(I,I+1,1) - ALFM * MATB.RE(I,I+1,1)
  122.             Y = BETM * MATA.RE(I+1,I,1)
  123.             Q = W * ZZ - X * Y
  124.             T = (X * S - ZZ * R) / Q
  125.             MATB.RE(I,EN,1) = T
  126.             IF ( ABS(X) . LE . ABS(ZZ) ) GOTO 650
  127.             MATB.RE(I+1,EN,1) = (-R - W * T) / X
  128.             GOTO 690
  129.   650       MATB.RE(I+1,EN,1) = (-S - Y * T) / ZZ
  130.   690       ISW = 3 - ISW
  131.   700    CONTINUE
  132. C
  133.          GOTO 800
  134. C
  135. C-----------------------------------------------------------------------
  136. C----------------------- VECTEUR PROPRE COMPLEXE -----------------------
  137. C-----------------------------------------------------------------------
  138. C
  139.   710    M = NA
  140.          ALMR = ALFR.PROG(M)
  141.          ALMI = ALFI.PROG(M)
  142.          BETM = BETA.PROG(M)
  143. C
  144.          Y = BETM * MATA.RE(EN,NA,1)
  145.          MATB.RE(NA,EN,1) = (ALMR * MATB.RE(EN,EN,1)
  146.      &                  - BETM * MATA.RE(EN,EN,1)) / Y
  147.          MATB.RE(NA,NA,1) = -ALMI * MATB.RE(EN,EN,1) / Y
  148.          MATB.RE(EN,NA,1) = 0.0D0
  149.          MATB.RE(EN,EN,1) = 1.0D0
  150.          ENM2 = NA - 1
  151.          IF ( ENM2 . EQ . 0 ) GOTO 795
  152. C
  153.          DO 790 II = 1, ENM2
  154.             I = NA - II
  155.             RA = 0.0D0
  156.             SA = 0.0D0
  157. C
  158.             W = BETM * MATA.RE(I,I,1) - ALMR * MATB.RE(I,I,1)
  159.             W1 = -ALMI * MATB.RE(I,I,1)
  160. C
  161.             DO 760 J = M, EN
  162.                X = BETM * MATA.RE(I,J,1) - ALMR * MATB.RE(I,J,1)
  163.                X1 = -ALMI * MATB.RE(I,J,1)
  164.                RA = RA + X * MATB.RE(J,NA,1) - X1 * MATB.RE(J,EN,1)
  165.                SA = SA + X * MATB.RE(J,EN,1) + X1 * MATB.RE(J,NA,1)
  166.   760       CONTINUE
  167. C
  168.             IF ( I . EQ . 1 . OR . ISW . EQ . 2 ) GOTO 770
  169.             SWAP = MATA.RE(I,I-1)
  170.             IF ( BETM * SWAP . EQ . 0.0D0 ) GOTO 770
  171.             ZZ = W
  172.             Z1 = W1
  173.             R = RA
  174.             S = SA
  175.             ISW = 2
  176.             GOTO 790
  177.   770       M = I
  178.             IF ( ISW . EQ . 2 ) GOTO 780
  179. C
  180.             TR = -RA
  181.             TI = -SA
  182.   773       DR = W
  183.             DI = W1
  184. C
  185.   775       IF ( ABS(DI) . GT . ABS(DR) ) GOTO 777
  186.             RR = DI / DR
  187.             D = DR + DI * RR
  188.             T1 = (TR + TI * RR) / D
  189.             T2 = (TI - TR * RR) / D
  190.             GOTO (787,782), ISW
  191.   777       RR = DR / DI
  192.             D = DR * RR + DI
  193.             T1 = (TR * RR + TI) / D
  194.             T2 = (TI * RR - TR) / D
  195.             GOTO (787,782), ISW
  196. C
  197.   780       continue
  198.             X = BETM * MATA.RE(I,I+1,1) - ALMR * MATB.RE(I,I+1,1)
  199.             Y = BETM * MATA.RE(I+1,I,1)
  200.             X1 = -ALMI * MATB.RE(I,I+1,1)
  201.             TR = Y * RA - W * R + W1 * S
  202.             TI = Y * SA - W * S - W1 * R
  203.             DR = W * ZZ - W1 * Z1 - X * Y
  204.             DI = W * Z1 + W1 * ZZ - X1 * Y
  205.             IF ( DR . EQ . 0.0D0 . AND . DI . EQ . 0.0D0 ) DR = EPSB
  206.             GOTO 775
  207.   782       continue
  208.             MATB.RE(I+1,NA,1) = T1
  209.             MATB.RE(I+1,EN,1) = T2
  210.             ISW = 1
  211.             IF ( ABS(Y) . GT . ABS(W) + ABS(W1) ) GOTO 785
  212.             TR = -RA - X * MATB.RE(I+1,NA,1) + X1 * MATB.RE(I+1,EN,1)
  213.             TI = -SA - X * MATB.RE(I+1,EN,1) - X1 * MATB.RE(I+1,NA,1)
  214.             GOTO 773
  215.   785       continue
  216.             T1 = (-R - ZZ * MATB.RE(I+1,NA,1)
  217.      &        + Z1 * MATB.RE(I+1,EN,1)) / Y
  218.             T2 = (-S - ZZ * MATB.RE(I+1,EN,1)
  219.      &        - Z1 * MATB.RE(I+1,NA,1)) / Y
  220.   787       continue
  221.             MATB.RE(I,NA,1) = T1
  222.             MATB.RE(I,EN,1) = T2
  223.   790    CONTINUE
  224. C
  225.   795    ISW = 3 - ISW
  226. C
  227. C-----------------------------------------------------------------------
  228. C------------------- FIN DE LA SUBSTITUTION INVERSE --------------------
  229. C-----------------------------------------------------------------------
  230. C
  231.   800 CONTINUE
  232. C
  233. C-----------------------------------------------------------------------
  234. C------- ON TRANSFORME DANS LE SYSTEME DE COORDONNEES D ORIGINE --------
  235. C-----------------------------------------------------------------------
  236. C
  237.       DO     JJ = 1, N
  238.          J = N + 1 - JJ
  239. C
  240.          DO     I = 1, N
  241.             ZZ = 0.0D0
  242. C
  243.             DO     K = 1, J
  244.             ZZ = ZZ + MATZ.RE(I,K,1) * MATB.RE(K,J,1)
  245.             enddo
  246. C
  247.             MATZ.RE(I,J,1) = ZZ
  248.          enddo  
  249.       enddo  
  250. C
  251. C-----------------------------------------------------------------------
  252. C----------------- NORMALISATION DES VECTEURS PROPRES ------------------
  253. C-----------------------------------------------------------------------
  254. C
  255.       DO 950 J = 1, N
  256.          D = 0.0D0
  257.          IF ( ISW . EQ . 2 ) GOTO 920
  258.          SWAP = ALFI.PROG(J)
  259.          IF ( SWAP . NE . 0.0D0 ) GOTO 945
  260. C
  261.          DO 890 I = 1, N
  262.             IF ( ABS(MATZ.RE(I,J,1)) . GT . D )
  263.      &       D = ABS(MATZ.RE(I,J,1))
  264.   890    CONTINUE
  265. C
  266.          DO     I = 1, N
  267.            MATZ.RE(I,J,1) = MATZ.RE(I,J,1) / D
  268.          enddo
  269. C
  270.          GOTO 950
  271. C
  272.   920       continue
  273.          DO 930 I = 1, N
  274.             R = ABS(MATZ.RE(I,J-1,1)) + ABS(MATZ.RE(I,J,1))
  275.             IF ( R . NE . 0.0D0 )
  276.      &        R =     SQRT((MATZ.RE(I,J-1,1)  )**2
  277.      &          + (MATZ.RE(I,J,1)  )**2)
  278.             IF ( R . GT . D ) D = R
  279.   930    CONTINUE
  280. C
  281.          DO 940 I = 1, N
  282.             MATZ.RE(I,J-1,1) = MATZ.RE(I,J-1,1) / D
  283.             MATZ.RE(I,J,1) = MATZ.RE(I,J,1) / D
  284.   940    CONTINUE
  285. C
  286.   945    ISW = 3 - ISW
  287.   950 CONTINUE
  288. C
  289. c       SEGDES , MATA,MATB,MATZ
  290. c       SEGDES , ALFR,ALFI,BETA
  291.       SEGSUP , MATA,MATB
  292. c
  293. C=======================================================================
  294. C========================== FIN DES CALCULS ============================
  295. C=======================================================================
  296. C
  297.       RETURN
  298.       END
  299.  
  300.  
  301.  
  302.  
  303.  
  304.  
  305.  
  306.  
  307.  

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