Engineering Mathematics for NM: Numerical Differentiation, and the Exact and Bernoulli Equations

The Naval Architecture and Marine Engineering (NM) paper’s Section 1 names two things the shared Engineering Mathematics chapters do not teach, and this chapter is written for exactly those two. The first is numerical methods for differentiation: the shared Numerical Methods chapter integrates by the trapezoidal and Simpson rules but never differentiates a table, so the forward, backward and central difference formulas, their truncation errors, the second-derivative formula, the Newton forward-difference derivative and Richardson extrapolation are here. The second is the non-linear half of "linear, non-linear, first and higher order ordinary differential equations": the borrowed first-order chapter teaches the separable, homogeneous and linear forms, so the exact equation, its integrating factors and the Bernoulli equation, which a substitution turns linear, are here too. Everything else in the NM Section 1 list is in the shared chapters. Both topics are small, quick and reliably examinable: a difference quotient on three tabulated values, or a test for exactness, is a one-mark question that costs a minute.

1. Difference formulas from the Taylor series

Every numerical-differentiation formula is a truncated Taylor series. Writing f(x + h) = f(x) + hf′(x) + (h²/2)f″(x) + (h³/6)f‴(x) + … and solving for f′ gives the forward difference f′(x) ≈ [f(x + h) − f(x)]/h, whose leading error is −(h/2)f″(ξ): first order, O(h). Expanding f(x − h) instead gives the backward difference [f(x) − f(x − h)]/h, also O(h). Subtracting the two expansions cancels the f″ terms and leaves the central difference f′(x) ≈ [f(x + h) − f(x − h)]/(2h), whose leading error is −(h²/6)f‴(ξ): second order, O(h²). Adding them cancels f′ and gives the second-derivative formula f″(x) ≈ [f(x + h) − 2f(x) + f(x − h)]/h², also O(h²), with leading error −(h²/12)f⁗.

The four formulas and their errors
FormulaExpressionTruncation error
Forward difference[f(x + h) − f(x)]/h−(h/2)f″, order h
Backward difference[f(x) − f(x − h)]/h+(h/2)f″, order h
Central difference[f(x + h) − f(x − h)]/(2h)−(h²/6)f‴, order h²
Second derivative[f(x + h) − 2f(x) + f(x − h)]/h²−(h²/12)f⁗, order h²
🧠 Check the formula on a function whose derivative you know
For f(x) = x³ at x = 1 with h = 0.1: forward gives (1.331 − 1)/0.1 = 3.31, central gives (1.331 − 0.729)/0.2 = 3.01, and the exact value is 3. The forward error 0.31 is close to (h/2)f″ = 0.05 × 6 = 0.3; the central error 0.01 is exactly (h²/6)f‴ = 0.01 × 6/6 = 0.01, because the Taylor series of a cubic stops at the third term.

2. Differentiating a table: Newton’s formula and Richardson extrapolation

When f is known only at equally spaced points x₀, x₀ + h, x₀ + 2h, …, build the forward-difference table Δf₀ = f₁ − f₀, Δ²f₀ = Δf₁ − Δf₀, and so on. Differentiating Newton’s forward interpolation polynomial at x₀ gives f′(x₀) ≈ (1/h)[Δf₀ − Δ²f₀/2 + Δ³f₀/3 − Δ⁴f₀/4 + …]. Keeping only the first term is the forward difference; each further term raises the order. If the data come from a polynomial of degree n, the nth differences are constant and the (n + 1)th vanish, so the series terminates and the derivative is exact.

Richardson extrapolation removes the leading error term by combining two step sizes. If D(h) has error ch², then D(h/2) has error ch²/4, and eliminating c gives the improved estimate D = [4D(h/2) − D(h)]/3, of order h⁴. For a first-order formula (error ch) the combination is 2D(h/2) − D(h). The same idea applied to the trapezoidal rule is Romberg integration.

