JacobiGaussQuadrature

This routine computes the n Gauss-Quadrature points.

They are n zeros of a jacobi polynomial defined with respect to the weight $(1-x)^{\alpha} (1+x)^{\beta}$.

All Gauss-Quadrature points are inside $(-1, 1)$

Interface

INTERFACE
  MODULE SUBROUTINE JacobiGaussQuadrature(n, alpha, beta, pt, wt)
    INTEGER(I4B), INTENT(IN) :: n
    !! It represents the order of Jacobi polynomial
    REAL(DFP), INTENT(IN) :: alpha
    REAL(DFP), INTENT(IN) :: beta
    REAL(DFP), INTENT(OUT) :: pt(:)
    !! the size is 1 to n
    REAL(DFP), OPTIONAL, INTENT(OUT) :: wt(:)
    !! the size is 1 to n
  END SUBROUTINE JacobiGaussQuadrature
END INTERFACE

Examples

️܀ See example

This example shows the usage of JacobiGaussQuadrature method. This routine returns the quadrature points for Jacobi weights.

In this example $\alpha=\beta=0.0$ (that is, Legendre polynomial)

program main
  use easifembase
  implicit none
  integer( i4b ) :: n
  real( dfp ), allocatable :: pt( : ), wt( : )
  real( dfp ), parameter :: alpha=0.0_DFP, beta=0.0_DFP
  type(string) :: msg, astr

order = 1

  n = 1; call callme

See results

Zeros of J(x), n = 1 alpha=0 beta=0

pt	wt
0	2

order = 2

  n = 2; call callme

See results

Zeros of J(x), n = 1 alpha=0 beta=0

pt	wt
-0.57735	1
0.57735	1

order = 5

  n = 5; call callme

See results

Zeros of J(x), n = 5 alpha=0 beta=0

pt	wt
-0.90618	0.23693
-0.53847	0.47863
-1.56541E-16	0.56889
0.53847	0.47863
0.90618	0.23693

  contains
  subroutine callme
    call reallocate( pt, n, wt, n )
    call JacobiGaussQuadrature( n=n, alpha=alpha, beta=beta, pt=pt, wt=wt )
    msg = "Zeros of J(x), n = " &
        & // tostring( n ) // " alpha="//tostring( alpha ) // &
        & " beta="//tostring( beta )
    call display(msg%chars())
    astr = MdEncode( pt .COLCONCAT. wt )
    call display( astr%chars(), "" )
  end subroutine callme

end program main

↢