step.M: Mathur's AutoDX-like automatic step selection

In pnd: Parallel Numerical Derivatives, Gradients, Jacobians, and Hessians of Arbitrary Accuracy Order

Mathur's AutoDX-like automatic step selection

Description

Mathur's AutoDX-like automatic step selection

Usage

step.M(
  FUN,
  x,
  h0 = NULL,
  deriv.order = 1,
  acc.order = 2,
  range = NULL,
  shrink.factor = 0.5,
  min.valid.slopes = 5L,
  seq.tol = 0.1,
  correction = TRUE,
  max.rel.error = .Machine$double.eps^(7/8),
  cores = 1,
  preschedule = getOption("pnd.preschedule", TRUE),
  cl = NULL,
  ...
)

Arguments

FUN	Function for which the optimal numerical derivative step size is needed.
x	Numeric scalar: the point at which the derivative is computed and the optimal step size is estimated.
h0	Numeric scalar: initial step size, defaulting to a relative step of slightly greater than .Machine$double.eps^(1/3) (or absolute step if x == 0).
deriv.order	Integer or vector of integers indicating the desired derivative order, \mathrm{d}^m / \mathrm{d}x^m, for each element of x.
acc.order	Integer or vector of integers specifying the desired accuracy order for each element of x. The final error will be of the order O(h^{\mathrm{acc.order}}).
range	Numeric vector of length 2 defining the valid search range for the step size.
shrink.factor	A scalar less than 1 that is used to create a sequence of step sizes. The recommended value is 0.5. Change to 0.25 for a faster search. This number should be a negative power of 2 for the most accurate representation.
min.valid.slopes	Positive integer: how many points must form a sequence with the correct slope with relative difference from 2 less than seq.tol. If shrink.factor is small (< 0.33), consider reducing this to 4.
seq.tol	Numeric scalar: maximum relative difference between old and new step sizes for declaring convergence.
correction	Logical: if TRUE, returns the corrected step size (last point in the sequence times a less-than-1 number to account for the possible continuation of the downwards slope of the total error); otherwise, returns the grid point that is is lowest in the increasing sequence of valid error estimates.
max.rel.error	Error bound for the relative function-evaluation error (\frac{\hat f(\hat x) - f(x)}{f(x)}). Measures how noisy a function is. If the function is relying on numerical optimisation routines, consider setting to sqrt(.Machine$double.eps). If the function has full precision to the last bit, set to .Machine$double.eps/2.
cores	Integer specifying the number of CPU cores used for parallel computation. Recommended to be set to the number of physical cores on the machine minus one.
preschedule	Logical: if TRUE, disables pre-scheduling for mclapply() or enables load balancing with parLapplyLB(). Recommended for functions that take less than 0.1 s per evaluation.
cl	An optional user-supplied cluster object (created by makeCluster or similar functions). If not NULL, the code uses parLapply() (if preschedule is TRUE) or parLapplyLB() on that cluster on Windows, and mclapply (fork cluster) on everything else.
...	Passed to FUN.

Details

This function computes the optimal step size for central differences using the \insertCitemathur2012analyticalpnd algorithm. It consists of the following steps.

1. Choose a reasonable large (but not too large) initial step size h_0 and a reduction factor (1/2 for fast, 1/4 for slow functions is a reasonable choice).

2. Compute a series of truncation error estimates via third derivatives or Richardson extrapolation.

3. Find the leftmost range of consecituve step sizes for which the slope of the trunctation error it approximately equal (within 10% tolerance) to the accuracy order and which is long enough (e.g. at least length 5).

4. Use the leftmost point of this range as the uncorrected optimal step size, or correct it by shrinking it by a small amount given in the article.

Value

A list similar to the one returned by optim():

• par – the optimal step size found.

• value – the estimated numerical first derivative (using central differences).

• counts – the number of iterations (each iteration includes two function evaluations).

• abs.error – an estimate of the truncation and rounding errors.

• exitcode – an integer code indicating the termination status:

  • 0 – Optimal termination due to a sequence of correct reductions.

  • 1 – Reductions are slightly outside the tolerance.

  • 2 – Tolerances are significantly violated; an approximate minimum is returned.

  • 3 – Not enough finite function values; a rule-of-thumb value is returned.

  • 4 – Step trimmed to 0.1|x| when |x| is not tiny and within range.

• message – A summary message of the exit status.

• iterations – A list including the step and argument grids, function values on those grids, estimated derivative values, and estimated error values.

References

\insertAllCited

Examples

f <- function(x) x^4  # The derivative at 1 is 4
step.M(x = 1, f)
step.M(x = 1, f, h0 = 1e-9) # Starting low
step.M(x = 1, f, h0 = 1000) # Starting high

f <- sin  # The derivative at pi/4 is sqrt(2)/2
plot(step.M(x = pi/2, f))
plot(step.M(x = pi/4, f))
step.M(x = pi/4, f, h0 = 1e-9) # Starting low
step.M(x = pi/4, f, h0 = 1000) # Starting high
# where the truncation error estimate is invalid

pnd documentation built on Sept. 9, 2025, 5:44 p.m.