⚠️ A smaller step is not always better
Truncation error falls as h shrinks, but the round-off error in f(x + h) − f(x), divided by h, grows as 1/h. The total error therefore has a minimum at an optimum step, below which the answer gets worse. A question asking what happens as h → 0 on a computer with finite precision wants this answer, not "the error goes to zero".

3. Exact equations and integrating factors

The first-order equation M(x, y) dx + N(x, y) dy = 0 is exact when its left side is the total differential of some function u(x, y), which happens exactly when ∂M/∂y = ∂N/∂x (on a simply connected region). The solution is u(x, y) = C, found by integrating M with respect to x (holding y constant), then adding the terms of N that contain no x. For (2xy + 3) dx + (x² + 4y) dy = 0: ∂M/∂y = 2x = ∂N/∂x, so it is exact; ∫M dx = x²y + 3x, and N contributes the y-only term 4y, whose integral is 2y². The solution is x²y + 3x + 2y² = C.

An equation that is not exact can often be made so by an integrating factor μ. If (∂M/∂y − ∂N/∂x)/N is a function of x alone, μ = exp∫[(∂M/∂y − ∂N/∂x)/N] dx; if (∂N/∂x − ∂M/∂y)/M is a function of y alone, μ = exp∫[(∂N/∂x − ∂M/∂y)/M] dy. Some forms are recognised on sight: y dx − x dy is made exact by 1/x², 1/y² or 1/(xy), giving −d(y/x), d(x/y) and d(ln x − ln y) respectively. The linear equation dy/dx + Py = Q is the special case whose integrating factor is e∫P dx.

⚠️ Test exactness with the cross-derivatives, not the same-variable ones
The test is ∂M/∂y against ∂N/∂x — M is differentiated by y and N by x. Comparing ∂M/∂x with ∂N/∂y is the commonest slip and proves nothing.

4. The Bernoulli equation

The Bernoulli equation dy/dx + P(x)y = Q(x)yⁿ is non-linear for n ≠ 0, 1, but the substitution v = y1−n makes it linear: dividing by yⁿ and using dv/dx = (1 − n)y−n dy/dx gives dv/dx + (1 − n)P v = (1 − n)Q, solved by the integrating factor e(1−n)∫P dx. For n = 0 it is already linear and for n = 1 it is separable, which is why those two are excluded.

Worked example, the logistic equation dy/dx = y − y², that is dy/dx − y = −y² (n = 2). With v = 1/y: dv/dx + v = 1, so v = 1 + Ce−x. If y(0) = 0.5 then v(0) = 2 and C = 1, so y = 1/(1 + e−x), which rises from 0.5 towards 1; at x = ln 2 it is 1/(1 + 0.5) = 0.667.

  • Identify n from the power of y on the right; the substitution is v = y1−n, not v = yⁿ.
  • After substituting, the coefficient of v is (1 − n)P and the right side is (1 − n)Q; forgetting the factor (1 − n) is the usual error.
  • Convert the initial condition to v before finding the constant, then convert back to y at the end.

5. Where the rest of NM Section 1 is taught

NM Section 1, mapped onto the chapters that teach it
Syllabus wordingWhere it is taught
Determinants and matrices, systems of linear equations, eigen values and eigen vectorsLinear Algebra chapter
Functions, definite and indefinite integralsThe borrowed single-variable calculus chapters
Chain rules, partial and directional derivatives, gradient, divergence, curl, line, surface and volume integrals, Stokes, Gauss and GreenMultivariable Calculus and Vector Calculus chapters
First-order linear equations; higher-order ODEs, PDEs, separation of variables, Laplace transformationThe borrowed first-order chapter and the Differential Equations chapter
Fourier seriesThe first section of the Transforms chapter only
Analytical functions of complex variables, complex analysisComplex Variables chapter
Numerical integration; probability and statisticsNumerical Methods (its quadrature section) and Probability and Statistics chapters
Numerical differentiation; exact and Bernoulli equationsThis chapter
ℹ️ Numerical integration in hydrostatics
Naval architecture integrates offsets with Simpson’s rules — waterplane areas, sectional areas, volumes and moments. The quadrature itself is taught in the shared Numerical Methods chapter; its use on half-breadths is in the NM naval architecture chapter.

