Boost Units Library

• Download source - 16.4 KB

Introduction

Some time ago, I began to write an engine tuning application. Development started smoothly, however as time went on, a number of annoying bugs occurred - air pressure in millibars passed into functions that expected the air pressure in pascals and the like. The inputs into the functions were not always fully documented (well ... the code was for my own use) and sometimes the documentation was wrong. It occurred to me that it would be much better if I could use the compiler to catch what are conceptually type-mismatch errors.

The boost units library seemed like just the thing to solve my problems.

Dimensional Analysis

Dimensional Analysis makes the assumption that quantities have an underlying dimensionality. It is taught in Engineering and Scientific classes at University as a way of checking the correctness of formulas. Length, time, mass, temperature are "base" dimensions. Other dimensions can be assembled from them. For example, speed is a derived dimension measured in meters/sec. Energy is also a derived dimension and measured in kilogram*meters²/second².

The choice of base dimensions is somewhat arbitrary. The boost units library defines the following base dimensions and units for the SI (metric) system:

• Electric current (Amperes)
• Light intensity (Candela)
• Temperature (Kelvin)
• Mass (kilograms)
• Time (seconds)
• Length (meters)
• Number of atoms/molecules (moles)

The SI system is not the only possible system. For example, it is possible to define a system that uses fewer dimensions and different units such as feet, seconds, pounds and degrees Fahrenheit.

Formulas must be dimensionally consistent - it does not make sense to add a mile to a kilogram. The boost units library defines addition, subtraction, multiplication and division for units with the same dimensionality. The following code works fine:

quantity<mass> m = 10.0 * kilograms;                          // dim(mass) = M
quantity<acceleration> a = 9.8 * meter_per_second_squared;    // dim(acceleration) = L.T^-2

quantity<force> f = m * a;                                    // dim(force) = M.L.T^-2

The following code is not dimensionally correct and will not compile.

quantity<mass> m = 10.0 * kilograms;
quantity<acceleration> a = 9.8 * meter_per_second_squared;

quantity<force> f = m * a * a;    // <== error! Won't compile. dim(LHS) does not equal dim(RHS)

So how does all this magic occur? The answer is the boost meta-programming library (mpl). Fortunately, like all good libraries, it is not necessary to understand all the details of its implementation to use it. However, like all good libraries, if you really want to use it to the full, you really need to get familiar with its innards.

The SI system defines a wide range of derived dimensions and associated units which can be accessed by simply including the si.hpp.

#include <boost/units/systems/si.hpp>

...

void example_declarations()
{
    quantity<length> length = 2.5 * meters;

    quantity<velocity> velocity = 2.5 * meters_per_second;

    quantity<mass> m = 2.5 * kilogram;

    quantity<si::time> t(2.5 * si::seconds);  // name clash with std::time, std::second.

    quantity<pressure> pressure(101300.0 * pascals);

    quantity<mass_density> density(50.0 * kilogrammes_per_cubic_metre);

    quantity<temperature> temp(285.15 * kelvin);
}

The stronger typing provided by units can also catch bugs like the one below.

double density(double T, double P, double humidity);  // prototype

double T = 20.0; // celsius
double P = 101300.0; // pascals

double air_density = density(P, T, humidity);   // <== T and P are swapped around!

Using units ...

quantity<mass_density> density(quantity<temperature> T, quantity<pressure> P, double humidity);  // prototype

quantity<temperature> T = 292.15 * kelvin;
quantity<pressure> P = 101300.0 * pascals;

quantity<mass_density> air_density = density(P, T, humidity);     // error! Wont compile. T and P are swapped around!

New Types

It is usually unnecessary to go outside the dimensions/units provided, however it is relatively straight-forward to do so. The code contains a definition of a type to support Avogadros numbers. The standard process is:

1. define the derived_dimension,
2. defined associated unit, and
3. define the type.

typedef derived_dimension<length_base_dimension, 2,
                         time_base_dimension, -2,
                         temperature_base_dimension, -1>::type Avogadros_constant_dim;

typedef unit<Avogadros_constant_dim, si::system> Avogadros_constant_unit;

typedef quantity<Avogadros_constant_unit> Avogadros_constant_type;

