The Differentiator

Fun With C# Expressions

Previously:

Differentiating f(x)^g(x)

Note: The latest code can be found at https://github.com/dogbiscuituk/Sid/.

Take a seat with Visual Studio and create a new solution - Console, Winforms, whatever - your choice. Now add a new class library FormulaBuilder, and add to that a static class Differentiator. This is going to let us manipulate objects of type System.Linq.Expression, representing real-valued functions of a single real variable, and ultimately to plot the graph of any such function, as well as its nth derivatives. It will use an easy-to-read fluent interface, while illustrating some key design tenets I'd like to highlight.

Back to the class, and here come our first three requirements. We will need ways to:

1. represent constants and variables, and combine these fluently into compound expressions;
2. display any such expression as an easy-to-read formula, and evaluate it for given values of x; and
3. verify correctness of these combinations and conversions.

Leaving aside the whole "combine these fluently into compound expressions" thing for just a moment, let's open with a ParameterExpression for x, and the following starter set of extension methods:

using System;
using System.Diagnostics;
using System.Linq.Expressions;

namespace FormulaBuilder
{
    public static class Differentiator
    {
        public static ParameterExpression x = Expression.Variable(typeof(double));

        public static ConstantExpression Constant(double c) => Expression.Constant(c);

        public static string AsString(this Expression e, string variableName = "x") =>
            e.ToString().Replace("Param_0", variableName);

        public static Func<double, double> AsFunction(this Expression e) =>
            Expression.Lambda<Func<double, double>>(e, x).Compile();

        public static double AsDouble(this Expression e, double x) => AsFunction(e)(x);

        public static void Check(double expected, double actual) =>
            Check(expected.ToString(), actual.ToString());

        public static void Check(string expected, string actual)
        {
            if (actual == expected)
                Debug.WriteLine("OK: " + actual);
            else
                Debug.WriteLine("*** Comparison failed. Expected: {0}, Actual: {1}.", expected, actual);
        }

On first look, the expression display requirement (1) seems to be met by the built-in ToString() method, available to every type; but that spits out an ugly Param_0 instead of a nice tidy x, so we supply a new AsString() method to tidy that up.

The expression evaluation/iteration requirement (2) means converting our expression into a LambdaExpression, compiling that into a function, and finally executing that function using our chosen values of x. For efficiency I've split this process into two steps. The conversion to a compiled LambdaExpression is performed by the ToFunction() method. Once we have this function, we can call it a million times with different x values, without having to perform a runtime lambda recompilation every time. The ToDouble() method is also provided for those out-of-loop test cases, where we just want to obtain a single value using a single method call.

As for testing requirement (3), the overloaded Check methods provide a line of debug output for each checked operation, beginning with either OK or *** Comparison failed.

Fluency First

Now back to those compounds. Actually, the System.Linq.Expression class already has some fine methods for compounding expressions into BinaryExpressions, but I'm going to deviate a little from these. The target is a fluent syntax for building expressions on the fly, and to that end, I find the names Plus, Minus, Times, Over a little easier to scan than Add, Subtract, Multiply, Divide. So let's add these to the Differentiator class:

        public static BinaryExpression Plus(this Expression f, Expression g) => Expression.Add(f, g);
        public static BinaryExpression Minus(this Expression f, Expression g) => Expression.Subtract(f, g);
        public static BinaryExpression Times(this Expression f, Expression g) => Expression.Multiply(f, g);
        public static BinaryExpression Over(this Expression f, Expression g) => Expression.Divide(f, g);
        public static BinaryExpression Power(this Expression f, Expression g) => Expression.Power(f, g);

Actually I also considered replacing Power with ToThe, but maybe that's just a colloquial idiom too far!

The Constant() syntax is quite clumsy, almost the opposite of the fluent design goal. It would be great to have an implicit operator for converting constants such as 2.0 or Math.PI into ConstantExpressions as required - just as it would be equally great to be able to overload the standard arithmetic operators for Expressions, and be finished with fluent method calls! Sadly, C# restrictions on overloaded operators make both of these dreams impossible, but at least in the Constant() case, what we can do instead is overload the calling methods at the point(s) of invocation. For example, we'd like to be able to write x.Power(2) or x.Power(2.0) instead of x.Power(Constant(2)). This can be achieved by overloading the Power() method, supplying a double in place of an Expression for its second argument, and similarly for the other BinaryExpression methods:

