Sources Esope | Cast3M

elpti7
C ELPTI7    SOURCE    CHAT      05/01/12    23:37:56     5004
        SUBROUTINE ELPTI7 (NMAX,NM,ENTIER)
      IMPLICIT INTEGER(I-N)
        IMPLICIT REAL*8 (A-H,O-Z)
-INC CCREEL
        SEGMENT STAB
           REAL*8 ICOS(INMAX), INTELL(INM,6)
        END SEGMENT
        LOGICAL ON_C
        INTEGER NMAX,NM
        INTEGER N,M,K,P,INT,ENTIER
        REAL*8 VALA,ALPHA,ALPHA0,THETA,VAL2,VAL3
        REAL*8 VAL4,VAL5,CO1,V1,V2,DV1,DV2,DK,L1
 
*        WRITE (6,*) '##discretisation integrale elliptique##'
        THETA= (45./NMAX)-90.
        ALPHA=0.
        ALPHA0=0.5*(1.-0.02)
        DK= (XPI/NMAX)
 
        INMAX=NMAX
        INM=NM+1
        SEGINI STAB
 
        DO 1 N=1,INMAX
        ICOS(N)= COS (THETA*XPI/180.)
        THETA=THETA+(90./INMAX)
  1     CONTINUE
 
        INTELL(1,1)=0.
        INTELL(1,2)=0.
 
        INT=1
  21    INT= INT + 1
        DO 2 M=INT,INM
        P=M-1
        VALA=0.
        ALPHA=0.5*P/(INM)
c        WRITE (6,*) 'ALPHA=',ALPHA
        INTELL(M,1)=ALPHA
 
              DO 3 K=1,INMAX
                VAL4=ABS(1.+(2.*ALPHA*ICOS(K)))
                VAL5= VAL4**0.5
                VALA=VALA-(ALPHA*ICOS(K)*XPI/(VAL5*INMAX))
  3           CONTINUE
 
c        WRITE (6,*) 'B1 M',P,'VALA',VALA
 
        K=0
        ON_C=.TRUE.
        DO WHILE (ON_C)
        K=K+1
        VAL2=ABS(1.-(2.*ALPHA*ICOS(K)))
        VAL3=VAL2**0.5
C       WRITE (6,*) 'VAL3',K,'= ',VAL3
        VALA=VALA +(ALPHA*ICOS(K)*XPI/(VAL3*INMAX))
                IF ((K.GE.INMAX).OR.((XPI/VAL3/INMAX).GE.(2*DK))) THEN
                ON_C=.FALSE.
                END IF
        END DO
c        WRITE (6,*) 'B2 M',P,'VALA',VALA,'K', K,(XPI/VAL3/INMAX),(2*DK)
 
        IF (K.GE.INMAX) THEN
        INTELL(M,2)=VALA
        GOTO 21
        END IF
 
        VAL2=0.
        V2=(K*XPI/2/INMAX)-(XPI/2)
        DV2= DK*((100.-(200.*ALPHA*(COS(V2))))**0.5)/10.
 
        ON_C=.TRUE.
        DO WHILE (ON_C)
        K=K+1
        V1=V2
        DV1=DV2
        V2=V1+(DV1/2.)
        CO1=COS(V2)
        VAL2=20.*ALPHA*CO1*DV1/((100.-(200.*ALPHA*CO1))**0.5)
        VALA=VALA+VAL2
        V2=V1+DV1
        DV2= DK*((100.-(200.*ALPHA*(COS(V2))))**0.5)/10.
                IF ((K.GE.6000).OR.(V2.GE.0.)) THEN
                 ON_C=.FALSE.
                END IF
        END DO
 
        L1=K
        VALA=VALA-VAL2
        CO1=COS((DV1-V2)/2)
        VAL2=20.*ALPHA*CO1*(DV1-V2)/((100.-(200.*ALPHA*CO1))**0.5)
        VALA=VALA+VAL2
        INTELL(M,2)=VALA
c        WRITE (6,*) 'FIN M',P,'VALA',VALA
  2     CONTINUE
 
c        WRITE (6,*) 'INTELL2= ',INTELL(INM,2),'L1= ', L1
        ENTIER=STAB
        SEGDES STAB
 
        RETURN
        END