C RAYLEI    SOURCE    CHAT      09/10/09    21:22:23     6519           
      SUBROUTINE RAYLEI (IPA,IPB,IV,IW,XRVP,XIVP,IPLMOX,IPLMOY,DELTA
     &      ,IPRX, IPIX)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H, O-Z)

*****************************************************************************
*                                                                           *
*                               RAYLEI                                      *
*                               ------                                      *
*****************************************************************************
*
*       Fonction:
*
*       Trouve les valeurs propres d'une matrice de Rayleigh.
*
*****************************************************************************
*
*      
*       Equation caracteristique: (1, -(LL11+LL22), (LL11*LL22-LL12*LL21))
*       Discriminant: Delta
*
*
*       PROCEDURE D APPEL: CALL RAYLEI (IPA, IPB, IV, IW, XRVP, XIVP,IPLMOX,
*                                      IPLMOY,DELTA, IPRX, IPIX)
*       
-INC CCOPTIO
-INC SMLMOTS
-INC SMCHPOI
-INC SMRIGID
-INC CCHAMP
        SEGMENT IDEMEN (0)

        REAL*8 LL(1:2,1:2)
        REAL*8 A, B, C, D, TO, TI, XRVP, XIVP
*       COORDONNEES DE LA MATRICE
        MCHPOI = 0
        MCHPO2 = 0
        CALL MUCPRI (IV, IPB, IPBV)
        SEGINI IDEMEN                                                             
        IDEMEN(**)=IPBV                                                       
        NOID=0                                                                    
        NOEN=1                                                                   
        CALL LDMT   (IPA,IDEMEN,NOID,NOEN,1.D-15)
        SEGACT IDEMEN                                   
        MCHPOI=IDEMEN(1)
        SEGSUP IDEMEN                                                          
        IF (IERR .NE. 0) RETURN
        CALL XTY1 (IV, MCHPOI, IPLMOX, IPLMOX, LL(1,1) )
        
        CALL MUCPRI (IW, IPB, IPBW)
        SEGINI IDEMEN                                                            
        IDEMEN(**)=IPBW                                                          
        NOID=0                                                                    
        NOEN=1                                                                    
        CALL LDMT   (IPA,IDEMEN,NOID,NOEN,1.D-15)
        SEGACT IDEMEN                                                             
        MCHPO2=IDEMEN(1)                                                           
        SEGSUP IDEMEN                                                             
        IF (IERR .NE. 0) RETURN
        CALL XTY1 (IV, MCHPO2, IPLMOX, IPLMOX, LL(1,2) )
        CALL XTY1 (IW, MCHPOI, IPLMOX, IPLMOX, LL(2,1) )
        CALL XTY1 (IW, MCHPO2, IPLMOX, IPLMOX, LL(2,2) )
        
        A=LL(1,1)
        B=LL(1,2)
        C=LL(2,1)
        D=LL(2,2)
        TO=(A+D)/2.
        
*       Discriminant
        DELTA=(A+D)**2-4.D0*(A*D-B*C)
        TI=ABS(DELTA)
        IF (TI .GE. 1.D-14) THEN 
                IF (DELTA .LT. 0.D0) THEN
                        GOTO 100
                ELSE
                    RETURN
                END IF
        ELSE
            GOTO 200
        END IF
        
        
        
 100    CONTINUE
        
        XRVP=TO
        XIVP=SQRT(TI)/2.D0
*
*     EXTRACTION DES VECTEURS PROPRES
*
        IF ( ABS( B ) .GE. ABS( C )) THEN
           X2R = (XRVP - A ) / B
           X2C = XIVP / B
           CALL ADCHPO (IV, IW, IPRX, 1.D0, X2R)
           CALL MUCHPO (IW, X2C, IPIX, 1)
        ELSE
           X1R = (XRVP - D ) / C
           X1C = XIVP / C
           CALL ADCHPO (IV, IW, IPRX, X1R, 1.D0)
           CALL MUCHPO (IV, X1C, IPIX, 1)
        END IF
*
*      COEFFICIENT DE RAYLEIGH
*      EXTRACTION DES VALEURS PROPRES
*      
        CALL XTMX(IPRX, IPB, RXBRX)
        CALL XTMX(IPIX, IPB, RIXBRIX)
        ALPHA = RXBRX + RIXBRIX
        CALL XTMX(IPRX, IPA, RXARX)
        CALL XTMX(IPIX, IPA, RIXARIX)
        BETA = RXARX + RIXARIX
        XRVP = BETA / ALPHA
        CALL MUCPRI(IPIX, IPA, IPAXI)
        CALL XTY1(IPRX, IPAXI, IPLMOX, IPLMOY, XRKXI)
        CALL MUCPRI(IPRX, IPA, IPAXR)
        CALL XTY1(IPIX, IPAXR, IPLMOX, IPLMOY, XIKXR)
        XIVP = (XIKXR - XRKXI) / ALPHA
        GOTO 300
 200    CONTINUE
        IPRX = IV
        CALL MUCPRI(IV, IPA, IPAX2)
        CALL XTY1(IV, IPAX2, IPLMOX, IPLMOY, X2AX2)
        IF (IERR .NE. 0 ) RETURN
        CALL MUCPRI(IV, IPB, IPBX2)
        CALL XTY1(IV, IPBX2, IPLMOX, IPLMOY, X2BX2)
        IF (IERR .NE. 0 ) RETURN
        XRVP = X2AX2 / X2BX2
        XIVP = 0.D0 
 300    CONTINUE   
        END