        public static BinaryExpression Plus(this Expression f, double g) => f.Plus(Constant(g));
        public static BinaryExpression Minus(this Expression f, double g) => f.Minus(Constant(g));
        public static BinaryExpression Times(this Expression f, double g) => f.Times(Constant(g));
        public static BinaryExpression Over(this Expression f, double g) => f.Over(Constant(g));
        public static BinaryExpression Power(this Expression f, double g) => f.Power(Constant(g));

This almost gets us where we want to go, but we still have to write e.g. Constant(3).Times(x), when we'd prefer just 3.Times(x). Here it's the first argument of the Times() extension method that needs overloading, and there's a spanner thrown in by C#'s rules of method resolution. Sure, we can supply alternate methods with a double in that position, but where previously we could rely on an int argument being promoted to double, no such luck this time. If we insist (and we do!) upon being able to write 3.Times(x), as opposed to workarounds like x.Times(3), 3.0.Times(x), or worst of all ((double)3).Times(x) - which seems hardly an improvement! - then we must supply separate double and int versions.

        public static BinaryExpression Plus(this double f, Expression g) => Constant(f).Plus(g);
        public static BinaryExpression Minus(this double f, Expression g) => Constant(f).Minus(g);
        public static BinaryExpression Times(this double f, Expression g) => Constant(f).Times(g);
        public static BinaryExpression Over(this double f, Expression g) => Constant(f).Over(g);
        public static BinaryExpression Power(this double f, Expression g) => Constant(f).Power(g);

        public static BinaryExpression Plus(this int f, Expression g) => Constant(f).Plus(g);
        public static BinaryExpression Minus(this int f, Expression g) => Constant(f).Minus(g);
        public static BinaryExpression Times(this int f, Expression g) => Constant(f).Times(g);
        public static BinaryExpression Over(this int f, Expression g) => Constant(f).Over(g);
        public static BinaryExpression Power(this int f, Expression g) => Constant(f).Power(g);

        public static UnaryExpression Negate(this Expression e) => Expression.Negate(e);
        public static BinaryExpression Reciprocal(this Expression e) => e.Power(-1);
        public static BinaryExpression Squared(this Expression e) => e.Power(2);
        public static BinaryExpression Cubed(this Expression e) => e.Power(3);

        public static void TestCompoundExpression()
        {
            var f = x.Squared().Plus(3.Times(x)).Minus(5);    // f(x) = x²+3x-5
            Check("(((x ^ 2) + (3 * x)) - 5)", f.AsString()); // Check the built expression formula
            Check(Math.Pow(7, 2) + 3 * 7 - 5, f.AsDouble(7)); // Check the expression value at x = 7 (should be 65)
        }

With the addition of these few minor utilities, we now have the means to build polynomials and other arithmetic expressions in x. Give your app a project reference to FormulaBuilder, add it to your using clauses, then add a call to Differentiator.TestCompoundExpression() somewhere in your main code and run it. You should see the following two lines of reassurance in your debug output:

OK: (((x ^ 2) + (3 * x)) - 5)
OK: 65

These show that the quadratic expression x² + 3x - 5 was built successfully, and that it correctly evaluates to 65 at x = 7.

Function Calls

Arithmetic operations are one way of building out our expressions, but what about invoking functions? What use a differentiator that can't tell you the derivative of sin is cos, while that of cos is -sin?

Elementary functions in the System.Linq.Expression hierarchy are represented by the MethodCallExpression class. The Function() code below generates such an expression from the correspondingly named method of the Math class, obtained via reflection. Below it are listed a generous dozen or so concrete functions, all which happen to be implemented in Math; and below these, a further group showing how functions "missing" from Math, such as the reciprocal trigonometrics and hyperbolics, as well as their inverses, can be built out of expressions in the elementary ones.

Naturally, the input to any functional expression can itself be an arbitrary expression tree, maybe containing further nested function calls, etc.

