Generic Numeric Computations in Rust

In this post, I will explain how to implement numeric computations in Rust that are generic over the actual numeric type used (e.g., f32/f64). I will do this by explaining step by step, which issues arise when going from a concrete implementation to a generic one, providing a sample implementation of the required traits. Step by step I will add more of the required features and explain why they are needed. I will only cover floating-point arithmetic and the basic floating-point data types f32 and f64. Integer maths and more complex (non-copy) data types will be covered in another post, if I get around to writing one.

Prerequisites

This article assumes you have basic familiarity with Rust, its numeric types, and how to define and use Traits to create generic methods. It helps to be familiar with the Copy and other basic traits.

The Problem

We want to implement a moderately complex mathematical computation, for example

fn double(value: f64) -> f64 {
	value + value 
}

Nothing fancy: we work with numbers (f64) and apply a basic mathematical operator (+).

But wait. Why are we using f64? I might want to use this for both f32 or f64 and I really don't want to have two functions for that. They would be identical except for the type! Sounds like we need to generalize. Maybe we can do something like this:

fn double<T: Number>(value: T) -> T {
	value + value 
}

Number here is a trait that we still need to define. But what do we require from Number? It has to support addition with another Number and produce another Number. Easy enough. Rust has a Trait that covers the requirements for doing addition. So our Number trait now looks like this:

trait Number: Add<Self, Output=Self> {}

We need two more things from our trait for this to work, though.

Since we use value twice, we have to ensure it is Copy. We could also relax that to Clone but that would mean we manually need to write the calls to clone() where needed, or restrict us to using references. In our current case both f32 and f64 are copy, so let's not worry about that, for now.


Aside: This is fine as long as we can assume copy but if we want to generalize over array like things, we need to put in the effort of thinking about memory management, where to clone, where to create new instances. Especially for large arrays this is very relevant for performance.


Secondly, we need to ensure that the type we use is Sized. That means the size of the values we operate on has to be known at compile time. Which makes sense: when we add two numbers, we want to store the result. That means there has to be a place for it to be stored, which has to have a certain size.

Okay, so with that cleared up we now have

trait Number: Sized + Copy + Add<Self, Output=Self> {}

And now this is the full working example, including the implementations for f64 and f32, which are pretty easy because we don't have to do anything except for stating that they fulfil the requirements of our trait.

A Basic Implementation

use core::ops::Add;

trait Number: Sized + Copy + Add<Self, Output=Self> {}

impl Number for f64{}
impl Number for f32{}

fn double<T: Number>(value: T) -> T {
	value + value 
}

pub fn main(){
    double(4.09_f32);
    double(5.01_f64);
}

Adding More Arithmetic Operations

To support more operations, we can just add the traits for division, multiplication, subtraction. That gives us all of basic arithmetic.

This is our trait by now

trait Number: Sized + Copy + Add<Self, Output=Self> + Mul<Self, Output=Self> + Sub<Self, Output=Self> + Div<Self, Output=Self>{}

So we are ready to do some real maths! Why don't we calculate the circumference of a circle based on its radius?

fn circle_circumference<T: Number>(r: T) -> T {
    2.0 * r * PI
}

No wait. We can't use a float constant. That's either a f32 or f64, but not a Number. So how do we get Number constants?

Working with Constants

Unfortunately, there is no easy way to define literals for our Number and thus we need to convert from the existing numeric types, e.g., by requiring Number to implement From<T> for something we can create literals for. In this case it makes the most sense to use f64 here as we get the best precision, but this does mean, that our computation might lose out on some of that precision if we downcast to f32.

There is another problem there. Where is PI coming from? That should be Number::PI. And that means we need to define it as part of our trait.

use std::f64;

trait Number: From<f64> + /* other traits */ {
    const PI: Self;
}

and then we need to provide the value when we implement the trait

impl Number for f64 {
    const PI: Self = std::f64::consts::PI;
}

impl Number for f32 {
    const PI: Self = std::f32::consts::PI;
}

Grand. That looks easy enough. And we can do the same things for every constant we need. We can even add a few quality of life constants like Zero and One that allow us to avoid writing T::from(1.0) or T::from(0.0) all the time.

