Units and scaling

Every worked problem in this book is posed on a domain like $[0, 1]$, with solutions of size about 1. That is not a coincidence or a simplification for teaching. It is a choice made before training, and a real engineering problem rarely arrives in that form. This section shows what goes wrong in real units and how to fix it on paper, before any code runs.

A steel rod in real units

A carbon-steel rod of length $L = 0.5$ m starts with a hump of temperature, $300\,{}^\circ\mathrm{C}$ in the middle, and its ends are held at $u_e = 20\,{}^\circ\mathrm{C}$. The thermal diffusivity of steel is of the order of $10^{-5}\ \mathrm{m^2/s}$; we take $a = 1.2\times10^{-5}\ \mathrm{m^2/s}$, typical of a high-carbon steel. We want the temperature for the first hour. The heat equation and its data are

$$u_t = a\,u_{xx}, \qquad u(0, t) = u(L, t) = u_e, \qquad u(x, 0) = u_e + \Delta T\sin(\pi x/L), \quad \Delta T = 280\ \mathrm{K}.$$

The numbers span many orders of magnitude: $x$ up to 0.5, $t$ up to 3600, $u$ around a few hundred, and $a$ of order $10^{-5}$. The equation itself suggests the natural unit of time:

import numpy as np
import torch

L, a = 0.5, 1.2e-5            # m, m^2/s
u_e, dT = 20.0, 280.0          # degC at the ends, K of rise in the middle at t = 0
t_end = 3600.0                 # s: one hour
tau = L**2 / a                 # the diffusion time of the rod
print(f"diffusion time L^2/a = {tau:.0f} s ({tau / 3600:.1f} h);  one hour is {t_end / tau:.3f} of it")

diffusion time L^2/a = 20833 s (5.8 h);  one hour is 0.173 of it

The rod needs hours to cool, and the hour we care about is about a sixth of its diffusion time. Keep that number, $0.173$, in mind.

What goes wrong in raw units

The inputs saturate the network. A freshly initialised layer (Xavier weights, Section 2.5) multiplies its inputs by numbers of order 0.2. Feed it $t$ in seconds and the pre-activations are in the hundreds, deep in the flat part of tanh, where the slope is zero to many digits (Section 2.3):

torch.manual_seed(0)
layer = torch.nn.Linear(2, 32)
torch.nn.init.xavier_normal_(layer.weight)
torch.nn.init.zeros_(layer.bias)

pts = torch.rand(1000, 2) * torch.tensor([L, t_end])              # x in [0, L], t in [0, t_end]
for name, inputs in [("raw (metres, seconds)", pts),
                     ("scaled to [-1, 1]", 2 * pts / torch.tensor([L, t_end]) - 1)]:
    z = layer(inputs)
    print(f"{name:22s} median |z| {z.abs().median():8.2f}   saturated (|z| > 3): {(z.abs() > 3).float().mean():.0%}")

raw (metres, seconds)  median |z|   174.38   saturated (|z| > 3): 98%
scaled to [-1, 1]      median |z|     0.08   saturated (|z| > 3): 0%

With raw inputs almost every neuron of the first layer is saturated before training begins, so almost no gradient reaches the weights. Scaled to $[-1, 1]$, none is. This part of the problem has a standard remedy, and the Studio applies it for you: the generated network rescales its inputs from the Domain's bounds to $[-1, 1]$ before the first layer, inside the network, so that derivatives with respect to $x$ and $t$ still come out in the original units (Section 5.3 shows why that matters).

The loss terms are of wildly different sizes. That part the Studio does not fix, because it depends on the equation and the data. Ask what the same relative mistake costs in each term: a network whose initial temperature rise is off by 1%, against one whose $u_t$ is off by 1% of the true cooling rate. Both can be computed from the exact solution, $u = u_e + \Delta T\,e^{-a\pi^2 t/L^2}\sin(\pi x/L)$:

x = np.linspace(0, L, 201)
X, T = np.meshgrid(x, np.linspace(0, t_end, 201))
u_t = -(a * np.pi**2 / L**2) * dT * np.exp(-a * np.pi**2 * T / L**2) * np.sin(np.pi * X / L)    # K per second
rise0 = dT * np.sin(np.pi * x / L)                                                              # K, at t = 0
pde, ic = np.mean((0.01 * u_t) ** 2), np.mean((0.01 * rise0) ** 2)
print(f"raw units:      residual term {pde:.1e}   initial-condition term {ic:.1e}   ratio {ic / pde:.0e}")