        public static MethodCallExpression Function(string functionName, Expression e) =>
            Expression.Call(typeof(Math).GetMethod(functionName, new[] { typeof(double) }), e);

        public static MethodCallExpression Abs(this Expression e) => Function("Abs", e);
        public static MethodCallExpression Sqrt(this Expression e) => Function("Sqrt", e);
        public static MethodCallExpression Exp(this Expression e) => Function("Exp", e);
        public static MethodCallExpression Log(this Expression e) => Function("Log", e);
        public static MethodCallExpression Log10(this Expression e) => Function("Log10", e);
        public static MethodCallExpression Sin(this Expression e) => Function("Sin", e);
        public static MethodCallExpression Cos(this Expression e) => Function("Cos", e);
        public static MethodCallExpression Tan(this Expression e) => Function("Tan", e);
        public static MethodCallExpression Asin(this Expression e) => Function("Asin", e);
        public static MethodCallExpression Acos(this Expression e) => Function("Acos", e);
        public static MethodCallExpression Atan(this Expression e) => Function("Atan", e);
        public static MethodCallExpression Sinh(this Expression e) => Function("Sinh", e);
        public static MethodCallExpression Cosh(this Expression e) => Function("Cosh", e);
        public static MethodCallExpression Tanh(this Expression e) => Function("Tanh", e);

        public static BinaryExpression Sec(this Expression e) => Reciprocal(Cos(e));
        public static BinaryExpression Csc(this Expression e) => Reciprocal(Sin(e));
        public static BinaryExpression Cot(this Expression e) => Reciprocal(Tan(e));
        public static MethodCallExpression Asec(this Expression e) => Acos(Reciprocal(e));
        public static MethodCallExpression Acsc(this Expression e) => Asin(Reciprocal(e));
        public static MethodCallExpression Acot(this Expression e) => Atan(Reciprocal(e));
        public static BinaryExpression Sech(this Expression e) => Reciprocal(Cosh(e));
        public static BinaryExpression Csch(this Expression e) => Reciprocal(Sinh(e));
        public static BinaryExpression Coth(this Expression e) => Reciprocal(Tanh(e));
        public static MethodCallExpression ArcSinh(this Expression e) => Log(e.Plus(Sqrt(e.Squared().Plus(1))));
        public static MethodCallExpression ArcCosh(this Expression e) => Log(e.Plus(Sqrt(e.Squared().Minus(1))));
        public static BinaryExpression ArcTanh(this Expression e) => Log(1.Plus(e).Over(1.Minus(e))).Over(2);
        public static MethodCallExpression ArcSech(this Expression e) => Log(1.Plus(Sqrt(1.Minus(e.Squared()))).Over(e));
        public static MethodCallExpression ArcCsch(this Expression e) => Log(1.Plus(Sqrt(1.Plus(e.Squared()))).Over(e));
        public static BinaryExpression ArcCoth(this Expression e) => Log(e.Plus(1).Over(e.Minus(1))).Over(2);

        public static void TestTrigonometricExpression()
        {
            var f = Sin(x).Squared().Plus(Cos(x).Squared());             // f(x) = sin²x+cos²x
            Check("((Sin(x) ^ 2) + (Cos(x) ^ 2))", f.AsString());        // Check the built expression formula
            double p3 = Math.PI / 3, s = Math.Sin(p3), c = Math.Cos(p3); // Check the expression value at x = π/3
            Check(Math.Pow(s, 2) + Math.Pow(c, 2), f.AsDouble(p3));      // (should be 1)
        }

In your app level code, add a call to Differentiator.TestTrigonometricExpression(), which builds the expression sin²x + cos²x: the LHS of the most famous trigonometric identity, whose RHS is 1 for all values of x. The new test method checks, first using AsString(), that the correct expression has in fact been built; then, using AsDouble(), that it evaluates to 1 as expected, in this case when supplied with π/3 as a parameter. At runtime, your debug window should display a further two comforting lines of verification:

OK: ((Sin(x) ^ 2) + (Cos(x) ^ 2))
OK: 1

Finally, A Derivative

Can't finish this article without a look at some actual differentiation, so let's add the function D() that differentiates expressions. It needs to handle just these four cases:

• a constant;
• x itself;
• a function call;
• an artithmetic operation (unary or binary).

The third and fourth cases here have potential for recursion. At long last, here is the code to differentiate any expression:

