program c2pack_test2d

  use ctpack
  implicit none
  integer, parameter :: im=8, km=5, nm=4
  real(8), parameter :: xmin=-1.0D0, xmax=1.0D0
  real(8), parameter :: eps=1.0D-14

  real(8), dimension(0:nm,0:im) :: y
  real(8), dimension(0:nm,0:km) :: t_y
  real(8), dimension(0:nm,0:im) :: yderiv

  integer :: k, i, m

  call ctinitial(im,km,xmin,xmax)

  y(0,:) = 1                          ! T_0
  y(1,:) = g_x                        ! T_1
  y(2,:) = 2*g_x**2 - 1               ! T_2
  y(3,:) = 4*g_x**3 - 3*g_x           ! T_3
  y(4,:) = 8*g_x**4 - 8*g_x**2 + 1
  yderiv(0,:) = 0.0D0                 ! T_0'= 0
  yderiv(1,:) = 1.0D0                 ! T_1'= T_0
  yderiv(2,:) = 4*g_x                 ! T_2'= 4*T_1
  yderiv(3,:) = 12*g_x**2 - 3         ! T_3'= 6*T_2 + 3*T_0
  yderiv(4,:) = 32*g_x**3 - 16*g_x    ! T_4'= 8*T_3 + 8*T_1

  t_y = t_dx_t(t_g(y))

  write(6,*) 'X grid points'
  write(6,*) g_x
  write(6,*)

  do m=0,nm
     write(6,*) 'n=',m
     write(6,*) t_y(m,:)
  enddo

  y=g_t(t_y)
  do m=0,nm
     do i=0,im
        if ( abs(y(m,i)-yderiv(m,i)) .gt. eps ) then
           write(6,*) m,i, y(m,i), yderiv(m,i), y(m,i)-yderiv(m,i)
        endif
     enddo
  enddo

end program c2pack_test2d
