Skip to main content

dualis_units/
lib.rs

1//! Dimensional analysis: physical quantities that refuse to be added wrongly.
2//!
3//! ```
4//! use dualis_units::{Area, Energy, Length, Mass, Power, SpecificHeat, Temperature, Time};
5//!
6//! // A unit-bearing constructor is the only place a factor of a thousand may appear.
7//! let side = Length::mm(10.0);
8//! let area: Area = side * side;                       // the dimension follows the product
9//! assert!((area.to_si() - 1e-4).abs() < 1e-18);
10//!
11//! // Absorbed power over a time is an energy, and the type says so without being told.
12//! let absorbed = Power::mw(96.0);
13//! let heat: Energy = absorbed * Time::s(1.0);
14//!
15//! // Divide it by a heat capacity and a temperature comes out.
16//! let capacity = Mass::g(2.0) * SpecificHeat::j_per_kg_k(858.0);
17//! let rise: Temperature = heat / capacity;
18//! assert!((rise.to_si() - 0.05594).abs() < 1e-4);
19//! ```
20//!
21//! And the mistake the whole crate exists to prevent does not compile:
22//!
23//! ```compile_fail
24//! use dualis_units::{Length, Time};
25//! let nonsense = Length::mm(3.0) + Time::s(1.0);
26//! ```
27//!
28//! One domain can get away with a convention. `dualis-core` began as optics and
29//! said "millimetres, nanometres and seconds, everywhere" in a doc comment, and
30//! that held because every number in the crate was a length, a wavelength or a
31//! fraction. It stops holding the moment a second domain arrives: a kelvin, a
32//! newton and a watt are all `f64`, they all add, and the compiler and the tests
33//! both stay green while the physics goes wrong.
34//!
35//! So dimension lives in the type. [`Qty`] carries the seven SI base exponents as
36//! const generic parameters, which makes `Length + Time` a compile error and
37//! `Force * Length` an [`Energy`] — and costs nothing at runtime, since a `Qty`
38//! is an `f64` and every operation on it is the `f64` operation.
39//!
40//! # Storage is always SI base units
41//!
42//! A `Qty` holds metres, kilograms, seconds, amperes, kelvin, moles, candela —
43//! never millimetres, never nanometres. Those are *entry and exit* forms:
44//!
45//! ```
46//! use dualis_units::{Length, Time, Velocity};
47//!
48//! let d = Length::mm(120.0);
49//! let t = Time::ms(4.0);
50//! let v: Velocity = d / t;
51//! assert!((v.to_si() - 30.0).abs() < 1e-12);   // 30 m/s
52//! assert!((d.in_nm() - 1.2e8).abs() < 1.0);
53//! ```
54//!
55//! That way there is exactly one representation to reason about, and the
56//! unit-bearing constructors are the only place a factor of 1000 can hide.
57//!
58//! # What this cannot do
59//!
60//! **Angles are dimensionless**, so [`Frequency`] and an angular velocity are the
61//! same type — SI says radians are m/m, and no dimensional system can separate
62//! them. Same for torque and energy. Where that distinction matters, it has to be
63//! carried by a newtype in the domain crate, not here.
64//!
65//! **Only declared products compose.** `Length * Length` is an [`Area`] because
66//! that pair is written down below. Deriving arbitrary products would need
67//! arithmetic on const generic parameters, which is unstable, so the alternative
68//! to a declared list is a dependency on `uom`. The list is cheap to extend, and
69//! anything undeclared can always go through [`Qty::from_si`].
70
71// Every public item carries a doc comment. Denied rather than warned: a public physics API
72// whose `Length::mm` shows a blank summary in rustdoc is documented in the sense that a
73// paragraph exists somewhere, and not in the sense a reader needs.
74#![deny(missing_docs)]
75#![forbid(unsafe_code)]
76
77use core::fmt;
78use core::ops::{Add, AddAssign, Div, Mul, Neg, Sub, SubAssign};
79
80use serde::{Deserialize, Deserializer, Serialize, Serializer};
81
82pub mod vector;
83pub use vector::{AccelerationVec, ForceVec, LengthVec, MomentumVec, QVec3, VelocityVec};
84
85/// A quantity, with the seven SI base dimensions in its type.
86///
87/// The parameters are the exponents of metre, kilogram, second, ampere, kelvin,
88/// mole and candela, in that order, so a velocity (m·s⁻¹) is `Qty<1,0,-1,0,0,0,0>`
89/// — which is what [`Velocity`] names.
90///
91/// Addition, subtraction, negation, comparison and scaling by a plain `f64` work
92/// for every dimension. Multiplication and division between two quantities work
93/// for the pairs declared in this module.
94#[derive(Clone, Copy, PartialEq, PartialOrd, Default)]
95pub struct Qty<
96    const L: i8,
97    const M: i8,
98    const T: i8,
99    const I: i8,
100    const K: i8,
101    const N: i8,
102    const J: i8,
103>(f64);
104
105impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
106    Qty<L, M, T, I, K, N, J>
107{
108    /// Zero, which is the one value every dimension shares.
109    pub const ZERO: Self = Qty(0.0);
110
111    /// Wrap a number already in SI base units. The escape hatch: use it when a
112    /// dimension has no name here, and name it if you use it twice.
113    ///
114    /// `const`, so a dimensioned constant can be written without a lazy static.
115    pub const fn from_si(value: f64) -> Self {
116        Qty(value)
117    }
118
119    /// The value in SI base units.
120    pub const fn to_si(self) -> f64 {
121        self.0
122    }
123
124    /// The seven exponents, for diagnostics and for a runtime dimension check at
125    /// a boundary the type system does not cross (deserialisation, FFI).
126    pub const fn dimension() -> [i8; 7] {
127        [L, M, T, I, K, N, J]
128    }
129
130    /// Magnitude without its sign, in the same dimension.
131    pub fn abs(self) -> Self {
132        Qty(self.0.abs())
133    }
134
135    /// The smaller of two quantities of the same dimension.
136    pub fn min(self, other: Self) -> Self {
137        Qty(self.0.min(other.0))
138    }
139
140    /// The larger of two quantities of the same dimension.
141    pub fn max(self, other: Self) -> Self {
142        Qty(self.0.max(other.0))
143    }
144
145    /// Whether the magnitude is neither infinite nor NaN.
146    ///
147    /// Worth checking where a limit is reported rather than computed: several methods here
148    /// return an infinity to mean "no limit", which is honest but arithmetic on it is not.
149    pub fn is_finite(self) -> bool {
150        self.0.is_finite()
151    }
152
153    /// Sign of the magnitude, as a plain number — a sign has no dimension.
154    pub fn signum(self) -> f64 {
155        self.0.signum()
156    }
157
158    /// Linear interpolation, which stays within the dimension.
159    pub fn lerp(self, other: Self, t: f64) -> Self {
160        Qty(self.0 + (other.0 - self.0) * t)
161    }
162}
163
164// ---------------------------------------------------------------------------
165// Dimension-preserving arithmetic: works for every dimension at once, because
166// none of it changes the exponents.
167// ---------------------------------------------------------------------------
168
169macro_rules! generic_op {
170    ($trait:ident, $method:ident, $op:tt) => {
171        impl<
172                const L: i8,
173                const M: i8,
174                const T: i8,
175                const I: i8,
176                const K: i8,
177                const N: i8,
178                const J: i8,
179            > $trait for Qty<L, M, T, I, K, N, J>
180        {
181            type Output = Self;
182            fn $method(self, rhs: Self) -> Self {
183                Qty(self.0 $op rhs.0)
184            }
185        }
186    };
187}
188
189generic_op!(Add, add, +);
190generic_op!(Sub, sub, -);
191
192impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
193    AddAssign for Qty<L, M, T, I, K, N, J>
194{
195    fn add_assign(&mut self, rhs: Self) {
196        self.0 += rhs.0;
197    }
198}
199
200impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
201    SubAssign for Qty<L, M, T, I, K, N, J>
202{
203    fn sub_assign(&mut self, rhs: Self) {
204        self.0 -= rhs.0;
205    }
206}
207
208impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8> Neg
209    for Qty<L, M, T, I, K, N, J>
210{
211    type Output = Self;
212    fn neg(self) -> Self {
213        Qty(-self.0)
214    }
215}
216
217impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
218    Mul<f64> for Qty<L, M, T, I, K, N, J>
219{
220    type Output = Self;
221    fn mul(self, k: f64) -> Self {
222        Qty(self.0 * k)
223    }
224}
225
226impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
227    Div<f64> for Qty<L, M, T, I, K, N, J>
228{
229    type Output = Self;
230    fn div(self, k: f64) -> Self {
231        Qty(self.0 / k)
232    }
233}
234
235impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
236    Mul<Qty<L, M, T, I, K, N, J>> for f64
237{
238    type Output = Qty<L, M, T, I, K, N, J>;
239    fn mul(self, q: Qty<L, M, T, I, K, N, J>) -> Qty<L, M, T, I, K, N, J> {
240        Qty(self * q.0)
241    }
242}
243
244/// Dividing two quantities of the *same* dimension gives a plain number — which
245/// is the one product rule that needs no exponent arithmetic, and the one every
246/// tolerance check uses.
247impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8> Div
248    for Qty<L, M, T, I, K, N, J>
249{
250    type Output = f64;
251    fn div(self, rhs: Self) -> f64 {
252        self.0 / rhs.0
253    }
254}
255
256impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
257    fmt::Debug for Qty<L, M, T, I, K, N, J>
258{
259    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
260        write!(f, "{}", self.0)?;
261        for (symbol, exponent) in [
262            ("m", L),
263            ("kg", M),
264            ("s", T),
265            ("A", I),
266            ("K", K),
267            ("mol", N),
268            ("cd", J),
269        ] {
270            match exponent {
271                0 => {}
272                1 => write!(f, "·{symbol}")?,
273                e => write!(f, "·{symbol}^{e}")?,
274            }
275        }
276        Ok(())
277    }
278}
279
280impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
281    fmt::Display for Qty<L, M, T, I, K, N, J>
282{
283    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
284        fmt::Debug::fmt(self, f)
285    }
286}
287
288// Serialised as the bare SI number: a scene file stays readable, and the
289// dimension is carried by the field's type rather than repeated in the data.
290impl<const L: i8, const M: i8, const T: i8, const I: i8, const K: i8, const N: i8, const J: i8>
291    Serialize for Qty<L, M, T, I, K, N, J>
292{
293    fn serialize<S: Serializer>(&self, s: S) -> Result<S::Ok, S::Error> {
294        self.0.serialize(s)
295    }
296}
297
298impl<
299        'de,
300        const L: i8,
301        const M: i8,
302        const T: i8,
303        const I: i8,
304        const K: i8,
305        const N: i8,
306        const J: i8,
307    > Deserialize<'de> for Qty<L, M, T, I, K, N, J>
308{
309    fn deserialize<D: Deserializer<'de>>(d: D) -> Result<Self, D::Error> {
310        f64::deserialize(d).map(Qty)
311    }
312}
313
314// ---------------------------------------------------------------------------
315// The dimensions themselves.
316// ---------------------------------------------------------------------------
317
318/// A pure ratio: reflectance, duty cycle, refractive index, Strehl.
319pub type Dimensionless = Qty<0, 0, 0, 0, 0, 0, 0>;
320
321/// Metres.
322pub type Length = Qty<1, 0, 0, 0, 0, 0, 0>;
323/// Kilograms.
324pub type Mass = Qty<0, 1, 0, 0, 0, 0, 0>;
325/// Seconds.
326pub type Time = Qty<0, 0, 1, 0, 0, 0, 0>;
327/// Amperes.
328pub type Current = Qty<0, 0, 0, 1, 0, 0, 0>;
329/// Absolute temperature. Kelvin only — see [`Temperature::celsius`].
330pub type Temperature = Qty<0, 0, 0, 0, 1, 0, 0>;
331/// Moles.
332pub type Amount = Qty<0, 0, 0, 0, 0, 1, 0>;
333/// Candelas.
334pub type LuminousIntensity = Qty<0, 0, 0, 0, 0, 0, 1>;
335
336/// Square metres.
337pub type Area = Qty<2, 0, 0, 0, 0, 0, 0>;
338/// Cubic metres.
339pub type Volume = Qty<3, 0, 0, 0, 0, 0, 0>;
340/// Metres per second.
341pub type Velocity = Qty<1, 0, -1, 0, 0, 0, 0>;
342/// Metres per second squared.
343pub type Acceleration = Qty<1, 0, -2, 0, 0, 0, 0>;
344/// kg·m·s⁻¹ — mass times velocity, and the thing a closed system conserves
345/// exactly rather than nearly.
346pub type Momentum = Qty<1, 1, -1, 0, 0, 0, 0>;
347/// Newtons.
348pub type Force = Qty<1, 1, -2, 0, 0, 0, 0>;
349/// Pascals. Also the unit of an energy density and of a stress, which are the same
350/// dimension and not a coincidence.
351pub type Pressure = Qty<-1, 1, -2, 0, 0, 0, 0>;
352/// Joules.
353pub type Energy = Qty<2, 1, -2, 0, 0, 0, 0>;
354/// Watts.
355pub type Power = Qty<2, 1, -3, 0, 0, 0, 0>;
356/// kg·m⁻³. Note that a glass catalogue quotes g/cm³, a factor of a thousand away —
357/// see [`Density::g_per_cm3`].
358pub type Density = Qty<-3, 1, 0, 0, 0, 0, 0>;
359/// Cycles per second. Dimensionally identical to an angular velocity, since a
360/// radian is m/m — the type system cannot and should not pretend otherwise.
361pub type Frequency = Qty<0, 0, -1, 0, 0, 0, 0>;
362/// Power per unit area, W·m⁻². What a detector face actually receives.
363pub type Irradiance = Qty<0, 1, -3, 0, 0, 0, 0>;
364/// W·m⁻¹·K⁻¹ — the `k` of Fourier's law.
365pub type ThermalConductivity = Qty<1, 1, -3, 0, -1, 0, 0>;
366/// J·kg⁻¹·K⁻¹ — the `c_p` that says how much heat a gram of glass can hide.
367pub type SpecificHeat = Qty<2, 0, -2, 0, -1, 0, 0>;
368/// m²·s⁻¹ — thermal diffusivity `α = k/(ρ c_p)`, and also mass diffusivity.
369pub type Diffusivity = Qty<2, 0, -1, 0, 0, 0, 0>;
370/// K⁻¹ — the coefficient that turns absorbed light into a focus shift.
371pub type ThermalExpansion = Qty<0, 0, 0, 0, -1, 0, 0>;
372/// kg·m² — how hard a body is to spin up about an axis.
373///
374/// The rotational counterpart of mass, and unlike mass it depends on the axis: a
375/// pencil is trivial to spin about its length and awkward about its middle. That
376/// direction-dependence is why it is a tensor and why a free body's rotation is
377/// interesting rather than uniform.
378pub type MomentOfInertia = Qty<2, 1, 0, 0, 0, 0, 0>;
379/// kg·m²·s⁻¹ — the rotational counterpart of momentum, and conserved for the same
380/// reason.
381pub type AngularMomentum = Qty<2, 1, -1, 0, 0, 0, 0>;
382/// N·m⁻¹ — a spring's `k`, and the penalty stiffness a contact is modelled with.
383///
384/// This is what sets a mechanical solver's stability limit: a mass on a spring
385/// oscillates with period `2π√(m/k)`, and an explicit integrator has to resolve that
386/// period whether or not anyone cares about it. Stiff contact is expensive for
387/// exactly this reason.
388pub type Stiffness = Qty<0, 1, -2, 0, 0, 0, 0>;
389/// N·s·m⁻¹ — a dashpot's `c`. Force proportional to velocity, and the only place a
390/// mechanical simulation loses energy on purpose.
391pub type Damping = Qty<0, 1, -1, 0, 0, 0, 0>;
392/// Coulombs.
393pub type Charge = Qty<0, 0, 1, 1, 0, 0, 0>;
394/// Volts.
395pub type Voltage = Qty<2, 1, -3, -1, 0, 0, 0>;
396/// J·K⁻¹ — mass times specific heat. How much heat a thing can hide before it
397/// shows up as a temperature.
398pub type HeatCapacity = Qty<2, 1, -2, 0, -1, 0, 0>;
399
400// ---------------------------------------------------------------------------
401// Declared products. Each line also gives the two divisions that undo it.
402// ---------------------------------------------------------------------------
403
404macro_rules! product {
405    ($a:ty, $b:ty => $c:ty) => {
406        impl Mul<$b> for $a {
407            type Output = $c;
408            fn mul(self, rhs: $b) -> $c {
409                Qty(self.0 * rhs.0)
410            }
411        }
412        impl Mul<$a> for $b {
413            type Output = $c;
414            fn mul(self, rhs: $a) -> $c {
415                Qty(self.0 * rhs.0)
416            }
417        }
418        impl Div<$b> for $c {
419            type Output = $a;
420            fn div(self, rhs: $b) -> $a {
421                Qty(self.0 / rhs.0)
422            }
423        }
424        impl Div<$a> for $c {
425            type Output = $b;
426            fn div(self, rhs: $a) -> $b {
427                Qty(self.0 / rhs.0)
428            }
429        }
430    };
431}
432
433macro_rules! square {
434    ($a:ty => $c:ty) => {
435        impl Mul<$a> for $a {
436            type Output = $c;
437            fn mul(self, rhs: $a) -> $c {
438                Qty(self.0 * rhs.0)
439            }
440        }
441        impl Div<$a> for $c {
442            type Output = $a;
443            fn div(self, rhs: $a) -> $a {
444                Qty(self.0 / rhs.0)
445            }
446        }
447    };
448}
449
450square!(Length => Area);
451product!(Area, Length => Volume);
452product!(Velocity, Time => Length);
453product!(Acceleration, Time => Velocity);
454product!(Mass, Acceleration => Force);
455product!(Mass, Velocity => Momentum);
456product!(Force, Length => Energy);
457product!(Force, Time => Momentum);
458product!(Pressure, Area => Force);
459product!(Power, Time => Energy);
460product!(Irradiance, Area => Power);
461product!(Density, Volume => Mass);
462product!(Current, Time => Charge);
463product!(Voltage, Current => Power);
464product!(Mass, Area => MomentOfInertia);
465product!(MomentOfInertia, Frequency => AngularMomentum);
466product!(Stiffness, Length => Force);
467product!(Damping, Velocity => Force);
468product!(Mass, SpecificHeat => HeatCapacity);
469product!(HeatCapacity, Temperature => Energy);
470product!(Frequency, Time => Dimensionless);
471
472impl Area {
473    /// The side of a square of this area. The one root worth naming, because it
474    /// is how a beam radius comes back out of a spot area.
475    pub fn sqrt(self) -> Length {
476        Qty(self.0.sqrt())
477    }
478}
479
480// ---------------------------------------------------------------------------
481// Unit-bearing entry and exit. The only place a factor of 1000 may appear.
482// ---------------------------------------------------------------------------
483
484impl Length {
485    /// Metres.
486    pub fn m(v: f64) -> Length {
487        Qty(v)
488    }
489    /// Millimetres.
490    pub fn mm(v: f64) -> Length {
491        Qty(v * 1e-3)
492    }
493    /// Micrometres.
494    pub fn um(v: f64) -> Length {
495        Qty(v * 1e-6)
496    }
497    /// Nanometres. The wavelength unit, and why every `Spectrum` field is named `_nm`.
498    pub fn nm(v: f64) -> Length {
499        Qty(v * 1e-9)
500    }
501    /// As millimetres.
502    pub fn in_mm(self) -> f64 {
503        self.0 * 1e3
504    }
505    /// As micrometres.
506    pub fn in_um(self) -> f64 {
507        self.0 * 1e6
508    }
509    /// As nanometres.
510    pub fn in_nm(self) -> f64 {
511        self.0 * 1e9
512    }
513}
514
515impl Time {
516    /// Seconds.
517    pub fn s(v: f64) -> Time {
518        Qty(v)
519    }
520    /// Milliseconds.
521    pub fn ms(v: f64) -> Time {
522        Qty(v * 1e-3)
523    }
524    /// Microseconds.
525    pub fn us(v: f64) -> Time {
526        Qty(v * 1e-6)
527    }
528    /// Nanoseconds.
529    pub fn ns(v: f64) -> Time {
530        Qty(v * 1e-9)
531    }
532    /// As milliseconds.
533    pub fn in_ms(self) -> f64 {
534        self.0 * 1e3
535    }
536    /// As microseconds.
537    pub fn in_us(self) -> f64 {
538        self.0 * 1e6
539    }
540}
541
542impl Temperature {
543    /// Kelvin, which is what is stored.
544    pub fn kelvin(v: f64) -> Temperature {
545        Qty(v)
546    }
547    /// Celsius is an *offset* scale, not a scaled one, which is why it gets a
548    /// named constructor rather than a factor: 20 °C is 293.15 K, and a
549    /// temperature *difference* of 20 K is a different thing entirely.
550    pub fn celsius(v: f64) -> Temperature {
551        Qty(v + 273.15)
552    }
553    /// As degrees Celsius. Subtracts the offset; see [`Temperature::celsius`].
554    pub fn in_celsius(self) -> f64 {
555        self.0 - 273.15
556    }
557}
558
559impl Mass {
560    /// Kilograms.
561    pub fn kg(v: f64) -> Mass {
562        Qty(v)
563    }
564    /// Grams.
565    pub fn g(v: f64) -> Mass {
566        Qty(v * 1e-3)
567    }
568}
569
570impl Density {
571    /// The way a glass catalogue quotes it: N-BK7 is 2.51 g/cm³.
572    pub fn g_per_cm3(v: f64) -> Density {
573        Qty(v * 1e3)
574    }
575    /// Kilograms per cubic metre, which is what is stored.
576    pub fn kg_per_m3(v: f64) -> Density {
577        Qty(v)
578    }
579}
580
581impl Power {
582    /// Watts.
583    pub fn w(v: f64) -> Power {
584        Qty(v)
585    }
586    /// Milliwatts.
587    pub fn mw(v: f64) -> Power {
588        Qty(v * 1e-3)
589    }
590    /// Microwatts.
591    pub fn uw(v: f64) -> Power {
592        Qty(v * 1e-6)
593    }
594    /// As milliwatts.
595    pub fn in_mw(self) -> f64 {
596        self.0 * 1e3
597    }
598}
599
600impl Energy {
601    /// Joules.
602    pub fn j(v: f64) -> Energy {
603        Qty(v)
604    }
605    /// Millijoules.
606    pub fn mj(v: f64) -> Energy {
607        Qty(v * 1e-3)
608    }
609}
610
611impl Frequency {
612    /// Hertz.
613    pub fn hz(v: f64) -> Frequency {
614        Qty(v)
615    }
616    /// Kilohertz.
617    pub fn khz(v: f64) -> Frequency {
618        Qty(v * 1e3)
619    }
620    /// Megahertz.
621    pub fn mhz(v: f64) -> Frequency {
622        Qty(v * 1e6)
623    }
624    /// Period: one over the frequency. Named because `1.0 / f` cannot typecheck.
625    pub fn period(self) -> Time {
626        Qty(1.0 / self.0)
627    }
628}
629
630impl Velocity {
631    /// Metres per second.
632    pub fn m_per_s(v: f64) -> Velocity {
633        Qty(v)
634    }
635    /// Millimetres per second.
636    pub fn mm_per_s(v: f64) -> Velocity {
637        Qty(v * 1e-3)
638    }
639}
640
641impl Irradiance {
642    /// Watts per square metre, which is what is stored.
643    pub fn w_per_m2(v: f64) -> Irradiance {
644        Qty(v)
645    }
646    /// How an illumination spec is usually written: mW/cm².
647    pub fn mw_per_cm2(v: f64) -> Irradiance {
648        Qty(v * 10.0)
649    }
650}
651
652impl ThermalConductivity {
653    /// W·m⁻¹·K⁻¹, the unit a materials table uses.
654    pub fn w_per_m_k(v: f64) -> ThermalConductivity {
655        Qty(v)
656    }
657}
658
659impl SpecificHeat {
660    /// J·kg⁻¹·K⁻¹, the unit a materials table uses.
661    pub fn j_per_kg_k(v: f64) -> SpecificHeat {
662        Qty(v)
663    }
664}
665
666impl ThermalExpansion {
667    /// Catalogues quote it in parts per million per kelvin: N-BK7 is 7.1.
668    pub fn ppm_per_k(v: f64) -> ThermalExpansion {
669        Qty(v * 1e-6)
670    }
671}
672
673impl Dimensionless {
674    /// A bare ratio, for the one case where a number genuinely has no dimension:
675    /// a reflectance, a duty cycle, a refractive index.
676    pub fn ratio(v: f64) -> Dimensionless {
677        Qty(v)
678    }
679}
680
681// ---------------------------------------------------------------------------
682// Physical constants, in SI base units, so that a formula written with them
683// carries its own dimensional proof.
684// ---------------------------------------------------------------------------
685
686/// Speed of light in vacuum, m·s⁻¹ (exact by definition).
687pub const C: Velocity = Qty(299_792_458.0);
688/// Planck constant, J·s (exact by definition).
689pub const PLANCK: Qty<2, 1, -1, 0, 0, 0, 0> = Qty(6.626_070_15e-34);
690/// Boltzmann constant, J·K⁻¹ (exact by definition).
691pub const BOLTZMANN: HeatCapacity = Qty(1.380_649e-23);
692/// Stefan-Boltzmann constant, W·m⁻²·K⁻⁴ — radiative exchange lives on this.
693pub const STEFAN_BOLTZMANN: Qty<0, 1, -3, 0, -4, 0, 0> = Qty(5.670_374_419e-8);
694/// Standard gravity, m·s⁻².
695pub const G0: Acceleration = Qty(9.806_65);
696
697/// Energy of one photon at a vacuum wavelength: `E = hc/λ`.
698///
699/// The bridge between a spectrum and a photon count, and the reason a detector's
700/// response is not the same shape as a lamp's output.
701pub fn photon_energy(wavelength: Length) -> Energy {
702    Qty(PLANCK.0 * C.0 / wavelength.0)
703}
704
705#[cfg(test)]
706mod tests {
707    use super::*;
708
709    /// The point of the crate: a product of dimensions lands on the type that
710    /// names it, whichever route it took there.
711    #[test]
712    fn products_land_on_the_named_dimension() {
713        let m = Mass::kg(2.0);
714        let a = Acceleration::from_si(3.0);
715        let f: Force = m * a;
716        assert!((f.to_si() - 6.0).abs() < 1e-12);
717
718        // Two different routes to the same energy, and they unify.
719        let by_work: Energy = f * Length::m(4.0);
720        let by_power: Energy = Power::w(24.0) * Time::s(1.0);
721        assert!((by_work - by_power).abs().to_si() < 1e-12);
722
723        // And multiplication commutes, as it must.
724        let swapped: Force = a * m;
725        assert_eq!(f, swapped);
726    }
727
728    /// Millimetres and nanometres are entry forms only; storage is metres. A
729    /// wavelength and a lens diameter therefore compare correctly without anyone
730    /// remembering which convention each was written in.
731    #[test]
732    fn unit_prefixes_are_only_a_doorway() {
733        let lens = Length::mm(25.4);
734        let green = Length::nm(550.0);
735        assert!(lens > green);
736        assert!((lens.to_si() - 0.0254).abs() < 1e-15);
737        assert!((green.in_nm() - 550.0).abs() < 1e-9);
738        // 25.4 mm is 46181.8... wavelengths of green light.
739        let waves = lens / green;
740        assert!((waves - 46_181.8).abs() < 0.1, "got {waves}");
741    }
742
743    /// Dividing like by like gives a plain number, which is what every
744    /// tolerance and every reflectance is.
745    #[test]
746    fn like_over_like_is_a_bare_number() {
747        let reflected = Power::mw(0.42);
748        let incident = Power::mw(10.0);
749        let r: f64 = reflected / incident;
750        assert!((r - 0.042).abs() < 1e-12);
751    }
752
753    /// Celsius offsets rather than scales, and getting that wrong is a 273 K
754    /// error that no dimensional check would ever catch.
755    #[test]
756    fn celsius_is_an_offset_not_a_factor() {
757        assert!((Temperature::celsius(20.0).to_si() - 293.15).abs() < 1e-12);
758        assert!((Temperature::kelvin(293.15).in_celsius() - 20.0).abs() < 1e-12);
759        // A *difference* of 20 K is not 293.15 K, and only one of these is a
760        // temperature you can put in Stefan-Boltzmann.
761        let rise = Temperature::kelvin(313.15) - Temperature::kelvin(293.15);
762        assert!((rise.to_si() - 20.0).abs() < 1e-12);
763    }
764
765    /// The optics-to-thermal chain this whole crate exists to make safe: a
766    /// surface absorbs a fraction of an irradiance over an area, and the watts
767    /// that result heat a mass with a known specific heat.
768    #[test]
769    fn absorbed_light_becomes_a_temperature_rise() {
770        let irradiance = Irradiance::mw_per_cm2(50.0); // 500 W/m²
771        let area: Area = Length::mm(10.0) * Length::mm(10.0); // 1e-4 m²
772        let absorptance = 0.02; // what SurfaceOptics::absorptance returns
773        let absorbed: Power = irradiance * area * absorptance;
774        assert!((absorbed.to_si() - 0.001).abs() < 1e-12, "{absorbed:?}");
775
776        // 1 mW into a 2 g piece of glass for 1 s.
777        let glass = Mass::g(2.0);
778        let c_p = SpecificHeat::j_per_kg_k(858.0); // N-BK7
779        let capacity: HeatCapacity = glass * c_p;
780        let heat: Energy = absorbed * Time::s(1.0);
781        let rise: Temperature = heat / capacity;
782        assert!(
783            (rise.to_si() - 0.000_582_7).abs() < 1e-7,
784            "expected about 0.58 mK, got {rise:?}"
785        );
786    }
787
788    /// A photon at 550 nm carries 3.6e-19 J, and the count per watt follows.
789    /// This is the number that separates radiometry from photon counting.
790    #[test]
791    fn photon_energy_matches_the_textbook_figure() {
792        let e = photon_energy(Length::nm(550.0));
793        assert!((e.to_si() - 3.612e-19).abs() < 1e-21, "{e:?}");
794        // 2.26 eV, and about 2.77e18 photons in a joule.
795        let per_joule = Energy::j(1.0) / e;
796        assert!((per_joule - 2.768e18).abs() < 1e15, "got {per_joule:e}");
797    }
798
799    /// Debug prints the dimension, so a mismatch found at a boundary can be
800    /// reported in a form a human recognises.
801    #[test]
802    fn debug_shows_the_dimension() {
803        assert_eq!(format!("{:?}", Force::from_si(6.0)), "6·m·kg·s^-2");
804        assert_eq!(format!("{:?}", Dimensionless::ratio(0.5)), "0.5");
805        assert_eq!(Force::dimension(), [1, 1, -2, 0, 0, 0, 0]);
806    }
807
808    /// Serialised as the bare SI number: the dimension is in the field's type,
809    /// not repeated in every scene file.
810    #[test]
811    fn serialises_as_a_bare_si_number() {
812        let json = serde_json::to_string(&Length::mm(25.4)).unwrap();
813        assert_eq!(json, "0.0254");
814        let back: Length = serde_json::from_str(&json).unwrap();
815        assert_eq!(back, Length::mm(25.4));
816    }
817}