        public static Expression D(this Expression e)
        {
            if (e is ConstantExpression) return Constant(0);  // d(c)/dx = 0
            if (e is ParameterExpression) return Constant(1); // d(x)/dx = 1
            if (e is MethodCallExpression m)
            {
                var f = m.Arguments[0];
                return DifferentiateFunction(m.Method.Name, f).Times(D(f)); // Chain Rule
            }
            if (e is UnaryExpression u)
            {
                var f = D(u.Operand);
                return u.NodeType == ExpressionType.UnaryPlus ? f : Negate(f);
            }
            if (e is BinaryExpression b)
            {
                Expression f = b.Left, g = b.Right;
                switch (b.NodeType)
                {
                    case ExpressionType.Add:      // (f+g)' = f'+g'
                        return D(f).Plus(D(g));
                    case ExpressionType.Subtract: // (f-g)' = f'-g'
                        return D(f).Minus(D(g));
                    case ExpressionType.Multiply: // (fg)' = f'g+fg'
                        return D(f).Times(g).Plus(f.Times(D(g)));
                    case ExpressionType.Divide:   // (f÷g)' = (f'g-fg')÷(g^2) = f'÷g-fg'÷(g^2)
                        return D(f).Over(g).Minus(f.Times(D(g)).Over(g.Squared()));
                    case ExpressionType.Power:    // (f^g)' = (f^g)*(f'g÷f+g'Log(f)) = (f'g+fg'Log(f))f^(g-1)
                        return D(f).Times(g).Plus(f.Times(D(g)).Times(Log(f))).Times(f.Power(g.Minus(1)));
                }
            }
            throw new InvalidOperationException();
        }

        private static Expression DifferentiateFunction(string methodName, Expression x)
        {
            switch (methodName)
            {
                case "Abs": return x.Over(Abs(x));                            // d(|x|)/dx = x/|x|
                case "Sqrt": return x.Power(-0.5).Over(2);                    // d(√x)/dx = 1/(2√x)
                case "Exp": return Exp(x);                                    // d(eˣ)/dx = eˣ
                case "Log": return Reciprocal(x);                             // d(ln(x))/dx = 1/x
                case "Log10": return Reciprocal(x).Times(Math.Log10(Math.E)); // d(log₁₀(x))/dx = log₁₀(e)/x
                case "Sin": return Cos(x);                                    // d(sin(x))/dx = cos(x)
                case "Cos": return Negate(Sin(x));                            // d(cos(x))/dx = -sin(x)
                case "Tan": return Cos(x).Power(-2);                          // d(tan(x))/dx = sec²(x)
                case "Asin": return 1.Minus(x.Squared()).Power(-0.5);         // d(arcsin(x))/dx = 1/√(1-x²)
                case "Acos": return Negate(1.Minus(x.Squared()).Power(-0.5)); // d(arccos(x))/dx = -1/√(1-x²)
                case "Atan": return 1.Over(1.Plus(x.Squared()));              // d(arctan(x))/dx = 1/(1+x²)
                case "Sinh": return Cosh(x);                                  // d(sinh(x))/dx = cosh(x)
                case "Cosh": return Sinh(x);                                  // d(cosh(x))/dx = sinh(x)
                case "Tanh": return 1.Minus(Tanh(x).Squared());               // d(tanh(x))/dx = 1-tanh²(x)
            }
            throw new InvalidOperationException();
        }

DifferentiateFunction() provides the basic template for a function's derivative, without actually applying the chain rule to its argument Expression. That's why it has been implemented as a private helper function; calling it out of context could result in incorrect results if chain rule application was omitted by mistake.

Algebraic Simplification

Still we're not quite ready to supply the final batch of unit tests, stick a fork in it and call it done. While the D() method may be returning strictly correct results, good enough for numeric applications such as graphing, algebraically they are a bit of a mess. For example, differentiating any candidate containing a function call, even one as simple as Sin(x), results in a clumsy *1 becoming appended to the output expression. Some level of simplification is clearly needed to provide visually acceptable results. For this reason, D() gets called not directly, but from within a wrapper function Differentiate() which applies a layer of simplification before returning the resultant Expression:

