!****************************************************
!*  Program to demonstrate Cheby_Ser() Subroutine   *
!* ------------------------------------------------ *
!* Reference: BASIC Scientific Subroutines, Vol. II *
!* By F.R. Ruckdeschel, BYTE/McGRAWW-HILL, 1981 [1].*
!*                                                  *
!*             F90 Version by J.-P. Moreau, Paris.  *
!*                     (www.jpmoreau.fr)            *
!* ------------------------------------------------ *
!* SAMPLE RUN:                                      *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 2   *
!*    A(0) = -1                                     *
!*    A(1) = 0                                      *
!*    A(2) = 2                                      *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 3   *
!*    A(0) = 0                                      *
!*    A(1) = -3                                     *
!*    A(2) = 0                                      *
!*    A(3) = 4                                      *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 4   *
!*    A(0) = 1                                      *
!*    A(1) = 0                                      *
!*    A(2) = -8                                     *
!*    A(3) = 0                                      *
!*    A(4) = 8                                      *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 5   *
!*    A(0) = 0                                      *
!*    A(1) = 5                                      *
!*    A(2) = 0                                      *
!*    A(3) = -20                                    *
!*    A(4) = 0                                      *
!*    A(5) = 16                                     *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 6   *
!*    A(0) = -1                                     *
!*    A(1) = 0                                      *
!*    A(2) = 18                                     *
!*    A(3) = 0                                      *
!*    A(4) = -48                                    *
!*    A(5) = 0                                      *
!*    A(6) = 32                                     *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 7   *
!*    A(0) = 0                                      *
!*    A(1) = -7                                     *
!*    A(2) = 0                                      *
!*    A(3) = 56                                     *
!*    A(4) = 0                                      *
!*    A(5) = -112                                   *
!*    A(6) = 0                                      *
!*    A(7) = 64                                     *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 8   *
!*    A(0) = 1                                      *
!*    A(1) = 0                                      *
!*    A(2) = -32                                    *
!*    A(3) = 0                                      *
!*    A(4) = 160                                    *
!*    A(5) = 0                                      *
!*    A(6) = -256                                   *
!*    A(7) = 0                                      *
!*    A(8) = 128                                    *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 9   *
!*    A(0) = 0                                      *
!*    A(1) = 9                                      *
!*    A(2) = 0                                      *
!*    A(3) = -120                                   *
!*    A(4) = 0                                      *
!*    A(5) = 432                                    *
!*    A(6) = 0                                      *
!*    A(7) = -576                                   *
!*    A(8) = 0                                      *
!*    A(9) = 256                                    *
!*                                                  *
!* Chebyshev polynomial coefficients for degree 10  *
!*    A(0) = -1                                     *
!*    A(1) = 0                                      *
!*    A(2) = 50                                     *
!*    A(3) = 0                                      *
!*    A(4) = -400                                   *
!*    A(5) = 0                                      *
!*    A(6) = 1120                                   *
!*    A(7) = 0                                      *
!*    A(8) = -1280                                  *
!*    A(9) = 0                                      *
!*    A(10) = 512                                   *
!*                                                  *
!****************************************************
PROGRAM CHEBYSER

parameter(SIZE=10)

real*8 B(SIZE+1,SIZE+1)
integer  i, n

   print *, ' '
   do n = 2, SIZE
     call Cheby_Ser(n,B)
     write(*,25) n
     do i = 0, n
       if (dabs(B(n,i)) < 1.d-10) B(n,i) = 0.d0
       write(*,50) i, B(n,i)
     end do   
     print *, ' '
     if (n < SIZE) pause ' <Enter> to continue...'
   end do
   print *, ' '

   stop
25 format(' Chebyshev polynomial coefficients for degree ',i2/)
50 format('   A(',i2,') = ',f8.0) 

end


!********************************************************
!* Chebyshev series coefficients evaluation subroutine  *
!* ---------------------------------------------------- *
!* The order of the polynomial is n. The coefficients   *
!* are returned in the array B[i,j), i is the degree of *
!* the polynomial, j is the coefficient order.          *
!********************************************************
Subroutine Cheby_Ser(n,B)
  parameter(SIZE=10)  
  integer i, j, n
  real*8 B(SIZE+1,SIZE+1)
  ! Establish t0 and t1 coefficients
  B(0,0) = 1; B(1,0) = 0; B(1,1) = 1
  ! Return if order is less than two
  if (n < 2) return
  do i = 2, n
    do j = 1, i
      ! Basic recursion relation
      B(i,j) = 2.d0 * B(i-1,j-1) - B(i-2,j)
    end do
    B(i,0) = -B(i-2,0)
  end do
  return
end

! End of file chebyser.f90