dqcheb Subroutine

public subroutine dqcheb(x, Fval, Cheb12, Cheb24)

chebyshev series expansion

this routine computes the chebyshev series expansion of degrees 12 and 24 of a function using a fast fourier transform method

  • f(x) = sum(k=1,..,13) (cheb12(k)*t(k-1,x))
  • f(x) = sum(k=1,..,25) (cheb24(k)*t(k-1,x))

where t(k,x) is the chebyshev polynomial of degree k.

See also

History

  • QUADPACK: revision date 830518 (yymmdd)

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: x(11)

vector of dimension 11 containing the values cos(k*pi/24), k = 1, ..., 11
real(kind=wp), intent(inout) :: Fval(25)

vector of dimension 25 containing the function values at the points (b+a+(b-a)*cos(k*pi/24))/2, k = 0, ...,24, where (a,b) is the approximation interval. fval(1) and fval(25) are divided by two (these values are destroyed at output).
real(kind=wp), intent(out) :: Cheb12(13)

vector of dimension 13 containing the chebyshev coefficients for degree 12
real(kind=wp), intent(out) :: Cheb24(25)

vector of dimension 25 containing the chebyshev coefficients for degree 24

Called by

proc~~dqcheb~~CalledByGraph proc~dqcheb quadpack_generic::dqcheb proc~dqc25c quadpack_generic::dqc25c proc~dqc25c->proc~dqcheb proc~dqc25f quadpack_generic::dqc25f proc~dqc25f->proc~dqcheb proc~dqc25s quadpack_generic::dqc25s proc~dqc25s->proc~dqcheb proc~dqawce quadpack_generic::dqawce proc~dqawce->proc~dqc25c proc~dqawoe quadpack_generic::dqawoe proc~dqawoe->proc~dqc25f proc~dqawse quadpack_generic::dqawse proc~dqawse->proc~dqc25s proc~dqawc quadpack_generic::dqawc proc~dqawc->proc~dqawce proc~dqawfe quadpack_generic::dqawfe proc~dqawfe->proc~dqawoe proc~dqawo quadpack_generic::dqawo proc~dqawo->proc~dqawoe proc~dqaws quadpack_generic::dqaws proc~dqaws->proc~dqawse proc~dqawf quadpack_generic::dqawf proc~dqawf->proc~dqawfe

