diff_module Module

Numerical differentiation of a 1D function f(x) using Neville's process.

Authors

  • J. Oliver, "An algorithm for numerical differentiation of a function of one real variable", Journal of Computational and Applied Mathematics 6 (2) (1980) 145–160. [Algol 60 source in original paper]
  • David Kahaner, Fortran 77 code from NIST
  • Jacob Williams : 2/17/2013 : Converted to modern Fortran. Some refactoring, addition of test cases.

Uses

  • module~~diff_module~~UsesGraph module~diff_module diff_module module~numdiff_kinds_module numdiff_kinds_module module~diff_module->module~numdiff_kinds_module iso_fortran_env iso_fortran_env module~numdiff_kinds_module->iso_fortran_env

Used by

  • module~~diff_module~~UsedByGraph module~diff_module diff_module module~numerical_differentiation_module numerical_differentiation_module module~numerical_differentiation_module->module~diff_module

Abstract Interfaces

abstract interface

  • private function func(me, x) result(fx)

    interface to function for diff

    Arguments

    Type IntentOptional Attributes Name
    class(diff_func), intent(inout) :: me
    real(kind=wp), intent(in) :: x

    Return Value real(kind=wp)


Derived Types

type, public ::  diff_func

class to define the function for diff

Components

Type Visibility Attributes Name Initial
logical, private :: stop = .false.
procedure(func), private, pointer :: f => null()

Type-Bound Procedures

procedure, private :: faccur
procedure, public :: set_function
procedure, public :: compute_derivative => diff
procedure, public :: terminate

Subroutines

private subroutine set_function(me, f)

Author
Jacob Williams
Date
12/27/2015

Set the function in a diff_func. Must be called before diff.

Arguments

Type IntentOptional Attributes Name
class(diff_func), intent(inout) :: me
procedure(func) :: f

private subroutine terminate(me)

Can be called by the user in the function to terminate the computation. This will set ifail=-1.

Arguments

Type IntentOptional Attributes Name
class(diff_func), intent(inout) :: me

private subroutine diff(me, iord, x0, xmin, xmax, eps, accr, deriv, error, ifail)

the procedure diff calculates the first, second or third order derivative of a function by using neville's process to extrapolate from a sequence of simple polynomial approximations based on interpolating points distributed symmetrically about x0 (or lying only on one side of x0 should this be necessary). if the specified tolerance is non-zero then the procedure attempts to satisfy this absolute or relative accuracy requirement, while if it is unsuccessful or if the tolerance is set to zero then the result having the minimum achievable estimated error is returned instead.

Arguments

Type IntentOptional Attributes Name
class(diff_func), intent(inout) :: me
integer, intent(in) :: iord

1, 2 or 3 specifies that the first, second or third order derivative,respectively, is required.
real(kind=wp), intent(in) :: x0

the point at which the derivative of the function is to be calculated.
real(kind=wp), intent(in) :: xmin

xmin, xmax restrict the interpolating points to lie in [xmin, xmax], which should be the largest interval including x0 in which the function is calculable and continuous.
real(kind=wp), intent(in) :: xmax

xmin, xmax restrict the interpolating points to lie in [xmin, xmax], which should be the largest interval including x0 in which the function is calculable and continuous.
real(kind=wp), intent(in) :: eps

denotes the tolerance, either absolute or relative. eps=0 specifies that the error is to be minimised, while eps>0 or eps<0 specifies that the absolute or relative error, respectively, must not exceed abs(eps) if possible. the accuracy requirement should not be made stricter than necessary, since the amount of computation tends to increase as the magnitude of eps decreases, and is particularly high when eps=0.
real(kind=wp), intent(in) :: accr

denotes that the absolute (accr>0) or relative (accr<0) errors in the computed values of the function are most unlikely to exceed abs(accr), which should be as small as possible. if the user cannot estimate accr with complete confidence, then it should be set to zero.
real(kind=wp), intent(out) :: deriv

the calculated value of the derivative
real(kind=wp), intent(out) :: error

an estimated upper bound on the magnitude of the absolute error in the calculated result. it should always be examined, since in extreme case may indicate that there are no correct significant digits in the value returned for derivative.
integer, intent(out) :: ifail

will have one of the following values on exit: 0 the procedure was successful. 1 the estimated error in the result exceeds the (non-zero) requested error, but the most accurate result possible has been returned. 2 input data incorrect (derivative and error will be undefined). 3 the interval [xmin, xmax] is too small (derivative and error will be undefined). -1 stopped by the user.

private subroutine faccur(me, h0, h1, facc, x0, twoinf, f0, f1, ifail)

This procedure attempts to estimate the level of rounding errors in the calculated function values near the point x0+h0 by fitting a least-squares straight-line approximation to the function at the six points x0+h0-j*h1, (j = 0,1,3,5,7,9), and then setting facc to twice the largest deviation of the function values from this line. hi is adjusted if necessary so that it is approximately 8 times the smallest spacing at which the function values are unequal near x0+h0.

Arguments

Type IntentOptional Attributes Name
class(diff_func), intent(inout) :: me
real(kind=wp), intent(inout) :: h0
real(kind=wp), intent(inout) :: h1
real(kind=wp), intent(out) :: facc
real(kind=wp), intent(in) :: x0
real(kind=wp), intent(in) :: twoinf
real(kind=wp), intent(in) :: f0
real(kind=wp), intent(in) :: f1
integer, intent(out) :: ifail

0 if no error, -1 if user termination.