Esope source files

   1 : C SJJCOI    SOURCE    BP208322  09/09/09    21:15:10     6495           
   2 : C
   3 :         SUBROUTINE SJjCOi(A,b,c,N,th)
   4 : C
   5 :       IMPLICIT REAL*8(A-H,o-Z)
   6 :       IMPLICIT INTEGER(I-N)
   7 :       dimension a(N,*),b(n,*),c(n,*)
   8 : C 
   9 : C
  10 : *        SEGACT ,IPK1*MOD
  11 : C
  12 :        do 200 j=1,n-1
  13 :         do 250 i=j+1,n
  14 :        CF1=SQRT (ABS(A(i,j)**2/(A(i,i)*A(j,j))))
  15 :        CF2=SQRT (ABS(b(i,j)**2/(b(i,i)*b(j,j))))
  16 : 
  17 :        IF((CF1.GT.TH).OR.(CF2.GT.TH)) THEN
  18 : 
  19 :      
  20 :            XKBII=A(I,I)*B(I,J)-B(I,I)*A(I,J)
  21 :            XKBJJ=A(J,J)*B(I,J)-B(J,J)*A(I,J)
  22 :            XKB=A(I,I)*B(J,J)-B(I,I)*A(J,J)
  23 :            IF (XKB.NE.0.D0) THEN
  24 :              X=XKB/2+(XKB/ABS(XKB))*SQRT(ABS((XKB/2)**2+XKBII*XKBJJ))
  25 :            ELSE
  26 :              X=SQRT(XKBII*XKBJJ)
  27 :            ENDIF
  28 : *          bp (septembre 2009): on ajoute le cas ou X = 0
  29 :            IF (X.NE.0.D0) THEN
  30 :              GAMMA=-XKBII/X
  31 :              ALPHA=XKBJJ/X
  32 :            ELSE
  33 :              GAMMA=-(A(I,J)/A(J,J))
  34 :              ALPHA=0.D0
  35 :            ENDIF
  36 : 
  37 :         DO 100 L=1,N
  38 : C
  39 :          xK1= A(I,L)
  40 :          xK2=A(J,L)
  41 :          A(I,L)=xK1+GAMMA*xK2
  42 :          A(J,L)=ALPHA*xK1+xK2
  43 : C             
  44 : 100     CONTINUE
  45 : C        
  46 :         DO 150 L=1,N
  47 :           xK1=A(L,I)
  48 :           xK2=A(L,J)
  49 :           A(L,I)=xK1+GAMMA*xK2
  50 :           A(L,J)=ALPHA*xK1+xK2    
  51 : 150     CONTINUE
  52 :         DO 101 L=1,N
  53 : C
  54 :          xK1= b(I,L)
  55 :          xK2=b(J,L)
  56 :          b(I,L)=xK1+GAMMA*xK2
  57 :          b(J,L)=ALPHA*xK1+xK2
  58 : C             
  59 : 101     CONTINUE
  60 : C        
  61 :         DO 151 L=1,N
  62 :           xK1=b(L,I)
  63 :           xK2=b(L,J)
  64 :           b(L,I)=xK1+GAMMA*xK2
  65 :           b(L,J)=ALPHA*xK1+xK2    
  66 : 151     CONTINUE
  67 : 
  68 :          DO 102 L=1,N
  69 :            xK1=c(L,I)
  70 :            xK2=c(L,J)
  71 :            c(L,I)=xK1+GAMMA*xK2
  72 :            c(L,J)=ALPHA*xK1+xK2    
  73 : 102     CONTINUE
  74 :         endif
  75 : 
  76 :  250  continue
  77 :  200  continue
  78 : *        SEGDES ,IPK1 
  79 :         RETURN
  80 :         END
  81 :  
  82 :  
  83 :  
  84 :