Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Legendre poly in Fotran 90

  1. Feb 27, 2009 #1
    Is it possible to evaluate a legendre polynomial p(n,x) using the recurrence relation

    p(n,x)=p(n-1,x)*x-p(n-2,x) in fotran 90


    {there are some other terms which i left out for brevity]
     
  2. jcsd
  3. Feb 28, 2009 #2
    Here is an example in F77 for N up to 100.
    Recursion is not necessary when the initial conditions are known, i.e. P(0) and P(1).
    P(2) was defined to avoid going back to P(0) which is not permitted in F77.

    Code (Text):
          FUNCTION FLEGENDRE(N,X)
          REAL*8 P(100)
          REAL*8 FI
          INTEGER I
    C     CHECK FOR VALID VALUES OF N AND X HERE BEFORE PROCEEDING
    C
    C     PROCEEDING WITH CALCULATIONS
          IF (N.EQ.0) THEN
            FLEGENDRE=1
          ELSEIF (N.EQ.1) THEN
            FLEGENDRE=X
          ELSEIF (N.EQ.2) THEN
            FLEGENDRE=(3*X*X-1)/2.0
          ELSE
    C       SYNTHETIC CALCULATIONS
            P(1)=X
            P(2)=(3*X*X-1)/2.0
            DO 20 I=3,N
            FI=I
            P(I)=((I+I-1)*X*P(I-1)-(I-1)*P(I-2))/FI
       20   CONTINUE
            FLEGENDRE=P(N)
          ENDIF
          END
    C  
    C     TEST PROGRAM
    C
          REAL*8 OUT
          DO 30 I=1,10
          OUT=FLEGENDRE(I,5.425)
          WRITE(6,999)I,OUT
       30 CONTINUE
      999 FORMAT(I3,F20.2)
          STOP
          END
     
     
  4. Feb 28, 2009 #3
    Here a computationally more efficient version that obviates the use of arrays and hence no limit on the size of N except for round-off errors.

    Code (Text):
          FUNCTION FLEGENDRE(N,X)
          REAL*8 PI,PIM1,PIM2
          REAL*8 FI
          INTEGER I
    C     CHECK FOR VALID VALUES OF N AND X HERE BEFORE PROCEEDING
    C
    C     PROCEEDING WITH CALCULATIONS
          IF (N.EQ.0) THEN
            FLEGENDRE=1
          ELSEIF (N.EQ.1) THEN
            FLEGENDRE=X
          ELSE
    C       SYNTHETIC CALCULATIONS
            PIM1=1
            PI=X
            DO 20 I=2,N
            FI=I
            PIM2=PIM1
            PIM1=PI
            PI=((I+I-1)*X*PIM1-(I-1)*PIM2)/FI
       20   CONTINUE
            FLEGENDRE=PI
          ENDIF
          END
    C  
    C     TEST PROGRAM
    C
          REAL*8 OUT
          DO 30 I=1,10
          OUT=FLEGENDRE(I,5.425)
          WRITE(6,999)I,OUT
       30 CONTINUE
      999 FORMAT(I3,F20.2)
          STOP
          END
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook




Similar Discussions: Legendre poly in Fotran 90
  1. Fortran 90 (Replies: 1)

Loading...