UltrasphericalGaussRadauQuadrature

This routine returns the $n+1$ Quadrature points and weights.

The Gauss-Radau quadrature points consists one of the end points denoted by $a$ .

So $a$ can be $\pm 1$ .

The remaining $n$ points are internal to $(-1, +1)$ , and they are n-zeros of Ultraspherical polynomial of order n.

Here n is the order of Ultraspherical polynomial.

If $a=1$ then n+1 quadrature point will be +1.
If $a=-1$ then 1st quadrature point will be -1.

Interface

܀ Interface
️܀ See example
↢

INTERFACE
  MODULE SUBROUTINE UltrasphericalGaussRadauQuadrature(a, n, lambda, pt, wt)
    REAL(DFP), INTENT(IN) :: a
    !! the value of one of the end points
    !! it should be either -1 or +1
    INTEGER(I4B), INTENT(IN) :: n
    !! order of Ultraspherical polynomial
    REAL(DFP), INTENT(IN) :: lambda
    !! lambda should be greater than -0.5
    REAL(DFP), INTENT(OUT) :: pt(:)
    !! n+1 quadrature points from 1 to n+1
    REAL(DFP), OPTIONAL, INTENT(OUT) :: wt(:)
    !! n+1 weights from 1 to n+1
  END SUBROUTINE UltrasphericalGaussRadauQuadrature
END INTERFACE

This example shows the usage of UltrasphericalGaussRadauQuadrature method.

Note that this routine returns n+1 quadrature points and one of the points will be end point.

program main
  use easifembase
  implicit none
  integer( i4b ) :: n
  real( dfp ), allocatable :: pt( : ), wt( : )
  real( dfp ), parameter :: left=-1.0, right=1.0, lambda=0.5_DFP
  type(string) :: msg, astr
  real( dfp ) :: a

  n = 1; a=left; call callme

See results

Ultraspherical Gauss Radau points, n+1 = 2


-1	0.5
0.33333	1.5

  n = 2; a=left; call callme

See results

Ultraspherical Gauss Radau points, n+1 = 3


-1	0.22222
-0.2899	1.025
0.6899	0.75281

  n = 3; a=left; call callme

See results

Ultraspherical Gauss Radau points, n+1 = 4


-1	0.125
-0.57532	0.65769
0.18107	0.77639
0.82282	0.44092

  n = 4; a=left; call callme

See results

Ultraspherical Gauss Radau points, n+1 = 5


-1	8E-02
-0.72048	0.44621
-0.16718	0.62365
0.44631	0.56271
0.88579	0.28743

  n = 4; a=right; call callme

See results

Ultraspherical Gauss Radau points, n+1 = 5


-0.88579	0.28743
-0.44631	0.56271
0.16718	0.62365
0.72048	0.44621
1	8E-02

  contains
  subroutine callme
    call reallocate( pt, n+1, wt, n+1 )
    call UltrasphericalGaussRadauQuadrature( a=a, n=n, &
      & pt=pt, wt=wt, lambda=lambda )
    msg="Ultraspherical Gauss Radau points, n+1 = " &
        & // tostring( n+1 )
    call display(msg%chars())
    astr = MdEncode( pt .COLCONCAT. wt )
    call display( astr%chars(), "" )
  end subroutine callme

end program main

UltrasphericalGaussRadauQuadrature

Interface​

Interface