        public static Expression Differentiate(this Expression e) => Simplify(D(e));

        public static Expression Simplify(this Expression e)
        {
            if (e is MethodCallExpression m) return SimplifyMethodCall(m);
            if (e is UnaryExpression u) return SimplifyUnary(u);
            if (e is BinaryExpression b) return SimplifyBinary(b);
            return e;
        }

        public static Expression SimplifyMethodCall(MethodCallExpression m)
        {
            var operand = Simplify(m.Arguments[0]);
            if (operand is ConstantExpression ce)
            {
                var c = (double)ce.Value;
                return Constant(Function(m.Method.Name, x).AsDouble(c));
            }
            return Function(m.Method.Name, operand);
        }

        public static Expression SimplifyUnary(UnaryExpression u)
        {
            var operand = Simplify(u.Operand);
            if (operand is ConstantExpression ce)
            {
                var c = (double)ce.Value;
                return Constant(u.NodeType == ExpressionType.UnaryPlus ? c : -c);
            }
            return u.NodeType == ExpressionType.UnaryPlus ? operand : Negate(operand);
        }

        public static Expression SimplifyBinary(BinaryExpression b)
        {
            b = Expression.MakeBinary(b.NodeType, Simplify(b.Left), Simplify(b.Right));
            if (b.Left is ConstantExpression ce1)
            {
                var c1 = (double)ce1.Value;
                if (b.Right is ConstantExpression ce2)
                    return Constant(b.NodeType, c1, (double)ce2.Value);
                switch (b.NodeType)
                {
                    case ExpressionType.Add:                   // c + x → x + c
                    case ExpressionType.Multiply:              // c * x → x * c
                        b = Expression.MakeBinary(b.NodeType, b.Right, b.Left);
                        break;
                    case ExpressionType.Subtract:
                        if (c1 == 0) return Negate(b.Right);   // 0 - x → -x
                        break;
                    case ExpressionType.Divide:
                        if (c1 == 0) return Constant(0);       // 0 ÷ x → 0
                        break;
                    case ExpressionType.Power:
                        if (c1 == 0) return Constant(0);       // 0 ^ x → 0
                        if (c1 == 1) return Constant(1);       // 1 ^ x → 1
                        break;
                }
            }
            if (b.Right is ConstantExpression ce3)
            {
                var c3 = (double)ce3.Value;
                switch (b.NodeType)
                {
                    case ExpressionType.Add:
                    case ExpressionType.Subtract:
                        if (c3 == 0) return b.Left;            // x ± 0 → x
                        break;
                    case ExpressionType.Multiply:
                        if (c3 == 0) return Constant(0);       // x * 0 → 0
                        if (c3 == 1) return b.Left;            // x * 1 → x
                        if (c3 == -1) return Negate(b.Left);   // x * -1 → -x
                        break;
                    case ExpressionType.Divide:
                        if (c3 == 0)
                            throw new DivideByZeroException();
                        if (c3 == 1) return b.Left;            // x ÷ 1 → x
                        if (c3 == -1) return Negate(b.Left);   // x ÷ -1 → -x
                        break;
                    case ExpressionType.Power:
                        if (c3 == 0) return Constant(1);       // x ^ 0 → 1
                        if (c3 == 1) return b.Left;            // x ^ 1 → x
                        break;
                }
                switch (b.NodeType)
                {
                    case ExpressionType.Add:                   // (x + c1) ± c2 → x + (c1 ± c2)
                    case ExpressionType.Subtract:              // (x - c1) ± c2 → x - (c1 ∓ c2)
                    case ExpressionType.Multiply:              // (x * c1) */ c2 → x * (c1 */ c2)
                    case ExpressionType.Divide:                // (x / c1) */ c2 → x / (c1 /* c2)
                        if (!(b.Left is BinaryExpression c && c.Right is ConstantExpression ce4))
                            break;
                        bool same = c.NodeType == b.NodeType, opps = c.NodeType == Invert(b.NodeType);
                        if (!(same || opps))
                            break;
                        var c4 = (double)ce4.Value;
                        var nodeType = c.NodeType == ExpressionType.Add || c.NodeType == ExpressionType.Subtract
                            ? ExpressionType.Add
                            : ExpressionType.Multiply;
                        if (opps) nodeType = Invert(nodeType);
                        return Expression.MakeBinary(c.NodeType, c.Left, Constant(nodeType, c4, c3));
                }
            }
            return b;
        }

