JacobiEval
Evaluate Jacobi polynomials of order n at single or several points.
Interface 1
INTERFACE
MODULE PURE FUNCTION JacobiEval(n, alpha, beta, x) RESULT(ans)
INTEGER(I4B), INTENT(IN) :: n
REAL(DFP), INTENT(IN) :: alpha
REAL(DFP), INTENT(IN) :: beta
REAL(DFP), INTENT(IN) :: x
REAL(DFP) :: ans
!! Evaluate Jacobi polynomial of order n at point x
END FUNCTION JacobiEval
END INTERFACE
Evaluate Jacobi polynomial of order n at single points.
- N, the order of polynomial to compute.
- alpha, beta are parameters
- x: the point at which the polynomials are to be evaluated.
Interface 2
INTERFACE
MODULE PURE FUNCTION JacobiEval(n, alpha, beta, x) RESULT(ans)
INTEGER(I4B), INTENT(IN) :: n
REAL(DFP), INTENT(IN) :: alpha
REAL(DFP), INTENT(IN) :: beta
REAL(DFP), INTENT(IN) :: x(:)
REAL(DFP) :: ans(SIZE(x))
!! Evaluate Jacobi polynomial of order n at point x
END FUNCTION JacobiEval
END INTERFACE
- N is order of polynomial to compute.
- alpha, beta are parameters
- x: the point at which the polynomials are to be evaluated.
- ans, the values of the Jacobi polynomials at the several points.
- ️܀ See example
- ↢
This example shows the usage of JacobiEval
method.
This routine evaluates Jacobi polynomial of order n, at single point
In this example (that is, Legendre polynomial)
program main
use easifembase
implicit none
integer( i4b ) :: n
real( dfp ) :: ans, x, exact
real( dfp ), parameter :: alpha=0.0_DFP, beta=0.0_DFP, tol=1.0E-10
n = 3
x = -1.0_DFP; call callme
exact = 0.5_DFP*(5.0 * x**3 - 3.0*x)
call ok( SOFTEQ(ans, exact, tol ))
x = 0.0_DFP; call callme
exact = 0.5_DFP*(5.0 * x**3 - 3.0*x)
call ok( SOFTEQ(ans, exact, tol ))
x = +1.0_DFP; call callme
exact = 0.5_DFP*(5.0 * x**3 - 3.0*x)
call ok( SOFTEQ(ans, exact, tol ))
contains
subroutine callme
ans= JacobiEval( n=n, alpha=alpha, beta=beta, x=x )
end subroutine callme
end program main
- ️܀ See example
- ↢
This example shows the usage of JacobiEval
method.
This routine evaluates Jacobi polynomial of order n, at several points.
In this example (that is, Legendre polynomial)
program main
use easifembase
implicit none
integer( i4b ) :: n
real( dfp ), allocatable :: ans(:), x(:), exact(:)
real( dfp ), parameter :: alpha=0.0_DFP, beta=0.0_DFP, tol=1.0E-10
type(string) :: astr
n = 3
x = [-1.0, 0.0, 0.25, 1.0]; call callme
exact = 0.5_DFP*(5.0 * x**3 - 3.0*x)
call ok( ALL(SOFTEQ(ans, exact, tol )))
contains
subroutine callme
ans= JacobiEval( n=n, alpha=alpha, beta=beta, x=x )
end subroutine callme
end program main