Symbolic calculation in F#. Part 2: simplifying algebraic expressions

In the previous post we have shown an F# program that calculates derivatives. It has a room for improvement in its presentation part: resulting expressions often can be rewritten in a simpler form, and F# suits very well to solve such task.

We will be using the same set of expression primitives and active pattern helpers “Op” and “Func”:

type Expression =
    | X
    | Const of float
    | Neg of Expression
    | Add of Expression * Expression
    | Sub of Expression * Expression
    | Mul of Expression * Expression
    | Div of Expression * Expression
    | Pow of Expression * Expression
    | Exp of Expression
    | Log of Expression
    | Sin of Expression
    | Cos of Expression

let (|Op|_|) (x : Expression) =
    match x with
    | Add(e1, e2) -> Some(Add, e1, e2)
    | Sub(e1, e2) -> Some(Sub, e1, e2)
    | Mul(e1, e2) -> Some(Mul, e1, e2)
    | Div(e1, e2) -> Some(Div, e1, e2)
    | Pow(e1, e2) -> Some(Pow, e1, e2)
    | _ -> None

let (|Func|_|) (x : Expression) =
    match x with
    | Exp(e) -> Some(Exp, e)
    | Log(e) -> Some(Log, e)
    | Sin(e) -> Some(Sin, e)
    | Cos(e) -> Some(Cos, e)
    | _ -> None

Now we define a recursive function Simplify that converts an input of Expression type to an output of the same type attempting to rewrite an expression in a simpler form. Like with Derivative function, the definition is a list of rules rather than a sequential algorithm:

let rec Simplify x : Expression =
    match x with
    | Add(Const(n1), Const(n2)) -> Const(n1 + n2)
    | Sub(Const(n1), Const(n2)) -> Const(n1 - n2)
    | Mul(Const(n1), Const(n2)) -> Const(n1 * n2)
    | Div(Const(n1), Const(n2)) -> Const(n1 / n2)
    | Neg(Const(0.)) -> Const(0.)
    | Neg(Neg(e)) -> e |> Simplify
    | Add(e, Const(0.)) -> e |> Simplify
    | Add(Const(0.), e) -> e |> Simplify
    | Add(Const(n), e) -> Add(e, Const(n)) |> Simplify
    | Add(e1, Neg(e2)) -> Sub(e1, e2) |> Simplify
    | Add(Neg(e1), e2) -> Sub(e2, e1) |> Simplify
    | Sub(e, Const(0.)) -> e |> Simplify
    | Sub(Const(0.), e) -> Neg(e) |> Simplify
    | Mul(e, Const(1.)) -> e |> Simplify
    | Mul(Const(1.), e) -> e |> Simplify
    | Mul(e, Const(0.)) -> Const(0.)
    | Mul(Const(0.), e) -> Const(0.)
    | Mul(e, Const(n)) -> Mul(Const(n), e) |> Simplify
    | Mul(Div(Const(n), e1), e2) -> Mul(Const(n), Div(e2, e1)) |> Simplify
    | Mul(e1, Div(Const(n), e2)) -> Mul(Const(n), Div(e1, e2)) |> Simplify
    | Mul(Neg(e1), e2) -> Neg(Mul(e1, e2)) |> Simplify
    | Mul(e1, Neg(e2)) -> Neg(Mul(e1, e2)) |> Simplify
    | Div(Const(0.), e) -> Const(0.)
    | Div(e, Const(1.)) -> e |> Simplify
    | Div(Neg(e1), e2) -> Neg(Div(e1, e2)) |> Simplify
    | Div(e1, Neg(e2)) -> Neg(Div(e1, e2)) |> Simplify
    | Pow(Const(0.), e) -> Const(0.)
    | Pow(Const(1.), e) -> Const(1.)
    | Pow(e, Const(0.)) -> Const(1.)
    | Pow(e, Const(1.)) -> e |> Simplify
    | Op (op, e1, e2) 
        ->
        let e1s = Simplify e1
        let e2s = Simplify e2
        if e1s <> e1 || e2s <> e2 then
            op(Simplify e1, Simplify e2) |> Simplify
        else
            op(e1, e2)
    | _ -> x

The set of rules is not complete. For example, there is no rule to rewrite “2 * (3 + X)” as “6 + 2*X”, but it should be enough to eliminate most obvious redundancies, such as multiplication by one and zero addition. So if we can fire F# interactive window, we can test how it works:

> let s = Simplify(Add(Mul(Const 0., X), Mul(Const 5., Const 1.)));;

val s : Expression = Const 5.0

What we can do now is extend the Derivative function written in the previous post, so it can take advantage of our new Simplify function:

let rec Derivative x : Expression =
    let y =
        match x with
        | X -> Const(1.)
        | Const(n) -> Const(0.)
        | Neg(e) -> Neg(Derivative(e))
        | Add(e1, e2) -> Add(Derivative(e1), Derivative(e2))
        | Sub(e1, e2) -> Sub(Derivative(e1), Derivative(e2))
        | Mul(e1, e2) -> Add(Mul(Derivative(e1), e2), Mul(e1, Derivative(e2)))
        | Pow(e, Const(n)) -> Mul(Const(n), Pow(e, Const(n-1.)))
        | Exp(X) -> Exp(X)
        | Log(X) -> Div(Const(1.), X)
        | Sin(X) -> Cos(X)
        | Cos(X) -> Neg(Sin(X))
        | Func(g, f) -> 
            let dg = Derivative(g(X))
            let df = Derivative(f)
            match dg with
            | Func(dgf, dge) -> Mul(dgf(f), df)
            | Op (op, e1, e2) -> Mul(op(e1, e2), df)
            | _ -> failwith(sprintf "Unable to match compound function [%A]" dg)
        | _ -> failwith(sprintf "Unable to match expression [%A]" x)
    Simplify y

Now let’s test what we have by calculating derivative for various functions:

> let d1 = Derivative(Add(Mul(Const(5.), X), Const(3.)))
let d2 = Derivative(Add(Pow(X, Const(3.)), Const(3.)))
let d3 = Derivative(Sin(Mul(Const(2.), X)))
let d4 = Derivative(Log(Mul(Const(2.), X)))
let d5 = Derivative(Exp(Mul(Const(2.), X)))
let d6 = Derivative(Exp(Pow(X, Const(2.))))
let d7 = Derivative(Log(Sin(X)))
let d8 = Derivative(Log(Cos(X)));;

> 
val d1 : Expression = Const 5.0
val d2 : Expression = Mul (Const 3.0,Pow (X,Const 2.0))
val d3 : Expression = Mul (Const 2.0,Cos (Mul (Const 2.0,X)))
val d4 : Expression = Div (Const 2.0,X)
val d5 : Expression = Mul (Const 2.0,Exp (Mul (Const 2.0,X)))
val d6 : Expression = Mul (Exp (Pow (X,Const 2.0)),Mul (Const 2.0,X))
val d7 : Expression = Div (Cos X,X)
val d8 : Expression = Neg (Div (Sin X,X))

Can we improve anything? Well, maybe… Since it’s getting so easy, why not set an ultimate goal: present input and output as a plain text, not as Expression items. So we can pass “log(cos(x))” and get back “-sin(x)/x”. I guess this will be slightly more work to achieve, but should be fun once we get there. In the next part.

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s