File DTLAP.FT (FORTRAN source file)

Directory of image this file is from
This file as a plain text file

C
C     ..................................................................
C
C        SUBROUTINE DTLAP
C
C        PURPOSE
C           A SERIES EXPANSION IN LAGUERRE POLYNOMIALS WITH INDEPENDENT
C           VARIABLE X IS TRANSFORMED TO A POLYNOMIAL WITH INDEPENDENT
C           VARIABLE Z, WHERE X=A*Z+B
C
C        USAGE
C           CALL DTLAP(A,B,POL,N,C,WORK)
C
C        DESCRIPTION OF PARAMETERS
C           A     - FACTOR OF LINEAR TERM IN GIVEN LINEAR TRANSFORMATION
C                   DOUBLE PRECISION VARIABLE
C           B     - CONSTANT TERM IN GIVEN LINEAR TRANSFORMATION
C                   DOUBLE PRECISION VARIABLE
C           POL   - COEFFICIENT VECTOR OF POLYNOMIAL (RESULTANT VALUE)
C                   COEFFICIENTS ARE ORDERED FROM LOW TO HIGH
C                   DOUBLE PRECISION VECTOR
C           N     - DIMENSION OF COEFFICIENT VECTORS POL AND C
C           C     - GIVEN COEFFICIENT VECTOR OF EXPANSION
C                   COEFFICIENTS ARE ORDERED FROM LOW TO HIGH
C                   POL AND C MAY BE IDENTICALLY LOCATED
C                   DOUBLE PRECISION VECTOR
C           WORK  - WORKING STORAGE OF DIMENSION 2*N
C                   DOUBLE PRECISION ARRAY
C
C        REMARKS
C           COEFFICIENT VECTOR C REMAINS UNCHANGED IF NOT COINCIDING
C           WITH COEFFICIENT VECTOR POL.
C           OPERATION IS BYPASSED IN CASE N LESS THAN 1.
C           THE LINEAR TRANSFORMATION X=A*Z+B OR Z=(1/A)(X-B) TRANSFORMS
C           THE RANGE (0,C) IN X TO THE RANGE (ZL,ZR) IN Z, WHERE
C           ZL=-B/A AND ZR=(C-B)/A.
C           FOR GIVEN ZL, ZR AND C WE HAVE A=C/(ZR-ZL) AND
C           B=-C*ZL/(ZR-ZL)
C
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C           NONE
C
C        METHOD
C           THE TRANSFORMATION IS BASED ON THE RECURRENCE EQUATION
C           FOR LAGUERRE POLYNOMIALS L(N,X)
C           L(N+1,X)=2*L(N,X)-L(N-1,X)-((1+X)*L(N,X)-L(N-1,X))/(N+1),
C           WHERE THE FIRST TERM IN BRACKETS IS THE INDEX,
C           THE SECOND IS THE ARGUMENT.
C           STARTING VALUES ARE L(0,X)=1, L(1,X)=1-X.
C           THE TRANSFORMATION IS IMPLICITLY DEFINED BY MEANS OF
C           X=A*Z+B TOGETHER WITH
C           SUM(POL(I)*Z**(I-1), SUMMED OVER I FROM 1 TO N)
C           =SUM(C(I)*L(I-1,X), SUMMED OVER I FROM 1 TO N).
C
C     ..................................................................
C
      SUBROUTINE DTLAP(A,B,POL,N,C,WORK)
C
      DIMENSION POL(1),C(1),WORK(1)
      DOUBLE PRECISION A,B,POL,C,WORK,H,P,Q,Q1,Q2,FI
C
C        TEST OF DIMENSION
      IF(N-1)2,1,3
C
C        DIMENSION LESS THAN 2
    1 POL(1)=C(1)
    2 RETURN
C
    3 POL(1)=C(1)+C(2)-B*C(2)
      POL(2)=-C(2)*A
      IF(N-2)2,2,4
C
C        INITIALIZATION
    4 WORK(1)=1.D0
      WORK(2)=1.D0-B
      WORK(3)=0.D0
      WORK(4)=-A
      FI=1.D0
C
C        CALCULATE COEFFICIENT VECTOR OF NEXT LAGUERRE POLYNOMIAL
C        AND ADD MULTIPLE OF THIS VECTOR TO POLYNOMIAL POL
      DO 6 J=3,N
      FI=FI+1.D0
      Q=1.D0/FI
      Q1=Q-1.D0
      Q2=1.D0-Q1-B*Q
      Q=Q*A
      P=0.D0
C
      DO 5 K=2,J
      H=-P*Q+WORK(2*K-2)*Q2+WORK(2*K-3)*Q1
      P=WORK(2*K-2)
      WORK(2*K-2)=H
      WORK(2*K-3)=P
    5 POL(K-1)=POL(K-1)+H*C(J)
      WORK(2*J-1)=0.D0
      WORK(2*J)=-Q*P
    6 POL(J)=C(J)*WORK(2*J)
      RETURN
      END



Feel free to contact me, David Gesswein djg@pdp8online.com with any questions, comments on the web site, or if you have related equipment, documentation, software etc. you are willing to part with.  I am interested in anything PDP-8 related, computers, peripherals used with them, DEC or third party, or documentation. 

PDP-8 Home Page   PDP-8 Site Map   PDP-8 Site Search