[!Aside] Doing this for a lot of constants is a bit repetitive with a fair amount of boilerplate. We could define a macro to simplify this:

macro_rules! impl_number_consts {
    ($ty:ty, $($const_name:ident),*) => {
        $(const $const_name: Self = std::$ty::consts::$const_name;)*
    };
}

impl Number for f64 {
    impl_number_consts!(f64, PI, E, SQRT_2);
}

This generates the constant definitions automatically for each type.

The working example then looks like this

use std::ops::*;

trait Number:
    Sized +
    Copy +
    Add<Self, Output=Self> +
    Sub<Self, Output=Self> +
    Mul<Self, Output=Self> +
    Div<Self, Output=Self> +
    From<f64> {

    const PI: Self;

}

impl Number for f64 {
    const PI: Self = std::f64::consts::PI;
}

impl Number for f32 {
    const PI: Self = std::f32::consts::PI;
}

<!-- SPELLING: "circumference" should be "circumference" -->
fn circle_circumference<T: Number>(r: T) -> T {
    T::from(2.0) * r * T::PI
}

Adding Comparisons

This leaves two important points. The first one is comparisons (<, >, ==) which is easily dealt with the same way as basic arithmetic operations: just using the respective traits from std:

trait Number: Sized + Copy + Add<Self, Output=Self> + Mul<Self, Output=Self> + Sub<Self, Output=Self> + Div<Self, Output=Self> + PartialEq + PartialOrd {}

[!Aside] We are working on floats, so no Eq for us. See the Rust doc for a more detailed explanation on why this is problematic: https://doc.rust-lang.org/std/cmp/trait.Eq.html

Done?

Mathematical Functions

Almost. We still need to cover all methods we need to use on our Number trait, like sin, exp etc. Unfortunately, for these, we don't have ready-made traits in the standard library that we can just add to our trait. Instead, now we have to add a bunch of boilerplate to delegate.

trait Number ...
    fn sin(self) -> Self;
    fn exp(self) -> Self;
    fn powf(self, other: Self) -> Self;
    // ... many more methods
}

and we need to provide an implementation for every type we want to use. We can use shortcuts like macros to ease the pain, but we cannot avoid it completely.

But once we have this, we can cover all basic mathematical operations without any new surprises. And that means we are done here!

We can now do pretty much everything and write something like this:

   fn fancy_maths<T: Number>(a: T, b: T) -> T {
       ( a.powf(T::from(3.0)) + b.powf(T::from(4.0)))
   }

Using the num_traits Crate

Great! Now we know what the pitfalls are, what to look out for and how it works under the hood. But hold on. Before you go off and implement all of this by hand, have a look at the num_traits crate. It implements a few things differently from explained here, but the general ideas hold up. So instead of implementing Number by hand, you can use the num_traits::Float and have access to all you would expect.

This will cover most of the common use cases, so in general there is no need to roll your own traits.

Here's how the examples from this post look using num_traits::Float:

use num_traits::Float;

fn double<T: Float>(value: T) -> T {
    value + value
}

<!-- SPELLING: "circumference" should be "circumference" -->
fn circle_circumference<T: Float>(r: T) -> T {
    T::from(2.0).unwrap() * r * T::PI()
}

fn fancy_maths<T: Float>(a: T, b: T) -> T {
    a.powf(T::from(3.0).unwrap()) + b.powf(T::from(4.0).unwrap())
}

pub fn main() {
    println!("{}", double(4.09_f32));
    println!("{}", circle_circumference(2.0_f64));
    println!("{}", fancy_maths(2.0_f32, 3.0_f32));
}

The num_traits::Float trait already provides all the arithmetic operations, comparison operators, mathematical functions (sin, cos, exp, powf, etc.), and constants you need—no manual implementation required!

Conclusion and Next Steps

This concludes this walkthrough on how to write code that is generic over f32/f64. The next step from here is to be compatible with even more interesting data structures: arrays of numbers (e.g., from ndarray or nalgebra) or data frames (e.g., polars). This step unlocks huge potential to speed up code by doing parallel computations with a relatively simple implementation of your core algorithm. But for that we need to consider memory management, as these data structures are not cheap to copy (and thus do not implement Copy). And that is a topic for another day.

Comments? Email me at thoughts@michaelmauderer.com