poly_legendre Subroutine

  1. pf_quadrature.f90
  2. pf_mod_quadrature
  3. poly_legendre

public subroutine poly_legendre(p, n)

Subroutine to compute Legendre polynomial coefficients using Bonnet's recursion formula.

Arguments

Type IntentOptional AttributesName
real(kind=pfqp), intent(out) :: p(0:n)
integer, intent(in), value:: n

Called by

proc~~poly_legendre~~CalledByGraph proc~poly_legendre poly_legendre proc~sdc_qnodes sdc_qnodes proc~sdc_qnodes->proc~poly_legendre proc~pf_quadrature pf_quadrature proc~pf_quadrature->proc~sdc_qnodes proc~pf_init_sdcmats pf_init_sdcmats proc~pf_init_sdcmats->proc~pf_quadrature proc~pf_level_setup pf_level_setup proc~pf_level_setup->proc~pf_init_sdcmats proc~pf_pfasst_setup pf_pfasst_setup proc~pf_pfasst_setup->proc~pf_level_setup

Contents

Source Code

poly_legendre

Source Code

  subroutine poly_legendre(p, n) 
    integer, intent(in), value :: n
    real(pfqp),       intent(out)       :: p(0:n)

    real(pfqp), dimension(0:n) :: p0, p1, p2
    integer :: j, m

    if (n == 0) then
       p = [ 1.0_pfqp ]
       return
    end if

    if (n == 1) then
       p = [ 0.0_pfqp, 1.0_pfqp ]
       return
    end if

    p0 = 0.0_pfqp;
    p1 = 0.0_pfqp;
    p2 = 0.0_pfqp

    p0(0) = 1.0_pfqp
    p1(1) = 1.0_pfqp

    ! (n + 1) P_{n+1} = (2n + 1) x P_{n} - n P_{n-1}
    do m = 1, n-1
       do j = 1, n
          p2(j) = ( (2*m + 1) * p1(j-1) - m * p0(j) ) / (m + 1)
       end do
       p2(0) = - m * p0(0) / (m + 1)

       p0 = p1
       p1 = p2
    end do

    p = p2
  end subroutine poly_legendre