Source Code

    subroutine dqcheb(x, Fval, Cheb12, Cheb24)
        implicit none

        real(wp), intent(in) :: x(11) !! vector of dimension 11 containing the
                                      !! values `cos(k*pi/24), k = 1, ..., 11`
        real(wp), intent(inout) :: Fval(25) !! vector of dimension 25 containing the
                                            !! function values at the points
                                            !! `(b+a+(b-a)*cos(k*pi/24))/2, k = 0, ...,24`,
                                            !! where `(a,b)` is the approximation interval.
                                            !! `fval(1)` and `fval(25)` are divided by two
                                            !! (these values are destroyed at output).
        real(wp), intent(out) :: Cheb12(13) !! vector of dimension 13 containing the
                                            !! chebyshev coefficients for degree 12
        real(wp), intent(out) :: Cheb24(25) !! vector of dimension 25 containing the
                                            !! chebyshev coefficients for degree 24

        real(wp) :: alam, alam1, alam2, part1, part2, part3, v(12)
        integer :: i, j

        do i = 1, 12
            j = 26 - i
            v(i) = Fval(i) - Fval(j)
            Fval(i) = Fval(i) + Fval(j)
        end do
        alam1 = v(1) - v(9)
        alam2 = x(6)*(v(3) - v(7) - v(11))
        Cheb12(4) = alam1 + alam2
        Cheb12(10) = alam1 - alam2
        alam1 = v(2) - v(8) - v(10)
        alam2 = v(4) - v(6) - v(12)
        alam = x(3)*alam1 + x(9)*alam2
        Cheb24(4) = Cheb12(4) + alam
        Cheb24(22) = Cheb12(4) - alam
        alam = x(9)*alam1 - x(3)*alam2
        Cheb24(10) = Cheb12(10) + alam
        Cheb24(16) = Cheb12(10) - alam
        part1 = x(4)*v(5)
        part2 = x(8)*v(9)
        part3 = x(6)*v(7)
        alam1 = v(1) + part1 + part2
        alam2 = x(2)*v(3) + part3 + x(10)*v(11)
        Cheb12(2) = alam1 + alam2
        Cheb12(12) = alam1 - alam2
        alam = x(1)*v(2) + x(3)*v(4) + x(5)*v(6) + x(7)*v(8) + x(9)*v(10) &
               + x(11)*v(12)
        Cheb24(2) = Cheb12(2) + alam
        Cheb24(24) = Cheb12(2) - alam
        alam = x(11)*v(2) - x(9)*v(4) + x(7)*v(6) - x(5)*v(8) + x(3)*v(10) &
               - x(1)*v(12)
        Cheb24(12) = Cheb12(12) + alam
        Cheb24(14) = Cheb12(12) - alam
        alam1 = v(1) - part1 + part2
        alam2 = x(10)*v(3) - part3 + x(2)*v(11)
        Cheb12(6) = alam1 + alam2
        Cheb12(8) = alam1 - alam2
        alam = x(5)*v(2) - x(9)*v(4) - x(1)*v(6) - x(11)*v(8) + x(3)*v(10) &
               + x(7)*v(12)
        Cheb24(6) = Cheb12(6) + alam
        Cheb24(20) = Cheb12(6) - alam
        alam = x(7)*v(2) - x(3)*v(4) - x(11)*v(6) + x(1)*v(8) - x(9)*v(10) &
               - x(5)*v(12)
        Cheb24(8) = Cheb12(8) + alam
        Cheb24(18) = Cheb12(8) - alam
        do i = 1, 6
            j = 14 - i
            v(i) = Fval(i) - Fval(j)
            Fval(i) = Fval(i) + Fval(j)
        end do
        alam1 = v(1) + x(8)*v(5)
        alam2 = x(4)*v(3)
        Cheb12(3) = alam1 + alam2
        Cheb12(11) = alam1 - alam2
        Cheb12(7) = v(1) - v(5)
        alam = x(2)*v(2) + x(6)*v(4) + x(10)*v(6)
        Cheb24(3) = Cheb12(3) + alam
        Cheb24(23) = Cheb12(3) - alam
        alam = x(6)*(v(2) - v(4) - v(6))
        Cheb24(7) = Cheb12(7) + alam
        Cheb24(19) = Cheb12(7) - alam
        alam = x(10)*v(2) - x(6)*v(4) + x(2)*v(6)
        Cheb24(11) = Cheb12(11) + alam
        Cheb24(15) = Cheb12(11) - alam
        do i = 1, 3
            j = 8 - i
            v(i) = Fval(i) - Fval(j)
            Fval(i) = Fval(i) + Fval(j)
        end do
        Cheb12(5) = v(1) + x(8)*v(3)
        Cheb12(9) = Fval(1) - x(8)*Fval(3)
        alam = x(4)*v(2)
        Cheb24(5) = Cheb12(5) + alam
        Cheb24(21) = Cheb12(5) - alam
        alam = x(8)*Fval(2) - Fval(4)
        Cheb24(9) = Cheb12(9) + alam
        Cheb24(17) = Cheb12(9) - alam
        Cheb12(1) = Fval(1) + Fval(3)
        alam = Fval(2) + Fval(4)
        Cheb24(1) = Cheb12(1) + alam
        Cheb24(25) = Cheb12(1) - alam
        Cheb12(13) = v(1) - v(3)
        Cheb24(13) = Cheb12(13)
        alam = 1.0_wp/6.0_wp
        do i = 2, 12
            Cheb12(i) = Cheb12(i)*alam
        end do
        alam = 0.5_wp*alam
        Cheb12(1) = Cheb12(1)*alam
        Cheb12(13) = Cheb12(13)*alam
        do i = 2, 24
            Cheb24(i) = Cheb24(i)*alam
        end do
        Cheb24(1) = 0.5_wp*alam*Cheb24(1)
        Cheb24(25) = 0.5_wp*alam*Cheb24(25)

    end subroutine dqcheb