Differentiation

type Expression =
    | Var of string
    | Pi
    | Const of float
    | Add of Expression * Expression
    | Sub of Expression * Expression
    | Mul of Expression * Expression
    | Div of Expression * Expression
    | Pow of Expression * Expression
    | Sqrt of Expression
    | Exp of Expression
    | Ln of Expression
    | Sin of Expression
    | Cos of Expression
    | Tan of Expression
    | Sec of Expression
    | Csc of Expression
    | Cot of Expression
    | Log10 of Expression
    | Log of Expression * Expression
    | Neg of Expression
    | Undefined
    | PositiveInfinity
    | NegativeInfinity

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)
    | Log(e1, e2) -> Some(Log, e1, e2)
    | _ -> None

let (|Func|_|) (x : Expression) =
    match x with
    | Neg(e) -> Some(Neg, e)
    | Sqrt(e) -> Some(Sqrt, e)
    | Exp(e) -> Some(Exp, e)
    | Ln(e) -> Some(Ln, e)
    | Log10(e) -> Some(Log10, e)
    | Sin(e) -> Some(Sin, e)
    | Cos(e) -> Some(Cos, e)
    | Tan(e) -> Some(Tan, e)
    | Csc(e) -> Some(Csc, e)
    | Sec(e) -> Some(Sec, e)
    | Cot(e) -> Some(Cot, e)
    | _ -> None

let rec Evaluate expression =
    match expression with
    | Neg(Const(v)) -> Const(-v)
    | Add(Const(v), Const(v2)) -> Const(v + v2)
    | Sub(Const(v), Const(v2)) -> Const(v - v2)
    | Mul(Const(v), Const(v2)) -> Const(v * v2)
    | Mul(Const(v), Mul(Const(v2), e2)) -> Mul(Const(v * v2), Evaluate e2)
    | Mul(Const(v), Mul(e2, Const(v2))) -> Mul(Const(v * v2), Evaluate e2)
    | Mul(Const(v), Mul(e, Mul(Const(v2), eq))) -> Mul(Mul(Const(v * v2), Evaluate e), Evaluate eq)
    //| Mul(Mul(Const(v), e2), Const(v2)) -> Mul(Const(v*v2), Evaluate e2)
    //| Mul(Mul(e2, Const(v)), Const(v2)) -> Mul(Const(v*v2), Evaluate e2)
    | Div(Const(v), Const(v2)) -> Const(v / v2)
    | Pow(Const(v), Const(v2)) -> Const(v ** v2)
    //| Sqrt(Const(v)) -> Const(sqrt(v))
    //| Log(Const(v), Const(v2)) -> Const(log10(v2) / log10(v))
    //| Log10(Const(v)) -> Const(log10(v))
    //| Ln(Const(v)) -> Const(log(v))
    //| Exp(Const(v)) -> Const(exp(v))
    //| Sin(Const(v)) -> Const(sin(v))
    //| Cos(Const(v)) -> Const(cos(v))
    //| Tan(Const(v)) -> Const(tan(v))
    //| Sec(Const(v)) -> Const(1. / cos(v))
    //| Csc(Const(v)) -> Const(1. / sin(v))
    //| Cot(Const(v)) -> Const(1. / tan(v))
    | Op(op, expr, expr2) -> op(Evaluate expr, Evaluate expr2)
    | Func(func, expr) -> func(Evaluate expr)
    | _ -> expression

let rec Simplify expression =
    match expression with
    | Neg(Const(0.)) -> Const(0.)
    | Neg(Neg(v)) -> v
    | Add(Const(0.), expr) -> Simplify expr
    | Add(expr, Const(0.)) -> Simplify expr
    | Sub(Const(0.), expr) -> Neg(Simplify expr)
    | Sub(expr, Const(0.)) -> Simplify expr
    | Mul(expr, Const(0.)) -> Const(0.)
    | Mul(Const(0.), expr) -> Const(0.)
    | Mul(expr, Const(1.)) -> Simplify expr
    | Mul(Const(1.), expr) -> Simplify expr
    | Div(expr, Const(0.)) -> Undefined
    | Div(Const(0.), expr) -> Const(0.)
    | Div(expr, Const(1.)) -> Simplify expr
    | Div(expr, expr2) when Simplify expr = Simplify expr2 -> Const(1.)
    | Pow(Const(0.), expr) -> Const(0.)
    | Pow(Const(1.), expr) -> Const(1.)
    | Pow(expr, Const(1.)) -> Simplify expr
    | Pow(expr, Const(0.)) -> Const(1.)
    | Op (op, e1, e2)
        ->
        let e1s = Simplify e1
        let e2s = Simplify e2
        if e1s <> e1 || e2s <> e2 then
            op(e1s, e2s) |> Simplify
        else
            op(e1, e2)
    | _ -> expression

