Numerical evaluation

Symbolica support efficient numerical evaluations of expressions. We will go through various alternatives and mention their relative advantages.

Evaluating once

To evaluate an expression once, use evaluate. In Python this evaluates with float and in Rust the type of the evaluation is inferred (one could evaluate with 4x4 simd floats or arbitrary-precision floats for example).

from symbolica import *
x, f, cos = Expression.symbol('x', 'f', 'cos')
e = cos(x)*3 + f(x,2)
print(e.evaluate({x: 1}, {f: lambda args: args[0]+args[1]}))
print(e.evaluate_complex({x: 1+2j}, {f: lambda args: args[0]+args[1]}))
use symbolica::{atom::Atom, state::State};

fn main() {
    let a = Atom::parse("cos(x)*3 + f(x,2)").unwrap();
    let mut h = HashMap::default();
    h.insert(State::get_symbol("x"), 1.);
    let mut fn_map = HashMap::default();
    fn_map.insert(State::get_symbol("f"), Box::new(|x, y| x + y));
    println!("{}", a.evaluate(&h, &fn_map));

    let mut h = HashMap::default();
    h.insert(State::get_symbol("x"), Complex::new(1., 2.));
    let mut fn_map = HashMap::default();
    fn_map.insert(State::get_symbol("f"), Box::new(|x, y| x + y));
    println!("{}", a.evaluate(&h, &fn_map));
}

The first argument of evaluate is a map from an atom to a constant and the second argument is a map of a function symbol to a closure that is called with the evaluated arguments.

The above example yields

4.620906917604419
(9.098169021058997-7.155693397455401j)

Evaluators

Naively evaluating an expression can be prohibitively slow for large expressions, especially when it has to be evaluated many times, such as for Monte Carlo simulations/integrations. Symbolica provides state-of-the-art ways to rewrite an expression in a form that makes it very fast to evaluate, often 10-100 times faster than a naive evaluation. Note that during these substitutions, “fast math” operations are applied.

It tries to find the best multivariate Horner scheme, performs common subexpression elimination and recycles intermediate variables and registers. It resembles what a programming language compiler such as gcc would do.

Nested expressions

Such an evaluator supports nested expressions; expression that call other expressions. For example, to evaluate \(e\) in

\[ \begin{align*} e &= x + \pi + \cos(x) + f(g(x+1),2x)\\ f(y,z) &=y^2 + z^2 y^2\\ g(y) &=x+y+5 \end{align*} \]

we use the following code:

from symbolica import *
x, y, z, pi, f, g = Expression.symbol(
    'x', 'y', 'z', 'pi', 'f', 'g')

e1 = Expression.parse("x + pi + cos(x) + f(g(x+1),x*2)")
fd = Expression.parse("y^2 + z^2*y^2")
gd = Expression.parse("x + y + 5")

ev = e1.evaluator({pi: Expression.num(22)/7},
             {(f, "f", (y, z)): fd,
              (g, "g", (y, )): gd},
            [x])
res = ev.evaluate([[5.]])
print(res)

Note that at the moment no floats may appear in the expressions or in the map. Use rationalize_coefficients() to convert them.

fn main() {
    let e1 = Atom::parse("x + pi + cos(x) + f(g(x+1),x*2)").unwrap();
    let f = Atom::parse("y^2 + z^2*y^2").unwrap();
    let g = Atom::parse("x+y+5").unwrap();

    let mut fn_map = FunctionMap::new();

    fn_map.add_constant(
        Atom::new_var(State::get_symbol("pi")),
        Rational::from((22, 7)).into(),
    );
    fn_map
        .add_function(
            State::get_symbol("f"),
            "f".to_string(),
            vec![State::get_symbol("y"), State::get_symbol("z")],
            f.as_view(),
        )
        .unwrap();
    fn_map
        .add_function(
            State::get_symbol("g"),
            "g".to_string(),
            vec![State::get_symbol("y")],
            g.as_view(),
        )
        .unwrap();

    let params = vec![Atom::parse("x").unwrap()];

    let evaluator = e1
        .evaluator(&fn_map, &params, OptimizationSettings::default())
        .unwrap();

    let mut e_f64 = evaluator.map_coeff(&|x| x.into());
    let r = e_f64.evaluate_single(&[5.]);
    println!("{}", r);
}

Which yields 25864.42651932832.

The evaluator takes a hashmap of functions with the following key:

(function symbol, function name, argument list)

where the function_name should be a C++-compatible name.

Multiple expressions

Multiple expressions can be turned into an evaluator at the same time, which means they will be optimized together. As a result, common subexpressions between the expressions will be computed only once.

from symbolica import *
x = Expression.symbol('x')

e1 = x * x + 5
e2 = x * x + 6

ev = Expression.evaluator_multiple([e1, e2], {}, {}, [x])
res = ev.evaluate([[5.]])
print(res)

fn main() {
    let x = Atom::parse("x").unwrap();
    let e1 = &x * &x + 5;
    let e2 = &x * &x + 6;

    let evaluator = Atom::evaluator_multiple(
        &[e1.as_view(), e2.as_view()],
        &FunctionMap::new(),
        &[x],
        OptimizationSettings::default(),
    )
    .unwrap();

    let mut e_f64 = evaluator.map_coeff(&|x| x.into());
    let mut out = vec![0., 0.];
    e_f64.evaluate(&[5.], &mut out);
    println!("{:?}", out);
}

Yields 30.0, 31.0.

Compiled evaluators

Evaluators can be turned into highly optimized C++ code, optionally with inline assembly (x86_64 for now). The latter may give better performance, but it will also seriously cut down on C++ compilation time, effectively making it independent on the size of the expression.

The compiled library accepts evaluations with double and complex<double> and can be directly loaded into Symbolica.

comp = eval.compile('my_fun', 'test.cpp', 'test.so')
res = comp.evaluate([[5.]])
print(res)
let mut compiled = e_f64
    .export_cpp("test.cpp", "my_fun", true, InlineASM::X64)
    .unwrap()
    .compile("test.so", CompileOptions::default())
    .unwrap()
    .load()
    .unwrap();

let mut out = vec![0.];
compiled.evaluate(&[5.], &mut out);
println!("{}", out[0]);