raw units:      residual term 2.5e-07   initial-condition term 3.9e+00   ratio 2e+07

A 1% error in the starting temperature costs ten million times more than a 1% error in the physics. An optimiser that minimises the sum will satisfy the initial and boundary data and treat the equation as an afterthought, and a loss curve that falls nicely will hide it. This is the imbalance that Chapter 10 is about, and in raw units it is built in before training starts.

Scaling the problem away

Introduce dimensionless variables, each measured in its natural unit:

$$\hat x = \frac{x}{L}, \qquad \hat t = \frac{t}{L^2/a}, \qquad \hat u = \frac{u - u_e}{\Delta T}.$$

By the chain rule, $\partial/\partial t = (a/L^2)\,\partial/\partial\hat t$ and $\partial/\partial x = (1/L)\,\partial/\partial\hat x$, so $u_t = (a\Delta T/L^2)\,\hat u_{\hat t}$ and $u_{xx} = (\Delta T/L^2)\,\hat u_{\hat x\hat x}$. Substituting,

$$\frac{a\,\Delta T}{L^2}\,\hat u_{\hat t} = a\,\frac{\Delta T}{L^2}\,\hat u_{\hat x\hat x} \quad\Longrightarrow\quad \hat u_{\hat t} = \hat u_{\hat x\hat x},$$

on $0 < \hat x < 1$, $0 < \hat t \le 0.173$, with $\hat u = 0$ at both ends and $\hat u(\hat x, 0) = \sin(\pi\hat x)$. Every constant has gone, the inputs are of order one, and so is the answer. The same 1% comparison now gives:

xh = x / L
Xh, Th = X / L, T / tau
uh_t = -np.pi**2 * np.exp(-np.pi**2 * Th) * np.sin(np.pi * Xh)
pde, ic = np.mean((0.01 * uh_t) ** 2), np.mean((0.01 * np.sin(np.pi * xh)) ** 2)
print(f"scaled units:   residual term {pde:.1e}   initial-condition term {ic:.1e}   ratio {ic / pde:.0e}")

scaled units:   residual term 1.4e-03   initial-condition term 5.0e-05   ratio 4e-02

The two terms are now within two orders of magnitude of each other, and it is the equation that weighs more. Once trained, the answer is turned back into degrees with $u = u_e + \Delta T\,\hat u(x/L,\ a t/L^2)$.

Dimensionless numbers tell you how hard it is

Scaling does more than tidy up. What is left after it is a handful of dimensionless numbers, and they measure the difficulty of the problem. The heat equation above had none left: every rod behaves the same in its own units. Other equations keep one:

The Burgers problem of Section 13.6 has $U = 1$, $L = 1$ and $\nu = 0.01/\pi$, so $\mathrm{Re} \approx 314$. The thickness of its shock is proportional to $\nu/U \approx 0.003$: at $t = 1$ the solution drops by about 1.4 over a distance of about 0.02, fifty times narrower than the domain. That single number explains why it needs a bigger network, ten thousand points and a second optimiser.

Exercises

1. Nondimensionalise the wave equation $u_{tt} = c^2 u_{xx}$ on $0 < x < L$. What is the natural time scale, and what equation is left?
2. Write Burgers' equation $u_t + u\,u_x = \nu u_{xx}$ with $\hat x = x/L$, $\hat t = tU/L$, $\hat u = u/U$. Show that only $\mathrm{Re} = UL/\nu$ remains.
3. An aluminium rod has $a \approx 9.7\times10^{-5}\ \mathrm{m^2/s}$. For $L = 0.5$ m, what is its diffusion time, and what range of $\hat t$ does one hour cover?

Recap

• In raw units the inputs saturate a tanh network and the loss terms differ by many orders of magnitude.
• Scale each variable by its natural unit ($L$, $L^2/a$, $\Delta T$, ...): the equation loses its constants, and the terms of the loss become comparable. The Studio scales the inputs; the output and the equation are yours to scale.
• What survives scaling is a set of dimensionless numbers, such as Péclet and Reynolds, and they predict how hard the problem will be.