Alternatively, a type can be defined from an instance:

const auto Rd = 287.05 * joules/(kilogram * kelvin);    // dry air
const auto Rv = 461.495 * joules/(kilogram * kelvin);   // water vapour

typedef decltype(Rd) Avogadros_constant_type;

New Units

It may on occassions be convenient to use non-SI units, which will nearly always differ from SI units by a constant factor. The library provides a simple template for defining scaled units.

// Definition of millibars (1 millibar = 100 pascals)
typedef make_scaled_unit<si::pressure, scale<10, static_rational<2> > >::type millibars;
const auto mbars = 100 * pascals;

quantity<millibars> mb(10 * mbars);
quantity<pressure> P(mb);

std::cout << "P = " << P << "\n";
std::cout << "P = " << mb << "\n";

assert(is_equal(P, mb));  // passes OK.

The code above produces the following output:

1000 m^-1 kg s^-2
10 h(m^-1 kg s^-2)

'h' indicates a heterogeneous system. Scaled units can be explicitly converted to base units, which makes it easy and safe to pass quantities to functions that expect quantities in different units.

quantity<mass_density> density(quantity<temperature> T,
                                     quantity<pressure> P,      // pascals
                                     double humidity);  // prototype

quantity<temperature> T = 292.15 * kelvin;
quantity<millibars> P = 1013.0 * mbars;

quantity<mass_density> density1 = density(T, P, humidity);  // error! Wont compile!
quantity<mass_density> density2 = density(T, quantity<pressure>(P), humidity); // converted.

The Mars Climate Orbiter might have benefited greatly from this sort of facility.

Relative vs Absolute Units

For most quantities, it is not necessary to make any distinction between relative and absolute measurements. The "zero" value for units is usually zero. For example, if there is no mass present, then the mass = 0 kilograms = 0 grams = 0 pounds. The one common exception is temperature. For historic reasons, the "zero" value for the Fahrenheit and Celsius scales is not zero. This results in ambiguity. If a temperature of 0 Kelvin is passed into a formula, does that mean a temperature difference of 0 Kelvin = 0 Fahrenheit = 0 Celsius? Or does it mean an absolute temperature of 0 Kelvin = -459.67 Fahrenheit = -273.15 Celsius?

The library handles this by defining 2 separate types as shown below.

typedef absolute<celsius::temperature> abs_celsius;   // saves typing...

quantity<celsius::temperature> boiling_pt = 100 * celsius::degrees; // relative temperature  - the default.
quantity<temperature> T_kelvin(boiling_pt);                         // convert to Kelvin
std::cout << T_kelvin << "\n";

quantity<abs_celsius> boiling_pt2(100 *abs_celsius());  // absolute temperature.
quantity<temperature> T2_kelvin(boiling_pt);            // convert to Kelvin
std::cout << T2_kelvin << "\n";

The code above produces the following output:

100 K
373.15 Absolute K

The Sample Code

The sample code consists of the implementation of a number of simple functions using both double(s) and quantity<double>. Test code is used to determine the time taken to execute the code in both cases.

The sample was built with NetBeans on Mint 17. A slightly different version of the code was tested using VS 2013 on Windows 8.

Caveats

Zero Overhead (Almost)

One of the great attributes of the boost unit library is that it is supposed to have zero overhead. The sample code bears that out ... almost. The one exception is the dew_point__Newton_Raphson function. The MSVC compiler is particularly bad: the boost units version runs 100+ time slower than the code which uses double values.

Why? It is not particularly clear. It could be that optimiser is not up to the job. It could be that the compiler added expensive object-unwinding code - a similar slow-down can be achieved by adding code to throw an exception. The effect is very compiler specific. Either way, the claim of zero overhead is not entirely true.

Usage

If you inspect the code for saturation_pressure, you may notice that it calls the polynomial function p takes a double. Not everything can be made dimensionally safe: polynomial functions are by definition of the form X + X² + X³ + ... and so do not have a defined dimensionality. It is probably a good objective to use the units library to make all public functions/methods as safe as possible.

Conclusion

The boost units library is an awesome demonstration of the power of templates and strong typing. Hopefully, this article will raise awareness and encourage others to use it.