        public static ExpressionType Invert(this ExpressionType nodeType)
        {
            switch (nodeType)
            {
                case ExpressionType.Add: return ExpressionType.Subtract;
                case ExpressionType.Subtract: return ExpressionType.Add;
                case ExpressionType.Multiply: return ExpressionType.Divide;
                case ExpressionType.Divide: return ExpressionType.Multiply;
            }
            throw new InvalidOperationException();
        }

        public static ConstantExpression Constant(ExpressionType nodeType, double c, double d) =>
            Constant(Apply(nodeType, c, d));

        public static double Apply(ExpressionType nodeType, double c, double d)
        {
            switch (nodeType)
            {
                case ExpressionType.Add: return c + d;
                case ExpressionType.Subtract: return c - d;
                case ExpressionType.Multiply: return c * d;
                case ExpressionType.Divide:
                    if (d == 0)
                        throw new DivideByZeroException();
                    return c / d;
                case ExpressionType.Power:
                    if (c == 0 && d == 0)
                        throw new InvalidOperationException();
                    return Math.Pow(c, d);
            }
            throw new InvalidOperationException();
        }

This Simplify() method does little more than apply some arithmetic identity, associativity and commutativity results, to remove some of the redundancy littering the differentiator's output. That's why I've implemented it here as simple, direct logic. There are many more simplification rules to consider, and as such rules are added, it becomes clear that a rule application engine is needed, allowing extension of the list without spaghettification of the logic.

Now we can at least start writing an acceptable suite of unit tests. Here are some examples:

        public static void TestSimplify(Expression e, string expected) => Check(expected, Simplify(e).AsString());

        public static void TestSimplifications()
        {
            TestSimplify(x.Plus(6).Plus(2), "(x + 8)");
            TestSimplify(6.Plus(x).Plus(2), "(x + 8)");
            TestSimplify(6.Plus(x.Plus(2)), "(x + 8)");
            TestSimplify(x.Plus(6).Minus(2), "(x + 4)");
            TestSimplify(6.Plus(x).Minus(2), "(x + 4)");
            TestSimplify(x.Minus(6).Plus(2), "(x - 4)");
            TestSimplify(x.Minus(6).Minus(2), "(x - 8)");
            TestSimplify(x.Times(6).Times(2), "(x * 12)");
            TestSimplify(6.Times(x).Times(2), "(x * 12)");
            TestSimplify(6.Times(x.Times(2)), "(x * 12)");
            TestSimplify(x.Times(6).Over(2), "(x * 3)");
            TestSimplify(6.Times(x).Over(2), "(x * 3)");
            TestSimplify(x.Over(2).Times(8), "(x / 0.25)");
            TestSimplify(x.Over(2).Over(8), "(x / 16)");
        }

        public static void TestDerivative(Expression e, string expected) => Check(expected, Differentiate(e).AsString());

