Between two steps: lerp¶
A quantity drawn between two steps, linear interpolation: x = before + (after - before) alpha, alpha in [0, 1] the share of the interval from the step before to the last one that the frame is drawn at (0 draws the step before, 1 the last).
A view's frames outnumber its steps, so a camera and the shapes it shows move at the frame rate. Drawing only: nothing interpolated goes back into the simulation. Any unit, and any length of vector. Source: the linear interpolation of numerical analysis (Burden and Faires, Numerical Analysis, section 3.1). Named lerp since between became a word of the language. Checked against a reference in numpy, which has its own tests.
- Source
- linear interpolation (Burden and Faires, Numerical Analysis, section 3.1).
- Proof
- reference, invariant, limit.
lerp.trait
namespace lusine.render;
vector<float<_>, _> lerp(vector<float<_>, _> before, vector<float<_>, _> after, float alpha)
{
return before + (after - before) * alpha;
}
A quantity drawn between the step before and the last, before + (after - before) alpha, the frame's alpha the host's. Any unit and any length of vector: the program here binds it at 3, a ball's centre. One of the kernels lusine.render.table draws a table's shapes with.
import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
import matplotlib.pyplot as plt
from lusine import trait_run
F32 = np.float32
U = float(np.finfo(F32).eps) / 2 # binary32's unit roundoff, 2^-24
rng = np.random.default_rng(2026)
The kernel¶
The law as the library states it.
print(trait_run.source('render/lerp.trait'))
// lusine.render.lerp: a quantity drawn between two steps, linear interpolation: x = before + (after - before)
// alpha, alpha in [0, 1] the share of the interval from the step before to the last one that the frame is drawn at
// (0 draws the step before, 1 the last). A view's frames outnumber its steps, so a camera and the shapes it shows
// move at the frame rate. Drawing only: nothing
// interpolated goes back into the simulation. Any unit, and any length of vector.
// Source: the linear interpolation of numerical analysis (Burden and Faires, Numerical Analysis, section 3.1).
// Named lerp since `between` became a word of the language.
// Checked against a reference in numpy, which has its own tests.
namespace lusine.render;
vector<float<_>, _> lerp(vector<float<_>, _> before, vector<float<_>, _> after, float alpha)
{
return before + (after - before) * alpha;
}
The program¶
Each section runs the kernel over a lattice of items, one a case: the host sets every parameter, and main() binds the kernel over the items.
print(trait_run.law_program('render/lerp.trait', 'lerp', 3))
from lusine.kernels.render.lerp import lerp as reference
namespace probe;
lusine.version version = { { number = 1 } };
struct item {
host vector<float<_>, 3> before;
host vector<float<_>, 3> after;
host float alpha;
vector<float<_>, 3> result_0;
};
lattice<item, 1> items(spacing = 1 m, boundary = wall);
void main() { items.result_0 = lusine.render.lerp(items.before, items.after, items.alpha); }
Reference: the line between the steps¶
Random centres at two steps and random alphas: the kernel is the binary32 sequence before + (after - before) * alpha, each operation rounded once, bit for bit, and the float64 reference within a few roundings.
def lerp(before, after, alpha):
return trait_run.law('render/lerp.trait', 'lerp', vector_length=3, before=before, after=after, alpha=alpha)
b = F32(rng.normal(size=(400, 3))).astype(float)
a = F32(rng.normal(size=(400, 3))).astype(float)
t = F32(rng.random(400)).astype(float)
got = lerp(b, a, t)
model = F32(b) + F32(F32(a) - F32(b)) * F32(t)[:, None]
assert np.array_equal(F32(got), model)
ref = b + (a - b) * t[:, None]
err = np.max(np.abs(got - ref) / (np.abs(a) + np.abs(b)))
assert err <= 4 * U, err / U
print(f'{len(b)} centres bit for bit with the binary32 sequence; the float64 reference within {err / U:.2f} x 2^-24')
400 centres bit for bit with the binary32 sequence; the float64 reference within 1.73 x 2^-24
Invariant: a line, not an ease¶
Equal steps of alpha move the quantity by equal steps (to binary32's rounding): a ball rolling at a steady speed is drawn at a steady speed between two steps.
alphas = np.linspace(0.0, 1.0, 101)
x = lerp(np.zeros((101, 3)), np.tile([2.0, 0.0, -1.0], (101, 1)), alphas)
steps = np.diff(x[:, 0])
assert np.allclose(steps, 0.02, rtol=0, atol=8 * U)
fig, ax = plt.subplots(figsize=(7, 3))
ax.plot(alphas, x[:, 0], '.', ms=3)
ax.set_xlabel('alpha'); ax.set_ylabel('x (m)')
plt.show()
Limit: the two steps themselves¶
At alpha 0 the frame draws the step before exactly, and at alpha 1 the last step to within a rounding (after - before is rounded once before it is added back).
assert np.array_equal(lerp(b, a, 0.0), b)
end = lerp(b, a, 1.0)
assert np.all(np.abs(end - a) <= 2 * U * (np.abs(a) + np.abs(b)))
print('alpha 0 is the step before exactly; alpha 1 the last within a rounding')
alpha 0 is the step before exactly; alpha 1 the last within a rounding