# Workbook 5 — Finite-difference approximations from Taylor series
> [!abstract] The problem
> Starting from the Taylor expansions of $f$ about $x_0$ evaluated at $x_0\pm\Delta x$,
> $f(x_0\pm\Delta x)=f(x_0)\pm f'(x_0)\,\Delta x+\tfrac{f''(x_0)}{2!}\Delta x^2\pm\tfrac{f'''(x_0)}{3!}\Delta x^3+\mathcal{O}(\Delta x^4),$
> derive (1) the **forward-difference** approximation to $f'(x_0)$, (2) the **centered-difference** approximation to $f'(x_0)$, and (3) the second-order **central approximation to $f''(x_0)$**. State each formula's [order of accuracy](https://en.wikipedia.org/wiki/Order_of_accuracy), and explain in one sentence why the centered formulas beat the forward one.
---
## 0. What are we actually being asked to do?
A computer never has $f'$ or $f''$ as formulas — it has function *values* on a grid. [Finite differences](https://en.wikipedia.org/wiki/Finite_difference) are the bridge: combine a few nearby samples of $f$ to estimate a derivative. This problem asks you to *derive* three such formulas honestly from Taylor series and to read off, from the leftover terms, **how fast each one converges** as the grid spacing $\Delta x$ shrinks.
The whole exercise is algebra on the two expansions above. Every formula, and every error estimate, falls out of **adding or subtracting** them.
> [!tip] The one idea to hold onto
> Add the two expansions and the *odd* powers of $\Delta x$ cancel; subtract them and the *even* powers cancel. That single symmetry is the reason centered formulas are more accurate — the cancellation removes the largest error term for free.
---
## 1. The two building blocks
Write both expansions about $x_0$ (all derivatives evaluated at $x_0$):
$f(x_0+\Delta x)=f+f'\Delta x+\tfrac{f''}{2}\Delta x^2+\tfrac{f'''}{6}\Delta x^3+\mathcal{O}(\Delta x^4)\qquad(\mathrm{A})$
$f(x_0-\Delta x)=f-f'\Delta x+\tfrac{f''}{2}\Delta x^2-\tfrac{f'''}{6}\Delta x^3+\mathcal{O}(\Delta x^4)\qquad(\mathrm{B})$
Notice the sign pattern: in (B) the **odd**-order terms flip sign while the **even**-order terms do not. That is the lever for everything below.
## 2. Forward difference for $f'$ — from (A) alone
Keep only $(\mathrm{A})$, isolate $f'\Delta x$, and divide by $\Delta x$:
$f(x_0+\Delta x)-f=f'\Delta x+\tfrac{f''}{2}\Delta x^2+\cdots
\;\Longrightarrow\;
\boxed{\,f'(x_0)\approx\dfrac{f(x_0+\Delta x)-f(x_0)}{\Delta x}\,}$
The discarded terms are
$\underbrace{\dfrac{f(x_0+\Delta x)-f(x_0)}{\Delta x}}_{\text{approximation}}-f'(x_0)=\tfrac{f''(x_0)}{2}\Delta x+\mathcal{O}(\Delta x^2).$
The leading error is proportional to $\Delta x^1$, so the forward difference is **first-order accurate**, $\mathcal{O}(\Delta x)$. Halve $\Delta x$ and the error roughly halves.
## 3. Centered difference for $f'$ — subtract (B) from (A)
Subtracting kills the even terms (the $f$ and $\tfrac{f''}{2}\Delta x^2$ pieces cancel):
$f(x_0+\Delta x)-f(x_0-\Delta x)=2f'\Delta x+\tfrac{f'''}{3}\Delta x^3+\mathcal{O}(\Delta x^5).$
Divide by $2\Delta x$:
$\boxed{\,f'(x_0)\approx\dfrac{f(x_0+\Delta x)-f(x_0-\Delta x)}{2\Delta x}\,},\qquad
\text{error}=-\tfrac{f'''(x_0)}{6}\Delta x^2+\mathcal{O}(\Delta x^4).$
The $\Delta x^1$ term is *gone*; the leading error is $\propto\Delta x^2$, so the centered difference is **second-order accurate**, $\mathcal{O}(\Delta x^2)$. Halve $\Delta x$ and the error drops by a factor of **four**.
## 4. Central second difference for $f''$ — add (A) and (B)
Adding kills the odd terms (the $f'\Delta x$ and $\tfrac{f'''}{6}\Delta x^3$ pieces cancel):
$f(x_0+\Delta x)+f(x_0-\Delta x)=2f+f''\Delta x^2+\tfrac{f''''}{12}\Delta x^4+\cdots.$
Solve for $f''$:
$\boxed{\,f''(x_0)\approx\dfrac{f(x_0+\Delta x)-2f(x_0)+f(x_0-\Delta x)}{\Delta x^2}\,},\qquad
\text{error}=\tfrac{f''''(x_0)}{12}\Delta x^2+\mathcal{O}(\Delta x^4).$
Again the leading error is $\propto\Delta x^2$: **second-order accurate**, $\mathcal{O}(\Delta x^2)$.
## 5. Why centered beats forward (the one-sentence answer)
**Because the symmetric stencil cancels the leading $\mathcal{O}(\Delta x)$ error term** — subtracting (A) and (B) annihilates the even-power terms, pushing the first surviving error from $\Delta x$ up to $\Delta x^2$, so the centered estimate improves four-fold (not two-fold) each time the grid is refined.
---
## 6. Numerical confirmation
Take $f(x)=e^{-x^2}$ at $x_0=0.5$, where $f'(x_0)=-0.7788$ and $f''(x_0)=-0.7788$. Halving $\Delta x$ repeatedly:
| $\Delta x$ | forward $f'$ error | centered $f'$ error | central $f''$ error |
|---|---|---|---|
| $0.1000$ | $3.24\times10^{-2}$ | $6.46\times10^{-3}$ | $6.56\times10^{-4}$ |
| $0.0500$ | $1.78\times10^{-2}$ | $1.62\times10^{-3}$ | $1.63\times10^{-4}$ |
| $0.0250$ | $9.33\times10^{-3}$ | $4.06\times10^{-4}$ | $4.06\times10^{-5}$ |
| $0.0125$ | $4.77\times10^{-3}$ | $1.01\times10^{-4}$ | $1.01\times10^{-5}$ |
Read the columns as ratios: forward error **halves** each row (order 1), while both centered columns drop by **≈4** each row (order 2). On a log–log plot the slopes are exactly $1$ and $2$:
The one-line experiment in each language just fills this table. Runnable files live in `code/wb5/`; Python is run-verified, the others follow the course porting conventions.
> [!example]- Python (reference)
> ```python
> import numpy as np
> f = lambda x: np.exp(-x**2)
> x0 = 0.5
> fp_true = -2*x0*np.exp(-x0**2)
> fpp_true = 2*(2*x0**2 - 1)*np.exp(-x0**2)
> print(f"{'dx':>8}{'fwd':>12}{'cen':>12}{'2nd':>12}")
> for k in range(7):
> h = 0.4 / 2**k
> fwd = (f(x0+h) - f(x0)) / h # O(dx) -> f'
> cen = (f(x0+h) - f(x0-h)) / (2*h) # O(dx^2) -> f'
> sec = (f(x0+h) - 2*f(x0) + f(x0-h)) / h**2 # O(dx^2) -> f''
> print(f"{h:8.4f}{abs(fwd-fp_true):12.2e}{abs(cen-fp_true):12.2e}{abs(sec-fpp_true):12.2e}")
> ```
> *Look for:* the forward column halving while the centered columns quarter.
> [!example]- MATLAB
> ```matlab
> f = @(x) exp(-x.^2); x0 = 0.5;
> fpTrue = -2*x0*exp(-x0^2);
> fppTrue = 2*(2*x0^2 - 1)*exp(-x0^2);
> for k = 0:6
> h = 0.4 / 2^k;
> fwd = (f(x0+h) - f(x0)) / h;
> cen = (f(x0+h) - f(x0-h)) / (2*h);
> sec = (f(x0+h) - 2*f(x0) + f(x0-h)) / h^2;
> fprintf('%8.4f%12.2e%12.2e%12.2e\n', h, abs(fwd-fpTrue), abs(cen-fpTrue), abs(sec-fppTrue))
> end
> ```
> *Look for:* `f(x0+h)` calls the handle at a scalar — no dots needed here.
> [!example]- R
> ```r
> f <- function(x) exp(-x^2); x0 <- 0.5
> fpTrue <- -2*x0*exp(-x0^2)
> fppTrue <- 2*(2*x0^2 - 1)*exp(-x0^2)
> for (k in 0:6) {
> h <- 0.4 / 2^k
> fwd <- (f(x0+h) - f(x0)) / h
> cen <- (f(x0+h) - f(x0-h)) / (2*h)
> sec <- (f(x0+h) - 2*f(x0) + f(x0-h)) / h^2
> cat(sprintf("%8.4f%12.2e%12.2e%12.2e\n", h, abs(fwd-fpTrue), abs(cen-fpTrue), abs(sec-fppTrue)))
> }
> ```
> *Look for:* `sprintf` with `%12.2e` lines the columns up like the table above.
> [!example]- Mathematica (numerical only)
> ```wolfram
> f = Exp[-#^2] &; x0 = 0.5;
> fpTrue = -2 x0 Exp[-x0^2];
> fppTrue = 2 (2 x0^2 - 1) Exp[-x0^2];
> Table[With[{h = 0.4/2^k},
> {h, Abs[(f[x0+h]-f[x0])/h - fpTrue],
> Abs[(f[x0+h]-f[x0-h])/(2 h) - fpTrue],
> Abs[(f[x0+h]-2 f[x0]+f[x0-h])/h^2 - fppTrue]}], {k, 0, 6}] // Grid
> ```
> *Look for:* keep it numerical — `x0 = 0.5` (a machine number), not exact `1/2`.
---
# Sanity checks and context
> [!note] What the rest of this document is
> The derivations (§§2–4), the orders of accuracy, and the one-sentence reason (§5) are the **complete answer**. Below: where these formulas come from historically, why they matter, common pitfalls, and self-test questions.
## 7. Where these come from — and the limit you already know
The forward difference is just the [difference quotient](https://en.wikipedia.org/wiki/Difference_quotient) from the [limit definition of the derivative](https://en.wikipedia.org/wiki/Derivative#Definition), *stopped before the limit is taken*:
$f'(x_0)=\lim_{\Delta x\to0}\frac{f(x_0+\Delta x)-f(x_0)}{\Delta x}.$
Calculus takes $\Delta x\to0$; computing keeps $\Delta x$ small but finite and pays for it with the $\mathcal{O}(\Delta x)$ error you just derived. Finite differences are "calculus with the limit left unfinished," and Taylor's theorem is the exact bookkeeping of the unfinished part.
## 8. Why this matters
These three stencils are the atoms of numerical [differentiation](https://en.wikipedia.org/wiki/Numerical_differentiation) and, later, of solving differential equations: replace every derivative in an equation by its difference stencil and an unsolvable calculus problem becomes solvable algebra on a grid. The second difference of §4 is exactly the numerical $f''$ that detects concavity and inflection points from data (Workbook 1's geometry, done by a computer), and it is the operator behind discretizing the heat and wave equations. Understanding *order of accuracy* is what lets you choose a grid: to cut a second-order error by $100\times$, refine by only $10\times$, whereas first-order would demand $100\times$ the points.
> [!warning] Round-off has the last word
> The error analysis above is the **truncation** error, which shrinks with $\Delta x$. In floating point, subtracting nearly equal values ($f(x_0+\Delta x)\approx f(x_0)$) and dividing by a tiny $\Delta x$ *amplifies* [round-off](https://en.wikipedia.org/wiki/Round-off_error). Past some sweet-spot $\Delta x$, shrinking further makes the estimate **worse**. Truncation error says "smaller $\Delta x$ is better"; round-off says "not too small." The optimum is a balance.
## 9. Common pitfalls
> [!warning] Watch out for these
> - **Forgetting to divide by the right power.** Forward/centered $f'$ divide by $\Delta x$ and $2\Delta x$; the second difference divides by $\Delta x^{2}$. A wrong power silently rescales the answer.
> - **Miscounting the order.** The order is the exponent of $\Delta x$ in the *leading surviving* error term, not the highest term you wrote down. Centered $f'$ is $\mathcal{O}(\Delta x^2)$ because the $\Delta x^1$ term cancelled.
> - **The centered $f'$ does not use $f(x_0)$.** Only the two neighbors enter; the center value drops out in the subtraction. (It *does* appear in the second difference.)
> - **Sign slips in (B).** Every odd derivative flips sign at $x_0-\Delta x$. Getting one sign wrong destroys the cancellation that gives second order.
> - **Chasing $\Delta x\to0$ numerically.** Because of §8's round-off, an ever-smaller $\Delta x$ eventually increases error — the opposite of the truncation prediction.
## 10. Checkpoint questions
1. Derive the **backward** difference $\dfrac{f(x_0)-f(x_0-\Delta x)}{\Delta x}$ from (B) and state its order. Why does it match the forward difference's order but not its error sign?
2. Add (A) and (B) but keep the $\Delta x^4$ term. What is the *exact* leading error coefficient of the second-difference formula for $f''$, in terms of $f''''(x_0)$?
3. For $f(x)=e^{-x^2}$ at $x_0=0.5$, predict the centered-$f'$ error at $\Delta x=0.05$ from the coefficient $-\tfrac{f'''(x_0)}{6}\Delta x^2$, and compare with the table in §6.
4. You want $f'$ accurate to $10^{-6}$ using the centered formula. Starting from the $\Delta x=0.1$ error in §6, roughly what $\Delta x$ do you need — and what would first-order forward differencing demand instead?