        public static void TestFunctionDerivatives()
        {
            TestDerivative(Abs(x), "(x / Abs(x))");                       // d(|x|)/dx = x/|x|
            TestDerivative(Sqrt(x), "((x ^ -0.5) / 2)");                  // d(√x)/dx = 1/(2√x)
            TestDerivative(Exp(x), "Exp(x)");                             // d(eˣ)/dx = eˣ
            TestDerivative(Log(x), "(x ^ -1)");                           // d(ln(x))/dx = 1/x
            TestDerivative(Log10(x), "((x ^ -1) * 0.434294481903252)");   // d(log₁₀(x))/dx = log₁₀(e)/x
            TestDerivative(Sin(x), "Cos(x)");                             // d(sin(x))/dx = cos(x)
            TestDerivative(Cos(x), "-Sin(x)");                            // d(cos(x))/dx = -sin(x)
            TestDerivative(Tan(x), "(Cos(x) ^ -2)");                      // d(tan(x))/dx = sec²(x)
            TestDerivative(Asin(x), "((1 - (x ^ 2)) ^ -0.5)");            // d(arcsin(x))/dx = 1/√(1-x²)
            TestDerivative(Acos(x), "-((1 - (x ^ 2)) ^ -0.5)");           // d(arccos(x))/dx = -1/√(1-x²)
            TestDerivative(Atan(x), "(1 / ((x ^ 2) + 1))");               // d(arctan(x))/dx = 1/(x²+1)
            TestDerivative(Sinh(x), "Cosh(x)");                           // d(sinh(x))/dx = cosh(x)
            TestDerivative(Cosh(x), "Sinh(x)");                           // d(cosh(x))/dx = sinh(x)
            TestDerivative(Tanh(x), "(1 - (Tanh(x) ^ 2))");               // d(tanh(x))/dx = 1-tanh²(x)
        }

        public static void TestPolynomialDerivative()
        {
            TestDerivative(x.Power(4).Minus(3.Times(x.Cubed())).Plus(6.Times(x.Squared())).Minus(3.Times(x)).Plus(1),
                "(((((x ^ 3) * 4) - ((x ^ 2) * 9)) + (x * 12)) - 3)");    // d(x⁴-3x³+6x²-3x+1)/dx = 4x³-9x²+12x-3
        }

        public static void TestChainRule()
        {
            TestDerivative(Exp(x.Squared()), "(Exp((x ^ 2)) * (x * 2))"); // d(exp(x²))/dx = exp(x²)*2x
            TestDerivative(Log(Sin(x)), "((Sin(x) ^ -1) * Cos(x))");      // d(ln(sin(x)))/dx = cot(x)
            TestDerivative(Tan(x.Cubed().Plus(8.Times(x))),               // d(tan(x³+8x))/dx = sec²(x³+8x)*(3x²+8)
                "((Cos(((x ^ 3) + (x * 8))) ^ -2) * (((x ^ 2) * 3) + 8))");
        }

        public static void TestAll()
        {
            TestSimplifications();
            TestCompoundExpression();
            TestTrigonometricExpression();
            TestFunctionDerivatives();
            TestPolynomialDerivative();
            TestChainRule();
        }
    }
}

Next time, by overwhelming public demand (one request has flooded in), I'll be exercising these extension methods in a simple differential graphing application.

Acknowledgements

• Thanks are due to my colleague Stuart Stein for his code review with useful comments, particularly the suggestion that I should embrace the expression body syntax.
• Thanks also to Pavel Vladov for the free use of his C# Syntax Highlighter to format the code snippets in this article.

Differentiating f(x)^g(x)

Black Pixel Red Pixel

My favourite mathematical "YouTutor" is Steve Chow, the guy behind the superlative blackpenredpen collection of mathematics videos on YouTube. Steve has developed a way of explaining proofs using two or more pen colours, hence the name of his YouTube channel. He holds multiple marker pens in one hand while writing on a whiteboard, deftly switching colours so as to highlight the essential differences between each pair of lines in the mathematical exposition. It's a very effective expository technique.

Until recently, my number one blackpenredpen video was this one, where Steve uses logarithmic differentiation and the chain rule to prove the power, product, and quotient rules of derivatives. The reason it's so high on my list is because of the additional bonus content at the end, kind of a "hidden track" type deal, where he goes on to solve for the derivative of a function raised to the power of another function.



Just recently Steve, a keen runner of physical marathons, published a truly marathon mathematical video - my new favourite - in which he performs a hundred and one integration proofs in a single take. I'll let you search for that 6 hour video yourself, and if you feel like sitting through it, may I suggest you limit your accompanying playlist to just these three songs: Steely Dan's "My Old School", Supertramp's "School", and Bowling For Soup's "High School Never Ends". This will help you digest the integration marathon in healthy 15 minute chunks.

Getting back to the bonus differentiation, I remember thinking, this procedure has the effect of subtracting one from the power of one of the functions. Yet it doesn't use the power rule in its own derivation. Aha, I bet this could be used to derive the power rule itself! And since I was looking for a suitable subject for my first foray into maths typesetting, using LaTeX in the Overleaf environment, that became my target.

