1. Not finding help here? Sign up for a free 30min tutor trial with Chegg Tutors
    Dismiss Notice
Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Understanding fortran77 syntax

  1. Jun 27, 2015 #1
    1. The problem statement, all variables and given/known data
    I have two subroutine syntax program in fotran77 code, and I am stuck to understand it.
    this is the first subroutine code
    Code (Fortran):

    subroutine LZHES(N,A,NA,B,NB,X,NX,WANTX)
    c
    implicit real*8 (a-h,o-z)
    c
    complex*16 Y,W,Z,A(NA,N),B(NB,N),X(NX,N)
    real*8 C,D
    logical WANTX
    c
    NM1=N-1
    DO 30 I=1,NM1
    D=0.0
    IP1=I+1
    DO 10 K=IP1,N
    Y=B(K,I)
    C=dabs(dreal(Y))+dabs(dimag(Y))
    IF (C.LE.D) GOTO 10
    D=C
    II=K
    10  CONTINUE
    IF (D.EQ.0.0) GOTO 30
    Y=B(I,I)
    IF(D.LE.dabs(dreal(Y))+dabs(dimag(Y))) GOTO 15
    DO 11 J=1,N
    Y=A(I,J)
    A(I,J)=A(II,J)
    A(II,J)=Y
    11  CONTINUE
    DO 12 J=I,N
    Y=B(I,J)
    B(I,J)=B(II,J)
    B(II,J)=Y
    12  CONTINUE
    15  DO 20 J=IP1,N
    Y=B(J,I)/B(I,I)
    IF (dreal(Y).EQ.0.0.AND.dimag(Y).EQ.0.0) GOTO 20
    DO 18 K=1,N
    A(J,K)=A(J,K)-Y*A(I,K)
    18  CONTINUE
    DO 19 K=IP1,N
    B(J,K)=B(J,K)-Y*B(I,K)
    19  CONTINUE
    20  CONTINUE
    B(IP1,I)=(0.0,0.0)
    30  CONTINUE
    IF (.NOT.WANTX) GOTO 40
    DO 38 I=1,N
    DO 37 J=1,N
    X(I,J)=(0.0,0.0)
    37  CONTINUE
    X(I,I)=(1.0,0.0)
    38  CONTINUE
    40  NM2=N-2
    IF (NM2.LT.1) GOTO 100
    DO 90 J=1,NM2
    JM2=NM1-J
    JP1=J+1
    DO 80 II=1,JM2
    I=N+1-II
    IM1=I-1
    IMJ=I-J
    W=A(I,J)
    Z=A(IM1,J)
    IF (dabs(dreal(W))+dabs(dimag(W)).LE.dabs(dreal(Z))+
    &  dabs(dimag(Z))) GOTO 50
    DO 45 K=J,N
    Y=A(I,K)
    A(I,K)=A(IM1,K)
    A(IM1,K)=Y
    45  CONTINUE
    DO 46 K=IM1,N
    Y=B(I,K)
    B(I,K)=B(IM1,K)
    B(IM1,K)=Y
    46  CONTINUE
    50  Z=A(I,J)
    IF (dreal(Z).EQ.0.0.AND.dimag(Z).EQ.0.0) GOTO 58
    Y=Z/A(IM1,J)
    DO 52 K=JP1,N
    A(I,K)=A(I,K)-Y*A(IM1,K)
    52  CONTINUE
    DO 54 K=IM1,N
    B(I,K)=B(I,K)-Y*B(IM1,K)
    54  CONTINUE
    58  W=B(I,IM1)
    Z=B(I,I)
    IF(dabs(dreal(W))+dabs(dimag(W)).LE.dabs(dreal(Z))+
    &  dabs(dimag(Z))) GOTO 70
    DO 60 K=1,I
    Y=B(K,I)
    B(K,I)=B(K,IM1)
    B(K,IM1)=Y
    60  CONTINUE
    DO 64 K=1,N
    Y=A(K,I)
    A(K,I)=A(K,IM1)
    A(K,IM1)=Y
    64  CONTINUE
    IF (.NOT.WANTX) GOTO 70
    DO 68 K=IMJ,N
    Y=X(K,I)
    X(K,I)=X(K,IM1)
    X(K,IM1)=Y
    68  CONTINUE
    70  Z=B(I,IM1)
    IF (dreal(Z).EQ.0.0.AND.dimag(Z).EQ.0.0) GOTO 80
    Y=Z/B(I,I)
    DO 72 K=1,IM1
    B(K,IM1)=B(K,IM1)-Y*B(K,I)
    72  CONTINUE
    B(I,IM1)=(0.0,0.0)
    DO 74 K=1,N
    A(K,IM1)=A(K,IM1)-Y*A(K,I)
    74  CONTINUE
    IF (.NOT.WANTX) GOTO 80
    DO 76 K=IMJ,N
    X(K,IM1)=X(K,IM1)-Y*X(K,I)
    76  CONTINUE
    80  CONTINUE
    A(JP1+1,J)=(0.0,0.0)
    90  CONTINUE
    100  RETURN
    END
     
    and this is the second subroutine code
    Code (Fortran):

    subroutine LZIT(N,A,NA,B,NB,X,NX,WANTX,ITER,EIGA,EIGB)
    c
    implicit real*8 (a-h,o-z)
    c
    integer ITER(N)
    real*8 EPSA,EPSB,SS,R,ANORM,BNORM,ANI,BNI,C
    real*8 D0,D1,D2,E0,E1
    complex*16 A(NA,N),B(NB,N),EIGA(N),EIGB(N)
    complex*16 X(NX,N)
    complex*16 S,W,Y,Z
    complex*16 ANNM1,ALFM,BETM,D,SL,DEN,NUM,ANM1M1
    logical WANTX
    NN=N
    ANORM=0.
    BNORM=0.
    DO 5 I=1,N
    ANI=0.
    IF (I.EQ.1) GOTO 2
    Y=A(I,I-1)
    ANI=ANI+dabs(dreal(Y))+dabs(dimag(Y))
    2  BNI=0.
    DO 3 J=I,N
    ANI=ANI+dabs(dreal(A(I,J)))+dabs(dimag(A(I,J)))
    BNI=BNI+dabs(dreal(B(I,J)))+dabs(dimag(B(I,J)))
    3  CONTINUE
    IF (ANI.GT.ANORM) ANORM=ANI
    IF (BNI.GT.BNORM) BNORM=BNI
    5  CONTINUE
    IF (ANORM.EQ.0.) ANORM=1.0
    IF (BNORM.EQ.0.) BNORM=1.0
    EPSB=BNORM
    EPSA=ANORM
    6  EPSA=EPSA/2.0
    EPSB=EPSB/2.0
    C=ANORM+EPSA
    IF (C.GT.ANORM) GOTO 6
    IF (N.LE.1) GOTO 100
    10  ITS=0
    NM1=NN-1
    11  D2=dabs(dreal(A(NN,NN)))+dabs(dimag(A(NN,NN)))
    DO 12 LB=2,NN
    L=NN+2-LB
    SS=D2
    Y=A(L-1,L-1)
    D2=dabs(dreal(Y))+dabs(dimag(Y))
    SS=SS+D2
    Y=A(L,L-1)
    R=SS+dabs(dreal(Y))+dabs(dimag(Y))
    IF (R.EQ.SS) GOTO 13
    12  CONTINUE
    L=1
    13  IF (L.EQ.NN) GOTO 100
    IF (ITS.LT.30) GOTO 20
    ITER(NN)=-1
    IF (dabs(dreal(A(NN,NM1)))+dabs(dimag(A(NN,NM1))).GT.0.8*
    1  dabs(dreal(ANNM1))+dabs(dimag(ANNM1))) RETURN
    20  IF(ITS.EQ.10.OR.ITS.EQ.20) GOTO 28
    ANNM1=A(NN,NM1)
    ANM1M1=A(NM1,NM1)
    S=A(NN,NN)*B(NM1,NM1)-ANNM1*B(NM1,NN)
    W=ANNM1*B(NN,NN)*(A(NM1,NN)*B(NM1,NM1)-B(NM1,NN)*ANM1M1)
    Y=(ANM1M1*B(NN,NN)-S)/2.
    Z=CDSQRT(Y*Y+W)
    IF (dreal(Z).EQ.0.0.AND.dimag(Z).EQ.0.0) GOTO 26
    D0=Y/Z
    IF (D0.LT.0.0) Z=-Z
    26  DEN=(Y+Z)*B(NM1,NM1)*B(NN,NN)
    IF (dreal(DEN).EQ.0.0.AND.dimag(DEN).EQ.0.0) DEN=EPSA
    NUM=(Y+Z)*S-W
    GOTO 30
    28  Y=A(NM1,NN-2)
    NUM=dCMPLX(dabs(dreal(ANNM1))+dabs(dimag(ANNM1)),dabs(dreal(Y))+
    1dabs(dimag(Y)))
    DEN=(1.0,0.0)
    30  IF (NN.EQ.L+1) GOTO 35
    D2=dabs(dreal(A(NM1,NM1)))+dabs(dimag(A(NM1,NM1)))
    E1=dabs(dreal(ANNM1))+dabs(dimag(ANNM1))
    D1=dabs(dreal(A(NN,NN)))+dabs(dimag(A(NN,NN)))
    NL=NN-(L+1)
    DO 34 MB=1,NL
    M=NN-MB
    E0=E1
    Y=A(M,M-1)
    E1=dabs(dreal(Y))+dabs(dimag(Y))
    D0=D1
    D1=D2
    Y=A(M-1,M-1)
    D2=dabs(dreal(Y))+dabs(dimag(Y))
    Y=A(M,M)*DEN-B(M,M)*NUM
    D0=(D0+D1+D2)*(dabs(dreal(Y))+dabs(dimag(Y)))
    E0=E0*E1*(dabs(dreal(DEN))+dabs(dimag(DEN)))+D0
    IF(E0.EQ.D0) GOTO 36
    34  CONTINUE
    35  M=L
    36  CONTINUE
    ITS=ITS+1
    W=A(M,M)*DEN-B(M,M)*NUM
    Z=A(M+1,M)*DEN
    D1=dabs(dreal(Z))+dabs(dimag(Z))
    D2=dabs(dreal(W))+dabs(dimag(W))
    NP1=N+1
    LOR1=L
    NNORN=NN
    IF (.NOT.WANTX) GOTO 42
    LOR1=1
    NNORN=N
    42  DO 90 I=M,NM1
    J=I+1
    IF (I.EQ.M) GOTO 50
    W=A(I,I-1)
    Z=A(J,I-1)
    D1=dabs(dreal(Z))+dabs(dimag(Z))
    D2=dabs(dreal(W))+dabs(dimag(W))
    IF (D1.EQ.0.0) GOTO 11
    50  IF (D2.GT.D1) GOTO 60
    DO 55 K=I,NNORN
    Y=A(I,K)
    A(I,K)=A(J,K)
    A(J,K)=Y
    Y=B(I,K)
    B(I,K)=B(J,K)
    B(J,K)=Y
    55  CONTINUE
    IF (I.GT.M) A(I,I-1)=A(J,I-1)
    IF(D2.EQ.0.0) GOTO 65
    Y=dCMPLX(dreal(W)/D1,dimag(W)/D1)/dCMPLX(dreal(Z)/D1,
    &  dimag(Z)/D1)
    GOTO 62
    60  Y=dCMPLX(dreal(Z)/D2,dimag(Z)/D2)/dCMPLX(dreal(W)/D2,dimag(W)/D2)
    62   DO 64 K=I,NNORN
    A(J,K)=A(J,K)-Y*A(I,K)
    B(J,K)=B(J,K)-Y*B(I,K)
    64  CONTINUE
    IF (I.GT.M) A(J,I-1)=(0.0,0.0)
    65  Z=B(J,I)
    W=B(J,J)
    D2=dabs(dreal(W))+dabs(dimag(W))
    D1=dabs(dreal(Z))+dabs(dimag(Z))
    IF (D1.EQ.0.0) GOTO 11
    IF (D2.GT.D1) GOTO 80
    DO 70 K=LOR1,J
    Y=A(K,J)
    A(K,J)=A(K,I)
    A(K,I)=Y
    Y=B(K,J)
    B(K,J)=B(K,I)
    B(K,I)=Y
    70  CONTINUE
    IF(I.EQ.NM1) GOTO 75
    Y=A(J+1,J)
    A(J+1,J)=A(J+1,I)
    A(J+1,I)=Y
    75  IF(.NOT.WANTX) GOTO 79
    DO 78 K=1,N
    Y=X(K,J)
    X(K,J)=X(K,I)
    X(K,I)=Y
    78  CONTINUE
    79  IF(D2.EQ.0.0) GOTO 90
    Z=dCMPLX(dreal(W)/D1,dimag(W)/D1)/dCMPLX(dreal(Z)/D1,dimag(Z)/D1)
    GOTO 81
    80  Z=dCMPLX(dreal(Z)/D2,dimag(Z)/D2)/dCMPLX(dreal(W)/D2,dimag(W)/D2)
    81  DO 82 K=LOR1,J
    A(K,I)=A(K,I)-Z*A(K,J)
    B(K,I)=B(K,I)-Z*B(K,J)
    82  CONTINUE
    B(J,I)=(0.0,0.0)
    IF (I.LT.NM1) A(I+2,I)=A(I+2,I)-Z*A(I+2,J)
    IF(.NOT.WANTX) GOTO 90
    DO 85 K=1,N
    X(K,I)=X(K,I)-Z*X(K,J)
    85  CONTINUE
    90  CONTINUE
    GOTO 11
    100  CONTINUE
    EIGA(NN)=A(NN,NN)
    EIGB(NN)=B(NN,NN)
    IF (NN.EQ.1) GOTO 110
    ITER(NN)=ITS
    NN=NM1
    IF (NN.GT.1) GOTO 10
    ITER(1)=0
    GOTO 100
    110  IF(.NOT.WANTX) RETURN
    M=N
    115  CONTINUE
    ALFM=A(M,M)
    BETM=B(M,M)
    B(M,M)=(1.0,0.0)
    L=M-1
    IF (L.EQ.0) GOTO 140
    120  CONTINUE
    L1=L+1
    SL=0.
    DO 130 J=L1,M
    SL=SL+(BETM*A(L,J)-ALFM*B(L,J))*B(J,M)
    130  CONTINUE
    Y=BETM*A(L,L)-ALFM*B(L,L)
    IF(dreal(Y).EQ.0.0.AND.dimag(Y).EQ.0.0) Y=(EPSA+EPSB)/2.0
    B(L,M)=-SL/Y
    L=L-1
    140  IF (L.GT.0) GOTO 120
    M=M-1
    IF (M.GT.0) GOTO 115
    M=N
    200  CONTINUE
    DO 220 I=1,N
    S=0.
    DO 210 J=1,M
    S=S+X(I,J)*B(J,M)
    210  CONTINUE
    X(I,M)=S
    220  CONTINUE
    M=M-1
    IF (M.GT.0) GOTO 200
    M=N
    230  CONTINUE
    SS=0.
    DO 235 I=1,N
    R=dabs(dreal(X(I,M)))+dabs(dimag(X(I,M)))
    IF (R.LT.SS) GOTO 235
    SS=R
    D=X(I,M)
    235  CONTINUE
    IF(SS.EQ.0.0) GOTO 245
    DO 240 I=1,N
    X(I,M)=X(I,M)/D
    240  CONTINUE
    245  M=M-1
    IF (M.GT.0) GOTO 230
    RETURN
    END
     
    and the question is
    I don't know exactly the subroutine doing what,
    and what the subroutine output and input

    2. Relevant equations

    3. The attempt at a solution

    I thought the first sub routine (LZHES) trying to diagonalization matrix B
    and I thought the first subroutine input is
    N,A,NA,B,NB,X,NX,WANTX
    N is the size of matrix
    A is the identity matrix
    B is the matrix will be diagonalize
    NA, NB is the size of matrix A and matrix B
    WANTX is the variable to decided to diagonalize matrix or not.
    and the subroutine output, I thought is B, the diagonalize matrix B.
    but this subroutine when I tried to created in MATLAB code version, matrix B is not diagonalized.
    and I thought this subroutine LZHES is connected to the subroutine LZIT.

    and the second subroutine trying to find the eigen values of the matrix B.
    I thought the first subroutine input is
    N,A,NA,B,NB,X,NX,WANTX,ITER,EIGA,EIGB
    N is the size of matrix
    A is the identity matrix
    B is the matrix will be diagonalize
    NA, NB is the size of matrix A and matrix B
    WANTX is the variable to decided to diagonalize matrix or not.
    ITER is the number of the iteration
    EIGA is the initialize eigen value
    EIGB is the initilaize eigen value

    and the output, i thought is the eigen values which is EIGA and EIGB
     
    Last edited by a moderator: Jun 27, 2015
  2. jcsd
  3. Jun 27, 2015 #2

    SteamKing

    User Avatar
    Staff Emeritus
    Science Advisor
    Homework Helper

    An article was published in 1975 in the ACM journal Transactions on Mathematical Software (TOMS), as cited here:

    http://dl.acm.org/citation.cfm?doid=355644.355652

    The original Fortran source code (with comments) for these two subroutines LZHES and LZIT can be recovered from the ACM Collected Algorithms page located here:

    http://netlib.org/toms/

    The subroutines are in a zip file named 496.gz on that page, which you can download for free and compare with the code from the OP.

    It appears that someone copied these subroutines and clumsily stripped out any comments which would serve to inform users of these routines what their purpose and function were.

    As to the other details of these algorithms, you'll have to track down the accompanying article from the TOMS journal from 1975. This journal should be available in most university libraries, or a .pdf copy can be purchased directly from the ACM using the citation above. There may be a copy floating around the internet somewhere, but I haven't looked.

    Good Luck!
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted



Similar Discussions: Understanding fortran77 syntax
  1. Python syntax (Replies: 5)

  2. Matlab syntax (Replies: 2)

Loading...