stdlib_quadrature.fypp Source File

#:include "common.fypp"
module stdlib_quadrature
    !! ([Specification](../page/specs/stdlib_quadrature.html#description))
    use stdlib_kinds, only: sp, dp, xdp, qp

    implicit none

    private

    ! array integration
    public :: trapz
    public :: trapz_weights
    public :: simps
    public :: simps_weights
    public :: gauss_legendre
    public :: gauss_legendre_lobatto


    interface trapz
        !! version: experimental
        !!
        !! Integrates sampled values using trapezoidal rule
        !! ([Specification](../page/specs/stdlib_quadrature.html#description))
        #:for k1, t1 in REAL_KINDS_TYPES
        pure module function trapz_dx_${k1}$(y, dx) result(integral)
          ${t1}$, dimension(:), intent(in) :: y
          ${t1}$, intent(in) :: dx
          ${t1}$ :: integral
        end function trapz_dx_${k1}$
        #:endfor
        #:for k1, t1 in REAL_KINDS_TYPES
        module function trapz_x_${k1}$(y, x) result(integral)
            ${t1}$, dimension(:), intent(in) :: y
            ${t1}$, dimension(:), intent(in) :: x
            ${t1}$ :: integral
        end function trapz_x_${k1}$
        #:endfor
    end interface trapz


    interface trapz_weights
        !! version: experimental
        !!
        !! Integrates sampled values using trapezoidal rule weights for given abscissas
        !! ([Specification](../page/specs/stdlib_quadrature.html#description_1))
        #:for k1, t1 in REAL_KINDS_TYPES
        pure module function trapz_weights_${k1}$(x) result(w)
            ${t1}$, dimension(:), intent(in) :: x
            ${t1}$, dimension(size(x)) :: w
        end function trapz_weights_${k1}$
        #:endfor
    end interface trapz_weights


    interface simps
        !! version: experimental
        !!
        !! Integrates sampled values using Simpson's rule
        !! ([Specification](../page/specs/stdlib_quadrature.html#description_3))
        ! "recursive" is an implementation detail
        #:for k1, t1 in REAL_KINDS_TYPES
        pure recursive module function simps_dx_${k1}$(y, dx, even) result(integral)
            ${t1}$, dimension(:), intent(in) :: y
            ${t1}$, intent(in) :: dx
            integer, intent(in), optional :: even
            ${t1}$ :: integral
        end function simps_dx_${k1}$
        #:endfor
        #:for k1, t1 in REAL_KINDS_TYPES
        recursive module function simps_x_${k1}$(y, x, even) result(integral)
            ${t1}$, dimension(:), intent(in) :: y
            ${t1}$, dimension(:), intent(in) :: x
            integer, intent(in), optional :: even
            ${t1}$ :: integral
        end function simps_x_${k1}$
        #:endfor
    end interface simps


    interface simps_weights
        !! version: experimental
        !!
        !! Integrates sampled values using trapezoidal rule weights for given abscissas
        !! ([Specification](../page/specs/stdlib_quadrature.html#description_3))
        #:for k1, t1 in REAL_KINDS_TYPES
        pure recursive module function simps_weights_${k1}$(x, even) result(w)
            ${t1}$, dimension(:), intent(in) :: x
            integer, intent(in), optional :: even
            ${t1}$, dimension(size(x)) :: w
        end function simps_weights_${k1}$
        #:endfor
    end interface simps_weights


    interface gauss_legendre
        !! version: experimental
        !!
        !! Computes Gauss-Legendre quadrature nodes and weights.
        pure module subroutine gauss_legendre_fp64 (x, w, interval)
            real(dp), intent(out) :: x(:), w(:)
            real(dp), intent(in), optional :: interval(2)
        end subroutine
    end interface gauss_legendre


    interface gauss_legendre_lobatto
        !! version: experimental
        !!
        !! Computes Gauss-Legendre-Lobatto quadrature nodes and weights.
        pure module subroutine gauss_legendre_lobatto_fp64 (x, w, interval)
            real(dp), intent(out) :: x(:), w(:)
            real(dp), intent(in), optional :: interval(2)
        end subroutine
    end interface gauss_legendre_lobatto


    ! Interface for a simple f(x)-style integrand function.
    ! Could become fancier as we learn about the performance
    ! ramifications of different ways to do callbacks.
    abstract interface
        #:for k1, t1 in REAL_KINDS_TYPES
        pure function integrand_${k1}$(x) result(f)
            import :: ${k1}$
            ${t1}$, intent(in) :: x
            ${t1}$ :: f
        end function integrand_${k1}$
        #:endfor
    end interface



end module stdlib_quadrature