Key takeaways

  • Forward and backward differences are first order (error −(h/2)f″ and +(h/2)f″); the central difference is second order (error −(h²/6)f‴).
  • f″ ≈ [f(x + h) − 2f(x) + f(x − h)]/h², second order; from a table, f′(x₀) ≈ (1/h)[Δ − Δ²/2 + Δ³/3 − …].
  • Richardson for an h² method: D = [4D(h/2) − D(h)]/3; round-off makes a very small h worse, not better.
  • M dx + N dy = 0 is exact iff ∂M/∂y = ∂N/∂x; y dx − x dy takes 1/x², 1/y² or 1/(xy) as an integrating factor.
  • Bernoulli dy/dx + Py = Qyⁿ: put v = y1−n to get dv/dx + (1 − n)Pv = (1 − n)Q.

Practice questions (12)

Attempt each one before opening the answer. Every explanation names the tempting wrong option as well as the right one, because that is where marks are lost.

  1. The truncation error of the central difference approximation f′(x) ≈ [f(x + h) − f(x − h)]/(2h) is of order

    1. h²
    2. h
    3. h³
    4. h⁴
    Show answer

    Answer: A — h²

    Subtracting the Taylor expansions of f(x + h) and f(x − h) cancels the even-order terms, so the leading error is −(h²/6)f‴: order h². The one-sided forward and backward differences are order h; order h⁴ needs Richardson extrapolation or a five-point formula.
  2. For f(x) = x³, the forward-difference estimate of f′(1) with h = 0.1, to two decimal places, is ____.

    Numerical answer — type the value.

    Show answer

    Answer: 3.31

    [f(1.1) − f(1)]/0.1 = (1.331 − 1)/0.1 = 3.31, against the exact 3. The error 0.31 is roughly (h/2)f″(1) = 0.05 × 6 = 0.3, as the first-order error term predicts. The central difference would give 3.01.
  3. For f(x) = x⁴, the central-difference estimate of the second derivative at x = 1 with h = 0.1, to two decimal places, is ____.

    Numerical answer — type the value.

    Show answer

    Answer: 12.02

    [f(1.1) − 2f(1) + f(0.9)]/h² = (1.4641 − 2 + 0.6561)/0.01 = 0.1202/0.01 = 12.02, against the exact f″(1) = 12. The error 0.02 equals (h²/12)f⁗ = (0.01/12) × 24 = 0.02, the formula’s leading error term.
  4. Central-difference estimates of f′(1) for f(x) = x³ are 3.04 with h = 0.2 and 3.01 with h = 0.1. The Richardson-extrapolated estimate, to two decimal places, is ____.

    Numerical answer — type the value.

    Show answer

    Answer: 3

    The central difference has error ch², so D = [4D(h/2) − D(h)]/3 = (4 × 3.01 − 3.04)/3 = (12.04 − 3.04)/3 = 3.00, the exact value. Averaging the two estimates gives 3.025; using the first-order combination 2D(h/2) − D(h) gives 2.98.
  5. A function is tabulated at x = 0, 1, 2, 3 as f = 2, 5, 10, 17. Using Newton’s forward-difference formula for the derivative, f′(0) is ____.

    Numerical answer — type the value.

    Show answer

    Answer: 2

    First differences 3, 5, 7; second differences 2, 2; third difference 0. With h = 1, f′(0) = Δf₀ − Δ²f₀/2 + Δ³f₀/3 = 3 − 1 + 0 = 2. The data are f = x² + 2x + 2, whose derivative at 0 is exactly 2. Stopping at the first difference gives 3.
  6. Which statements about numerical differentiation are correct?

    1. halving h roughly halves the truncation error of the forward difference
    2. halving h roughly quarters the truncation error of the central difference
    3. round-off error grows as h becomes very small
    4. the central difference is exact for every cubic polynomial
    Show answer

    Answer: A — halving h roughly halves the truncation error of the forward difference; B — halving h roughly quarters the truncation error of the central difference; C — round-off error grows as h becomes very small

    Forward error is proportional to h and central error to h², so halving h halves one and quarters the other; the subtraction f(x + h) − f(x) loses significant digits as h shrinks, so round-off grows like 1/h. The central difference is exact only up to quadratics: for a cubic its error is (h²/6)f‴, a non-zero constant.
  7. The general solution of the exact equation (2xy + 3) dx + (x² + 4y) dy = 0 is

    1. x²y + 3x + 2y² = C
    2. x²y + 3x + 4y² = C
    3. 2xy + 3x + 2y² = C
    4. x²y + 3y + 2x² = C
    Show answer

    Answer: A — x²y + 3x + 2y² = C

    ∂M/∂y = 2x = ∂N/∂x, so it is exact. Integrating M in x gives x²y + 3x; the only term of N free of x is 4y, which integrates to 2y². Writing 4y² forgets to integrate, and 2xy + … keeps M instead of its integral.
  8. The equation (3x²y + ky²) dx + (x³ + 4xy) dy = 0 is exact for k = ____.

    Numerical answer — type the value.

    Show answer

    Answer: 2

    ∂M/∂y = 3x² + 2ky and ∂N/∂x = 3x² + 4y; equal for all x and y when 2k = 4, so k = 2. Comparing ∂M/∂x with ∂N/∂y instead compares 6xy with 4x and gives no constant k at all.
  9. Which of the following are integrating factors of y dx − x dy = 0?

    1. 1/x²
    2. 1/y²
    3. 1/(xy)
    4. 1/(x + y)
    Show answer

    Answer: A — 1/x²; B — 1/y²; C — 1/(xy)

    Multiplying by 1/x² gives (y/x²) dx − (1/x) dy = −d(y/x); by 1/y², (1/y) dx − (x/y²) dy = d(x/y); by 1/(xy), dx/x − dy/y = d(ln x − ln y). With 1/(x + y), ∂M/∂y = x/(x + y)² but ∂N/∂x = −y/(x + y)², which differ, so it is not an integrating factor.
  10. The Bernoulli equation dy/dx + y = xy³ becomes linear under the substitution

    1. v = y⁻²
    2. v = y³
    3. v = y⁻³
    4. v = y²
    Show answer

    Answer: A — v = y⁻²

    With n = 3 the substitution is v = y1−n = y−2, giving dv/dx − 2v = −2x, which is linear. Using v = yⁿ = y³ or v = y−n = y−3 does not cancel the non-linearity; v = y² has the right magnitude of exponent but the wrong sign.
  11. The solution of dy/dx = y − y² with y(0) = 0.5, evaluated at x = ln 2, to three decimal places, is ____.

    Numerical answer — type the value.

    Show answer

    Answer: 0.667

    Bernoulli with n = 2: v = 1/y gives dv/dx + v = 1, so v = 1 + Ce−x; v(0) = 2 gives C = 1. At x = ln 2, v = 1 + 1/2 = 1.5 and y = 1/1.5 = 0.667. Leaving the answer as v = 1.5 forgets to convert back to y.
  12. For the Bernoulli equation dy/dx + Py = Qyⁿ, the values of n for which it is already linear or separable, and so needs no substitution, are

    1. n = 0 and n = 1
    2. n = 1 and n = 2
    3. n = −1 and n = 1
    4. n = 0 and n = 2
    Show answer

    Answer: A — n = 0 and n = 1

    With n = 0 the equation is dy/dx + Py = Q, linear; with n = 1 it is dy/dx = (Q − P)y, separable. Every other n, including 2 and −1, needs v = y1−n.