let rec Derive expression =
    let deriv =
        match expression with
        | Const(v) -> Const(0.)
        | Var(v) -> Const(1.)
        | Neg(e) -> Neg(Derive(e))
        | Add(left, right) -> Add(Derive(left), Derive(right))
        | Sub(left, right) -> Sub(Derive(left), Derive(right))
        | Mul(Const(v), f) -> Mul(Const(v), Derive(f))
        | Mul(expr, Exp(Mul(expr2, Ln(expr1)))) -> Mul(expr, Pow(expr1, expr2))
        | Mul(Exp(Mul(expr2, Ln(expr1))), expr) -> Mul(expr, Pow(expr1, expr2))
        | Mul(f, g) -> Add(Mul(Derive(f), g), Mul(f, Derive(g)))
        | Div(Const(1.), e) -> Div(Derive(e), Pow(e, Const(2.)))
        | Div(f, g) -> Div(Sub(Mul(Derive(f), g), Mul(Derive(g), f)), Pow(g, Const(2.)))
        | Sqrt(Var(x)) -> Pow(Var(x), Div(Const(1.), Const(2.))) |> Derive
        | Pow(Var(x), Const(g)) -> Mul(Const(g), Pow(Var(x), Const(g - 1.)))
        | Pow(Const(f), g) -> Mul(Ln(Const(f)), Pow(Const(f), g))
        | Ln(Var(x)) -> Div(Const(1.), Var(x))
        | Sin(Var(x)) -> Cos(Var(x))
        | Cos(Var(x)) -> Neg(Sin(Var(x)))
        | Tan(Var(x)) -> Pow(Sec(Var(x)), Const(2.))
        | Sec(Var(x)) -> Mul(Sec(Var(x)), Tan(Var(x)))
        | Csc(Var(x)) -> Neg(Mul(Csc(Var(x)), Cot(Var(x))))
        | Cot(Var(x)) -> Neg(Pow(Csc(Var(x)), Const(2.)))
        | Exp(Var(x)) -> Exp(Var(x))
                | Op(f, g, e) ->
            let f' = Derive(f(Var("x"), e))
            let g' = Derive(g)
            match f' with
            | Func(f'x, expr) -> Mul(g', f'x(g))
            | Op (op, e1, e2) -> Mul(g', op(e1, g))
            | _ -> failwith(sprintf "Unable to match expression [%A]" f')
        | Func(g, f) ->
            let dg = Derive(g(Var("x")))
            let df = Derive(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 "Unsupported expression"
    Simplify deriv

let formatOp (e: Expression) : string =
    match e with
    | Add(e1, e2) -> "+"
    | Sub(e1, e2) -> "-"
    | Mul(e1, e2) -> "*"
    | Div(e1, e2) -> "/"
    | Pow(e1, e2) -> "^"
    | _ -> failwith(sprintf "Unrecognized operator [%A]" e)

let formatFunc (e: Expression) (a : string) : string =
    match e with
    | Sqrt(x) -> sprintf "sqrt(%s)" a
    | Exp(x) -> sprintf "e^(%s)" a
    | Ln(x) -> sprintf "ln(%s)" a
    | Log10(x) -> sprintf "log10(%s)" a
    | Sin(x) -> sprintf "sin(%s)" a
    | Cos(x) -> sprintf "cos(%s)" a
    | Tan(x) -> sprintf "tan(%s)" a
    | Csc(x) -> sprintf "csc(%s)" a
    | Sec(x) -> sprintf "sec(%s)" a
    | Cot(x) -> sprintf "cot(%s)" a
    | _ -> failwith(sprintf "Unrecognized function [%A]" e)

let rec format (inner : Expression) : string =
    match inner with
    | Var(n) -> n
    | Const(n) -> sprintf "%g" n
    | Pi -> "pi"
    | Neg x -> x |> format |> sprintf "-%s"
    | Mul(Pi, Const(v)) -> sprintf "%gpi" v
    | Mul(Const(v), Pi) -> sprintf "%gpi" v
    | Mul(Const(v), Var(n)) -> sprintf "%g%s" v n
    | Pow(expr, expr2) -> sprintf "%s^(%s)" (expr |> format) (expr2 |> format)
    | Op(op, e1, e2) -> sprintf "(%s %s %s)" (e1 |> format) (inner |> formatOp) (e2 |> format)
    | Func(f, e) -> (inner |> formatFunc) (e |> format)
    | _ -> failwith "Could not format"

let pi =
    System.Math.PI

[<EntryPoint>]
let main argv =
    let fx = Add(Pow(Var("x"), Const(3.)), Mul(Const(3.), Pow(Var("x"), Const(2.)))) |> Simplify |> Evaluate
    let fx' = fx |> Derive |> Simplify |> Evaluate |> format
    // cos(sqrt(sin(tan(pi*x))))
    //let fx2 = Cos(Pow(Sin(Tan(Mul(Const(pi), Var("x")))), Div(Const(1.), Const(2.))))
    //let fx2 = Pow(Var("x"), Pow(Sin(Var("x")), Const(2.)))
    let fx2 = Pow(Add(Const(1.), Pow(Var("x"), Const(2.))), Const(100.)) |> Simplify |> Evaluate
    let fx2' = fx2 |> Derive |> Simplify |> Evaluate |> format

    printfn "f(x) = %A\nf(x)' = %A\ng(x) = %A\ng(x)' = %A" (fx |> format) fx' (fx2 |> format) fx2'
    System.Console.ReadLine() |> ignore
    0