Here is the resultant PDF. You'll notice the homage to blackpenredpen in some of the text colouring.



The actual Overleaf project which produced this PDF can be viewed here, if you're interested in examining the source (LaTeX is a markup language, just as HTML is a markup language).

This is the first in a series of articles looking at aspects of differentiation. Next time I'll be presenting a neat & simple differentiator, utilising the Expression Trees feature of Microsoft C#.

Kedgeree

My Curry Here

Smoked haddock kedgeree, with added petits pois and raisins why. Because they seem to go nicely. This recipe has been refined for more than a quarter century (I made it once in 1993, then suddenly again today). Wife says it's terrific, and she's a fantastic cook, so there's always that.

Preparation time: 15 minutes.
Cook time: 40 minutes.
Serves 4.
Main Ingredients:

• 300g smoked fillet of haddock
• 4 eggs
• 100g frozen petits pois
• 300 ml milk
• bunch of coriander
• 2 bay leaves

Rice Ingredients:

• 300g long grain easy cook rice (well rinsed)
• 1 large onion (finely chopped)
• 30g raisins
• 1 tsp ground coriander
• 1 tsp ground turmeric
• 2 tsp roasted curry powder
• 1 chicken or vegetable stock cube
• 2 tbsp rapeseed oil

Method:

1. Very important: partition your ingredients into rice and non-rice work areas!

2. Fry the onion in rapeseed oil until translucent (about 5 minutes).

3. Add the spices & salt seasoning; stir and continue frying for a couple of minutes.
4. Add the rice, raisins, crumbled stock cube, and 300ml of water.
5. Cover and simmer for 5 minutes, while...
6. ...poaching the haddock and bay in 300ml of milk, until flaky.

7. Drain the milk from the fish pan into the rice, then roughly flake the fish.
8. Cover the rice pan and remove from the heat.
9. Add the petits pois to the fish (they'll defrost almost instantly).
10. Boil the eggs for 4½ minutes, then plunge them into cold water.
11. After 5 minutes, peel and quarter the eggs, then add to the fish along with the rice and coriander.

12. Mix gently, season with black pepper, then serve with any additional herbage
(e.g. parsley, if you like that sort of thing).

Smiddy Shank Soup

A Photo Essay

A very hearty and rustic Ham Hough Red Lentil Soup. Seriously, if they had this recipe in Zelda: Breath of the Wild, it would break the cookery system. It brings health, stamina, attack and defence boosts, while simultaneously providing resistance to heat, cold and electricity. One plate of it would fetch at least two hundred and eighty rupees. Pro tip: for an additional speed boost, replace the carrots with "swift carrots".

Preparation time: 20 minutes.
Cook time: 3 hours 15 minutes.
Serves 8.
Ingredients:

• 800g Blair Drummond Smiddy smoked ham hough (shank)
• 2 baking potatoes
• 1 large onion
• 3 celery stalks with leaves
• 6 carrots
• 250g red split lentils (well rinsed)
• 4 garlic cloves
• 2 stock cubes (vegetable or chicken)
• 6-8 black peppercorns
• 1 teaspoon dried oregano
• 1 bouquet garni
• 2-3 bay leaves

Method:

1. To 3 litres of boiling water in a pot, add the ham, bouquet garni, peppercorns and bay.

2. Boil for 2½ hours.

3. Meanwhile, back at the vegetable station...

...peel the potatoes, onion, celery, carrots & garlic.

4. Dice the potatoes, onion, celery, and four of the carrots.

5.Crush the garlic and grate the remaining two carrots.

6. At the 2½ hour mark, remove the ham from the pot.

7. Skim the water, then add the potatoes...

...the chopped carrots...

...the celery...

...the onion...

...the grated carrots...

...the crushed garlic...

...the oregano and stock cubes...

...and finally the rinsed lentils. Simmer for 45 minutes.

8. Strip the ham off the shank and add back into the soup before serving.

9. It sticks tae yer ribs. Serve with a twist of pepper.

Acknowledgements

I am gratefully indebted to my beautiful assistant and arm model @catsjamas for her help with the visual production aspects of this project.