LegendreEval
Evaluate Legendre polynomials of order n at single or several points.
Interface 1
- ܀ Interface
- ️܀ See example
- ↢
INTERFACE
MODULE PURE FUNCTION LegendreEval(n, x) RESULT(ans)
INTEGER(I4B), INTENT(IN) :: n
!! order of polynomial
REAL(DFP), INTENT(IN) :: x
!! point of evaluation, it should be between -1 and 1
REAL(DFP) :: ans
!! Evaluate Legendre polynomial of order n at point x
END FUNCTION LegendreEval
END INTERFACE
This example shows the usage of LegendreEval
method.
This routine evaluates Legendre polynomial of order n, at single point
program main
use easifembase
implicit none
integer( i4b ) :: n
real( dfp ) :: ans, x, exact
real( dfp ), parameter :: 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= LegendreEval( n=n, x=x )
end subroutine callme
end program main
Interface 2
- ܀ Interface
- ️܀ See example
- ↢
INTERFACE
MODULE PURE FUNCTION LegendreEval(n, x) RESULT(ans)
INTEGER(I4B), INTENT(IN) :: n
!! order of polynomial
REAL(DFP), INTENT(IN) :: x(:)
!! several points of evaluation
REAL(DFP) :: ans(SIZE(x))
!! Evaluate Legendre polynomial of order n at points x
END FUNCTION LegendreEval
END INTERFACE
This example shows the usage of LegendreEval
method.
This routine evaluates Legendre polynomial of order n, at several points.
program main
use easifembase
implicit none
integer( i4b ) :: n
real( dfp ), allocatable :: ans(:), x(:), exact(:)
real( dfp ), parameter :: tol=1.0E-10
type(string) :: astr
n = 3
x = [-1.0, 0.0, 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= LegendreEval( n=n, x=x )
end subroutine callme
end program main