To estimate the derivative of sin(x) at x=1, a centred finite difference with step h=0.00001 gives an error of roughly 10⁻¹¹ in an ordinary binary64 test. Shrinking h all the way to 10⁻¹² makes the error roughly a million times larger. The calculus limit says smaller h is better in exact arithmetic—but a floating-point computer must subtract rounded function values before dividing by h.
Numerical differentiation exposes two competing error sources: approximation error from finite h and rounding error magnified as h becomes tiny. This article derives forward and central formulas, computes reproducible examples, explains the dependence on function smoothness and input precision, and shows where automatic differentiation or complex-step methods are better suited.
Find the question you need
Use the reading routes to reach a formula, worked comparison or method check. All teaching chapters remain open; the grouped index below is expandable.
Derive the finite differences
Chapters 1–3. Start here
Watch accuracy reverse
Chapters 4–6. Start here
Balance two error sources
Chapters 7–9. Start here
Choose an alternative method
Chapters 10–10. Start here
Practise and report
Chapters 11–12. Start here
Full chapter index · Worked practice · Scientific Computing subject root
Chapter index
Derive the finite differences · Chapters 1–3
Watch accuracy reverse · Chapters 4–6
Balance two error sources · Chapters 7–9
Choose an alternative method · Chapters 10–10
Practise and report · Chapters 11–12
The mathematical derivative at x is f′(x)=lim as h→0 of [f(x+h)−f(x)]/h, provided the limit exists. This is an exact statement about a function and infinitely accurate real-number operations. It does not instruct a finite-precision computer to divide two almost equal rounded values with an arbitrarily tiny h.
A simple numerical approximation chooses a nonzero step h, evaluates the function twice and forms the quotient. This is called a forward difference. Its value depends on h, so a result should never be presented as though the step size were irrelevant metadata.
We will test f(x)=sin x at x=1 radian, where the exact derivative is f′(1)=cos1≈0.540302305868. Choosing a function with a known analytic derivative lets us compare computed results with a trustworthy reference.
The point x=1 is expressed in radians. Using degrees inside the sine function while comparing against the radian derivative would invalidate the experiment before floating-point precision became relevant.
Contents · Next · Subject root
Expand f(x+h) by Taylor’s theorem for a sufficiently smooth function: f(x+h)=f(x)+hf′(x)+(h²/2)f″(x)+higher-order terms. Subtract f(x) and divide by h. The resulting forward difference is D₊(h)=[f(x+h)−f(x)]/h=f′(x)+(h/2)f″(x)+O(h²).
The leading discretisation error is proportional to h when the relevant derivatives are bounded. A step ten times smaller often reduces the leading truncation error by about ten for a smooth function in an arithmetic regime where roundoff does not dominate.
For f(x)=sin x, f″(x)=−sin x. At x=1, the leading forward-difference bias is approximately −h sin1/2. With h=0.1, that predicts about −0.0421, reasonably close to the exact observed error of roughly −0.04294.
The O(h) statement is an asymptotic mathematical relation under smoothness assumptions. It does not say that an actual binary64 result decreases its error monotonically as h tends toward the smallest representable number.
Contents · Previous · Next · Subject root
A centred approximation uses evaluations on both sides: D₀(h)=[f(x+h)−f(x−h)]/(2h). Taylor-expanding both functions and subtracting cancels the even-power second-derivative contributions in the numerator, leaving f′(x)+(h²/6)f‴(x)+O(h⁴).
The central formula therefore has second-order truncation error, O(h²), rather than the forward formula’s O(h). For the same moderately small h and a smooth function, this can yield substantially better accuracy at the cost of one evaluation on each side.
For sin(x), f‴(x)=−cos x, so the leading centred bias at x=1 is −h² cos1/6. At h=0.1 this is about −0.0009005, matching the observed central-difference error of about −0.0009001.
The method still computes a difference of two similar floating-point function values. Its higher order improves truncation error but does not remove cancellation or the magnification caused by division by a very small denominator.
Contents · Previous · Next · Subject root
Use h=0.1 and x=1. The forward quotient is [sin(1.1)−sin(1)]/0.1≈0.4973637525. Since cos1≈0.5403023059, the forward absolute error is about 0.04294.
The centred quotient is [sin(1.1)−sin(0.9)]/0.2≈0.5394022522, with absolute error roughly 0.0009001. This much smaller error comes from the cancellation of the leading Taylor-series error terms, not from an extra decimal setting on the calculator.
If a student rounds both sine values to three decimals before subtracting, the resulting quotient can differ substantially from the binary64 calculation. Keep intermediate precision appropriate to the method and distinguish displayed rounding from arithmetic precision.
This first test intentionally uses a step large enough that truncation error is visible. It allows the order-of-accuracy argument to be checked before the more subtle roundoff effect enters.
Contents · Previous · Next · Subject root
At h=0.01, the forward estimate is approximately 0.5360859810, with error 0.0042163. The centred result is approximately 0.5402933009, with error about 9.0×10⁻⁶. Both improve when h falls from 0.1 to 0.01.
At h=0.0001, the forward result is approximately 0.5402602314, with error 4.21×10⁻⁵. The centred result is approximately 0.54030230497, whose error is around 9.0×10⁻¹⁰. This is a particularly clear demonstration of the extra order of accuracy.
At h=0.00001, the centred estimate becomes roughly 0.540302305857 and is within about 1.1×10⁻¹¹ of cos1 for this evaluation. The approximation is very accurate, but that does not mean going six more orders of magnitude smaller must help.
These numbers are illustrative results for ordinary binary64 sine evaluations. Library implementations, extended intermediate precision and compiler optimisation can change last digits, so the claimed lesson should be the error mechanism and scale rather than one immutable machine bit pattern.
Contents · Previous · Next · Subject root
At h=10⁻¹⁰, the centred result in one binary64 evaluation is about 0.5403022474, roughly 5.8×10⁻⁸ away from cos1. That is much worse than the error around 10⁻¹¹ obtained at h=10⁻⁵.
At h=10⁻¹², the centred difference can yield about 0.5402900349, with absolute error around 1.23×10⁻⁵. The function values sin(1+h) and sin(1−h) are extremely close and subject to rounding, while the denominator 2h is tiny.
For h=10⁻¹⁶, the forward calculation in binary64 can return zero because 1+h rounds to exactly 1 before sine is evaluated. The numerator becomes zero, not the tiny nonzero difference that mathematical calculus intended. Another formulation or arithmetic precision is required.
The central formula may not fail in exactly the same way at the same h, because 1−h can round to a distinct representable value. But it is still unreliable there. The exact numerical symptom depends on representable spacing around the point x, not merely on the decimal digits written for h.
| Step h | Forward difference | Forward absolute error | Central difference | Central absolute error |
|---|---|---|---|---|
| 0.1 | 0.4973637525 | 4.294×10⁻² | 0.5394022522 | 9.001×10⁻⁴ |
| 0.01 | 0.5360859810 | 4.216×10⁻³ | 0.5402933009 | 9.005×10⁻⁶ |
| 10⁻⁴ | 0.5402602314 | 4.207×10⁻⁵ | 0.5403023050 | 9.004×10⁻¹⁰ |
| 10⁻⁵ | 0.5402980985 | 4.207×10⁻⁶ | 0.5403023059 | 1.114×10⁻¹¹ |
| 10⁻¹⁰ | 0.5403022474 | 5.848×10⁻⁸ | 0.5403022474 | 5.848×10⁻⁸ |
| 10⁻¹² | 0.5403455461 | 4.324×10⁻⁵ | 0.5402900349 | 1.227×10⁻⁵ |
Contents · Previous · Next · Subject root
A schematic forward-difference error model is C₁h+C₂ε/h, where ε represents a floating-point rounding scale and the constants depend on function magnitude and derivatives. The first term decreases as h shrinks; the second increases because cancellation errors are divided by h.
For a centred difference, the corresponding conceptual model is C₃h²+C₄ε/h. Balancing its two terms gives a characteristic h of order ε^(1/3) when the function scales and derivative magnitudes are of order one. For forward difference the analogous characteristic scale is ε^(1/2).
With binary64 ε around 10⁻¹⁶, these simple scales suggest forward steps around 10⁻⁸ and central steps around 10⁻⁵, consistent with our sin(1) experiment. They are not universally optimal constants and should not be copied to every input without considering units, cancellation and sensitivity.
The most important shape is U-like error as a function of decreasing step: approximation error first falls, then arithmetic noise increases. This is why taking h as small as possible is not a principled numerical differentiation strategy.
Contents · Previous · Next · Subject root
Suppose a model uses x measured in kilometres and another uses metres. A numerical step h=10⁻⁵ represents different physical changes in the two coordinate systems. A step-selection rule should account for the scale of x, expected derivative behaviour and data resolution, not a globally fixed decimal string.
If f(x) is very large while its derivative is small, subtracting f(x+h)−f(x) can be particularly vulnerable to losing meaningful digits. Conversely, a function with steep local curvature may need a smaller step to control truncation error. Both affect the balance of the error model.
A robust implementation often considers a relative step based on |x| and a minimum absolute scale, then compares estimates over several nearby h values. That is a diagnostic, not proof of convergence in every case.
When comparing two derivates from different studies, ask about units, coordinate transformations, scale and derivative definitions. A number without those details may look precise while answering a different physical question.
Contents · Previous · Next · Subject root
So far the sine function was evaluated to floating-point accuracy. Real data can carry measurement noise much larger than 10⁻¹⁶. Suppose two measurements of f have possible errors of ±0.001 units, and the forward difference uses h=0.001 in the corresponding x unit.
The difference of two errors can be as large as 0.002 units, so division by h=0.001 can inject an error of magnitude up to about two units into the derivative estimate, even before considering discretisation. The true derivative could be much smaller than that.
Making h ten times smaller would magnify this worst-case noise contribution by another factor of ten. A small derivative interval that is excellent for an analytic computer function can be disastrous for independently noisy sensor readings.
For measured data, suitable fitting, smoothing, differentiation-aware regularisation or a model-based derivative may be more appropriate. Any smoothing choice introduces assumptions and possible bias, so uncertainty and scale must be reported rather than hidden behind a smooth curve.
Contents · Previous · Next · Subject root
Automatic differentiation propagates derivatives through a sequence of elementary operations according to the chain rule. For differentiable implemented models, it can avoid the finite-difference step-size trade-off, though nondifferentiable branches, memory demands and software semantics still matter.
For suitable analytic functions with a valid complex extension, the complex-step method estimates a derivative from the imaginary part of f(x+ih) divided by h. It avoids subtracting two nearly equal real function values, allowing very small steps in favourable settings.
That trick does not automatically apply to functions involving arbitrary absolute values, branching, clipping, non-holomorphic operations or code that discards imaginary components. A successful complex-step derivative requires the function implementation to preserve the relevant analytic structure.
For a black-box model available only through real-valued evaluations, finite differences may remain the practical choice. In that case, use clear error tests and step-size studies rather than assuming a more sophisticated method can be called without verifying its requirements.
Contents · Previous · Next · Subject root
Exercise A: For f(x)=sin x at x=1, which is typically more accurate at h=0.1, forward or central difference? The central method, because its leading truncation error is O(h²), and the computed errors are about 0.000900 and 0.04294, respectively.
Exercise B: A student reduces h from 10⁻⁵ to 10⁻¹² in the central formula and obtains a much worse answer. Does that contradict the derivative limit? No. The limit presumes exact arithmetic; the machine evaluation introduces cancellation and rounding error.
Exercise C: The forward derivative at x=1 with h=10⁻¹⁶ becomes zero. Explain one cause. In binary64, the represented value of 1+h can equal 1, making the function evaluations identical and their difference exactly zero.
Exercise D: The function is a noisy physical measurement rather than an exact mathematical sin(x). Is h=10⁻⁵ automatically best? No. The appropriate scale depends on sensor precision, noise correlations, smoothness assumptions and how the derivative will be used.
Contents · Previous · Next · Subject root
First identify the desired derivative, its units and whether f is a symbolic function, differentiable program or set of noisy observations. Choose a method whose assumptions fit that object. For real finite differences, evaluate several steps on a log-spaced scale and look for an accuracy plateau rather than one lucky decimal.
Second, if a reference derivative is available, compute absolute errors explicitly. For a black-box problem, compare methods or use consistency checks that do not silently presume one numerical estimate is exact. Treat a sudden deterioration at tiny h as evidence of numerical limits.
Third, document the arithmetic precision, finite-difference formula, function-evaluation implementation, step-size rule and any smoothing or scaling. Do not equate a printed row of digits with accuracy unless the method and error checks support it.
Continue from The Gold Standard of Scientific Computing to stable quadratic-root evaluation and compensated summation. All three problems expose different ways information can disappear during an otherwise valid-looking calculation.
Applied numerical methods lab
A numerical table as a debugging tool
At h=0.1, the forward and central estimates differ by about four hundredths; the difference mostly reflects truncation. At h=0.0001, central difference is extremely close to cos1. At h=10⁻¹², central difference becomes less reliable, illustrating the roundoff branch of the error curve.
Always store the chosen h alongside the estimate and an independent reference when one exists. A table of results without step sizes conceals the experiment’s main variable. This is especially important when a later scientific article cites a derivative as if it were obtained directly from a mathematical expression.
Do not infer that the best-looking row is an exact optimum for every computer. Different implementations of sin and higher-precision arithmetic may change the last few rows; the structural pattern and need for a suitable h remain.
Derive a forward-difference warning from Taylor expansion
Expand sin(x+h) as sinx+hcosx−(h²/2)sinx plus higher powers. Dividing the difference by h leaves cosx−(h/2)sinx+O(h²). This explicitly predicts that the forward estimate at x=1 should be below cos1 for positive small h.
The leading sign is a useful quick check. If code returns a forward estimate well above cos1 at a moderate positive h, investigate angle units, argument ordering and how the function is evaluated before blaming fine-scale rounding.
The central difference cancels the first-order bias. This advantage comes from symmetry, not from averaging two ordinary forward differences after they were rounded in an arbitrary way.
Do not differentiate an instrument’s quantisation noise
Suppose an instrument rounds a physical measurement to the nearest 0.01 unit. If you choose changes in x so small that f(x+h) and f(x) frequently round to the same recorded level, the forward difference may report zero over a whole region even though the underlying physical function has nonzero slope.
When the reading finally jumps one quantisation level, the resulting quotient can spike. Shrinking h further amplifies these staircase effects instead of recovering the true smooth derivative.
Model-based fitting over a physically defensible window can improve such estimates, but only if the fitted assumptions are appropriate. Always report the measurement resolution and estimation window together with a derivative inferred from quantised data.
A practical classroom investigation
Provide students with the formula f(x)=sin x and an analytic reference cos1. Ask them to compute forward and central differences at h=0.1,0.01,0.0001,0.00001 and 0.000000000001, retaining the order of evaluations.
Students should plot absolute error against h on logarithmic axes, identify the decreasing truncation-error region and the rising rounding-error region, and mark where the central method performs best in their environment.
Then repeat the experiment after rounding every function value to four decimal places. Their task is to explain why the best h changes when input evaluation accuracy is degraded. The lesson transfers from textbook calculus to real measurement systems without demanding an unsupported claim of universal machine precision.
Frequently asked questions
Is smaller h always better?
No. Truncation error generally falls as h shrinks, but floating-point roundoff and measured-data noise can be amplified by division by very small h.
Is central difference always best?
Not universally. It has smaller leading truncation error for smooth functions, but boundary constraints, evaluation cost, noise and rounding can make another method more appropriate.
Can automatic differentiation avoid step-size choice?
For supported differentiable programs it can propagate derivatives without a finite-difference h, though it comes with implementation and differentiability requirements.
Sources and connected reading
All example values are computations using the stated analytic test function and arithmetic model; real measured data require their own uncertainty checks.
- Numerical Analysis Notes — Differentiation and truncation versus roundoff
- NumPy Reference — numpy.gradient and finite differences
- The Gold Standard of Scientific Computing — subject root
Advanced numerical lab: representable steps, extrapolation and noisy observations
When the requested step is not representable at the input’s scale
Floating-point step choice depends on the magnitude of x, not just the number of zeros in h. In binary64 near x=10¹², adjacent representable values are roughly 0.000122 apart. The mathematical quantity 10¹²+0.00001 therefore rounds to the same represented number as 10¹².
If a program attempts [f(x+h)−f(x)]/h with x=10¹² and h=10⁻⁵, the actual represented arguments can be identical even before f is called. A zero numerator follows, and the apparent derivative may be zero even when the true derivative is not.
This problem illustrates the difference between a nominal h entered by the programmer and the actual change in the represented input. A more careful numerical implementation may choose a step relative to the input’s scale and sometimes compute the effective representable increment rather than assume it exactly equals the decimal h.
Changing coordinate units can reveal the issue: x measured in metres and x measured in kilometres have different magnitudes, so a constant numerical h has different physical meaning. Choose step rules with dimensional scales and function sensitivity in mind.
Richardson extrapolation can cancel an error term
For a sufficiently smooth function, the centred difference D(h) has leading error proportional to h². Two centred estimates at h and h/2 can be combined as [4D(h/2)−D(h)]/3. This removes the leading h² contribution in ideal Taylor-series arithmetic.
At x=1 for sin(x), the centred estimate using h=0.1 is about 0.53940225217 and the estimate using h=0.05 is about 0.54007720805. Their Richardson combination is about 0.54030219334, close to cos1≈0.54030230587 and much better than either original estimate.
The absolute error of this extrapolated value is about 1.13×10⁻⁷ in the representative evaluation. That does not guarantee a higher-order combination is always more accurate at extremely small steps, because subtracting estimates can amplify their existing roundoff or noise.
Extrapolation is a reasoned tool when its assumed error expansion applies. If f is nonsmooth, the evaluation is noisy or the difference quotient is already roundoff-dominated, a formally higher order formula can provide a less reliable number.
Why the complex-step trick can work for sin(x)
For an analytic complex extension of sine, the identity sin(x+ih)=sin x cosh h+i cos x sinh h gives an imaginary part equal to cos x sinh h. Divide by h to obtain cos x ×[sinh h/h], which approaches cos x as h approaches zero.
Crucially, the imaginary component is extracted without subtracting two nearly equal real values. In a suitable implementation it can therefore permit far smaller steps than ordinary finite differences before cancellation becomes limiting. The method still has truncation effects and eventual representability or underflow limits.
This approach relies on a valid complex-analytic extension through the code path. An absolute value, comparison-driven branch, real-only function call, clipping operation or explicit discard of imaginary parts can invalidate the method or make its result misleading.
Before using a complex-step derivative in an engineering program, verify that the implementation accepts and preserves complex numbers and that the mathematical operation is appropriate for the desired derivative. It is not a universal workaround for arbitrary software.
A practical noise-resolution investigation
Imagine an instrument whose readings are rounded to the nearest 0.01 measurement unit. A smooth true f(x) may still produce an identical recorded value at x and x+h for many small h choices. The forward quotient then returns zero, not because the underlying derivative vanished, but because quantisation hid the difference.
When h becomes just large enough to cross a quantisation threshold, the recorded numerator can jump by 0.01. Divide that jump by a tiny h and the estimated derivative can spike dramatically. Making h smaller does not recover the lost smooth information; it may make the staircase artefact worse.
A useful classroom experiment records f(x)=sin x after rounding every value to four decimal places. Repeat the forward and centred calculations over progressively smaller h. Students can compare that error curve with the full-precision experiment and explain why the best step shifts.
For measured data, fit a physically defensible local model or use an appropriate regularisation strategy, but report its window, assumptions and uncertainty. Numerical sophistication should not conceal that the derivative came from limited-resolution observations.