fastmath.calculus — generateme/fastmath 2.4.0

Integration and derivatives

Integrate univariate and multivariate functions.

• VEGAS / VEGAS+ - Monte Carlo integration of multivariate function
• h-Cubature - h-adaptive integration of multivariate function
• Guass-Kronrod and Gauss-Legendre - quadrature integration of univariate functions
• Romberg, Simpson, MidPoint and Trapezoid

Integrant is substituted in case of improper integration bounds.

Derivatives (finite differences method):

• derivatives of any degree and any order of accuracy
• gradient and hessian for multivariate functions

cubature clj

(cubature f lower upper)

(cubature f
          lower
          upper
          {:keys [rel abs max-evals max-iters initdiv info?]
           :or {rel 1.0E-7
                abs 1.0E-7
                max-evals Integer/MAX_VALUE
                max-iters 64
                initdiv 2
                info? false}})

Cubature - h-adaptive integration of multivariate function, n>1 dimensions.

Algorithm uses Genz Malik method.

In each iteration a box with biggest error is subdivided and reevaluated.

Improper integrals with infinite bounds are handled by a substitution.

Arguments:

• f - integrant
• lower - seq of lower bounds
• upper - seq of upper bounds

Options:

• :initvid - initial subdivision per dimension, default: 2.
• :max-evals - maximum number of evaluations, default: max integer value.
• :max-iters - maximum number of iterations, default: 64.
• :rel - relative error, 1.0e-7
• :abs - absolute error, 1.0e-7
• :info? - return full information about integration, default: false

Function returns a map containing (if info? is true, returns result otherwise):

• :result - integration value
• :error - integration error
• :iterations - number of iterations
• :evaluations - number of evaluations
• :subdivisions - final number of boxes
• :fail? - set to :max-evals or :max-iters when one of the limits has been reached without the convergence.

derivative clj

(derivative f)

(derivative f n)

(derivative f
            n
            {:keys [acc h method extrapolate?]
             :or {acc 2 h 0.0 method :central}})

Create nth derivative of f using finite difference method for given accuracy :acc and step :h.

Returns function.

Arguments:

• n - derivative
• :acc - order of accuracy (default: 2)
• :h - step, (default: 0.0, automatic)
• :method - :central (default), :forward or :backward
• :extrapolate? - creates extrapolated derivative if set to true or a map with extrapolate function options

extrapolate clj

(extrapolate g)

(extrapolate g
             {:keys [contract power init-h rel abs max-evals tol]
              :or {contract 0.5
                   power 1.0
                   init-h 0.5
                   abs m/MACHINE-EPSILON
                   rel (m/sqrt (m/ulp init-h))
                   max-evals Integer/MAX_VALUE
                   tol 2.0}})

Richardson extrapolation for given function g=g(x,h). Returns extrapolated function f(x).

Options:

• :contract - shrinkage factor, default=1/2
• :power - set to 2.0 for even functions around x0, default 1.0
• :init-h - initial step h, default=1/2
• :abs - absolute error, default: machine epsilon
• :rel - relative error, default: ulp for init-h
• :tol - tolerance for error, default: 2.0
• :max-evals - maximum evaluations, default: maximum integer

f' clj

(f' f)

First central derivative with order of accuracy 2.

f'' clj

(f'' f)

Second central derivative with order of accuracy 2.

f''' clj

(f''' f)

Third central derivative with order of accuracy 2.

fx->gx+h clj

(fx->gx+h f)

Convert f(x) to g(x,h)=f(x+h)

fx->gx-h clj

(fx->gx-h f)

Convert f(x) to g(x,h)=f(x-h)

gradient clj

(gradient f)

(gradient f {:keys [h acc] :or {h 1.0E-6 acc 2}})

Create first partial derivatives of multivariate function for given accuracy :acc and step :h.

Returns function.

Options:

• :acc - order of accuracy, 2 (default) or 4.
• :h - step, default 1.0e-6

hessian clj

(hessian f)

(hessian f {:keys [h] :or {h 0.005}})

Creates function returning Hessian matrix for mulitvariate function f and given :h step (default: 5.0e-3).

integrate clj

(integrate f)

(integrate f lower upper)

(integrate f
           lower
           upper
           {:keys [rel abs max-iters min-iters max-evals info? integrator
                   integration-points]
            :or {rel BaseAbstractUnivariateIntegrator/DEFAULT_RELATIVE_ACCURACY
                 abs BaseAbstractUnivariateIntegrator/DEFAULT_ABSOLUTE_ACCURACY
                 min-iters
                   BaseAbstractUnivariateIntegrator/DEFAULT_MIN_ITERATIONS_COUNT
                 max-evals Integer/MAX_VALUE
                 integration-points 7
                 integrator :gauss-kronrod
                 info? false}
            :as options})

Univariate integration.

Improper integrals with infinite bounds are handled by a substitution.

Arguments:

• f - integrant
• lower - lower bound
• upper - upper bound

Options:

• :integrator - integration algorithm, one of: :romberg, :trapezoid, :midpoint, :simpson, :gauss-legendre and :gauss-kronrod (default).
• :min-iters - minimum number of iterations (default: 3), not used in :gauss-kronrod
• :max-iters - maximum number of iterations (default: 32 or 64)
• :max-evals - maximum number of evaluations, (default: maximum integer)
• :rel - relative error
• :abs - absolute error
• :integration-points - number of integration (quadrature) points for :gauss-legendre and :gauss-kronrod, default 7
• :initdiv - initial number of subdivisions for :gauss-kronrod, default: 1
• :info? - return full information about integration, default: false

:gauss-kronrod is h-adaptive implementation

Function returns a map containing (if info? is true, returns result otherwise):

• :result - integration value
• :error - integration error (:gauss-kronrod only)
• :iterations - number of iterations
• :evaluations - number of evaluations
• :subdivisions - final number of boxes (:gauss-kronrod only)
• :fail? - set to :max-evals or :max-iters when one of the limits has been reached without the convergence.

vegas clj

(vegas f lower upper)

(vegas f
       lower
       upper
       {:keys [max-iters rel abs nevals alpha beta warmup info? record-data?]
        :or {max-iters 10
             rel 5.0E-4
             abs 5.0E-4
             nevals 10000
             beta 0.75
             warmup 0
             info? false
             record-data? false}
        :as options})

VEGAS+ - Monte Carlo integration of multivariate function, n>1 dimensions.

Improper integrals with infinite bounds are handled by a substitution.

Arguments:

• f - integrant
• lower - seq of lower bounds
• upper - seq of upper bounds

Additional options:

• :max-iters - maximum number of iterations, default: 10
• :nevals - number of evaluations per iteration, default: 10000
• :nintervals - number of grid intervals per dimension (default: 1000)
• :nstrats - number of stratifications per dimension (calculated)
• :warmup - number of warmup iterations (results are used to train stratification and grid spacings, default: 0
• :alpha - grid refinement parameter, 0.5 slow (default for vegas+), 1.5 moderate/fast (defatult for vegas)
• :beta - stratification damping parameter for startification adaptation, default: 0.75
• :rel - relative accuracy, default: 5.0e-4
• :abs - absolute accuracy, default: 5.0e-4
• :random-sequence - random sequence used for generating samples: :uniform (default), low-discrepancy sequences: :r2, :sobol and :halton.
• :jitter - jittering factor for low-discrepancy random sequence, default: 0.75
• :info? - return full information about integration, default: false
• :record-data? - stores samples, number of strata, x and dx, default: false (requires, :info? to be set to true)

For original VEGAS algorithm set :nstrats to 1.

:nstrats can be also a list, then each dimension is divided independently according to a given number. If list is lower then number of dimensions, then it's cycled.

Function returns a map with following keys (if info? is true, returns result otherwise):

• :result - value of integral
• :iterations - number of iterations (excluding warmup)
• :sd - standard deviation of results
• :nintervals - actual grid size
• :nstrats - number of stratitfications per dimension
• :nhcubes - number of hypercubes
• :evaluations - number of function calls
• :chi2-avg - average of chi2
• :dof - degrees of freedom
• :Q - goodness of fit indicator, 1 - very good, <0.25 very poor
• :data - recorded data (if available)