Symbolic manipulation

Arjen Markus (14 august 2003) It was just an idea: symbolic manipulation or manipulation of mathematical expressions via Tcl. This could be handy for some mathematical playthings I have already. I did not want to spend the time and effort (a true, full-sized system would be as difficult to build as, say, hunting a snark, so I wondered if someone else had. The newsgroup did not bring any suggestions and with Google, I saw naught but implementations in Java, C++, LISP.

I pondered upon how to proceed with such system. I gave myself about an hour to actually implement and test it. Here is the result. It lacks a front-end, that will be the most difficult thing to get right. And if you run the test cases, you will see extra, unnecessary parentheses, but I did give myself one hour for the whole thing :).

As RS mentioned, the technique used in Expression parsing will come in handy for building the front-end.

AM (25 april 2005) This morning I came across a completely different technique to compute the derivative of a function - rather than manipulate the expression by which the function is computed you use the rules of differentiation and go from there (a floating-point value must become a value and "its" derivative, but it can be done easily ...). If you are interested: the page Automatic differentiation is meant to clarify the idea.


 # symbolic.tcl --
 #
 #    Experiment with symbolic (algebraic) manipulation:
 #    determine the derivative of an expression
 #
 #    Note:
 #    Parsing an expression like {$x/(1+$x)} is not implemented yet.
 #    So you will have to do that yourself:
 #    {/ $x {+ 1 $x}}
 #    Derivative w.r.t. x:
 #    {1/(1+$x)-$x/((1+$x)*(1+$x)) - not the simplest form, but
 #    [expr] won't mind
 #

 # toexpr --
 #    Reconstruct the expression from the parsed and manipulated
 #    form
 # Arguments:
 #    parse_tree     The parsed and manipulated expression
 # Result:
 #    String representing the expression as [expr] would take it
 # Note:
 #    The result is safe, as far as brackets are concerned, but
 #    very conservative
 #
 proc toexpr { parse_tree } {
    if { [llength $parse_tree] == 1 } {
       return $parse_tree
    } else {
       foreach {op operand1 operand2} $parse_tree {break}
       switch -- $op {
       "umin"  {return "-[toexpr2 - $operand1]"}
       "-"     {return "[toexpr2 - $operand1]-[toexpr2 - $operand2]"}
       "+"     {set add1 [toexpr2 + $operand1]
                set add2 [toexpr2 + $operand2]
                if { $add1 == "0" || $add1 == "(0)" } {
                   return $add2
                }
                if { $add2 == "0" || $add2 == "(0)" } {
                   return $add1
                }
                return "$add1+$add2"
               }
       "*"     {set mult1 [toexpr2 * $operand1]
                set mult2 [toexpr2 * $operand2]
                if { $mult1 == "1" || $mult1 == "(1)" } {
                   return $mult2
                }
                if { $mult2 == "1" || $mult2 == "(1)" } {
                   return $mult1
                }
                if { $mult1 == "0"   || $mult2 == "0"   ||
                     $mult1 == "(0)" || $mult2 == "(0)"    } {
                   return 0
                }
                return "$mult1*$mult2"
               }
       "/"     {return "[toexpr2 / $operand1]/[toexpr2 / $operand2]"}
       "npow"  {return "pow([toexpr $operand1],[toexpr $operand2])"}
       "atan2" {return "atan2([toexpr $operand1],[toexpr $operand2])"}
       "hypot" {return "hypot([toexpr $operand1],[toexpr $operand2])"}
       default {return "$op\([toexpr $operand1])"}
       }
    }
 }

 # toexpr2 --
 #    Reconstruct the expression from the parsed and manipulated
 #    form, add brackets if necessary given the context
 # Arguments:
 #    context        Operator context
 #    parse_tree     The parsed and manipulated expression
 # Result:
 #    String representing the expression as [expr] would take it
 #
 proc toexpr2 { context parse_tree } {
    if { [llength $parse_tree] == 1 } {
       return $parse_tree
    } else {
       set op [lindex $parse_tree 0]
       if { $op == $context || [llength $parse_tree] == 2 } {
          return [toexpr $parse_tree]
       } else {
          return "([toexpr $parse_tree])"
       }
    }
 }

 # deriv --
 #    Construct the derivative w.r.t. a given variable
 #
 # Arguments:
 #    var            Name of the variable
 #    parse_tree     The parsed and manipulated expression
 # Result:
 #    Parsed expression representing the derivative
 #
 proc deriv { var parse_tree } {
    #
    # Two cases:
    # - The parse tree consists of the expression "$var" only
    # - The parse tree is more complicated, then delegate the
    #   task to the subexpressions and assemble
    #
    if { [llength $parse_tree] == 1 } {
       if { "$parse_tree" == "\$$var" } {
          return 1
       } else {
          return 0
       }
    } else {
       foreach {op operand1 operand2} $parse_tree {break}
       switch -- $op {
       "umin"  {return [list umin [deriv $var $operand1]]}
       "-"     {return [list - [deriv $var $operand1] [deriv $var $operand2]]}
       "+"     {return [list + [deriv $var $operand1] [deriv $var $operand2]]}
       "*"     {return [list + \
                          [list * [deriv $var $operand1] $operand2 ] \
                          [list * $operand1 [deriv $var $operand2] ] ]}
       "/"     {return [list / \
                          [list - \
                             [list * [deriv $var $operand1] $operand2 ] \
                             [list * $operand1 [deriv $var $operand2] ] ] \
                          [list * $operand2 $operand2] ]}
       "npow"  {return [list * \
                          [list * $operand2 [deriv $var $operand1]] \
                          [list npow $operand1 [expr {$operand2-1}]]   ]}
       "sin"   {return [list * \
                          [deriv $var $operand1] [list cos $operand1] ]}
       "cos"   {return [list umin \
                          [list * \
                             [deriv $var $operand1] [list sin $operand1]] ]}
       "exp"   {return [list * \
                          [deriv $var $operand1] [list exp $operand1] ]}

       default {error "Derivative for '$op' not implemented"}
       }
    }
 }

 #
 # Simple test
 #
 set x {$x}
 puts "Original expression: {$x/(1+$x)}"
 puts "Reconstructed: [toexpr {/ $x {+ 1 $x}}]"
 puts "Derivative: [toexpr [deriv x {/ $x {+ 1 $x}}]]"
 puts "Original expression: {sin($x)*exp($x)}"
 puts "Reconstructed: [toexpr {* {sin $x} {exp $x}}]"
 puts "Derivative: [toexpr [deriv x {* {sin $x} {exp $x}}]]"

If you run this, the result is:

 Original expression: {$x/(1+$x)}
 Reconstructed: $x/(1+$x)
 Derivative: (((1+$x))-($x))/((1+$x)*(1+$x))
 Original expression: {sin($x)*exp($x)}
 Reconstructed: sin($x)*exp($x)
 Derivative: (cos($x)*exp($x))+(sin($x)*exp($x))

Not bad, eh? For one hour's work.

There, I have said it thrice, and what I say thrice, is true!

(Note: I just had get Lewis Carroll in there somewhere :D)


SYMDIFF is a computer algebra tool capable of taking symbolic derivatives. Using a natural syntax, it is possible to manipulate symbolic equations to aid derivation of equations for a variety of applications. Additional commands provide the means to simplify results, create common subexpressions, and order expressions for use as source code in a computer program. With its Python and Tcl interpreters, you have the ability to create algorithms to generate equations programatically.

The source code is available under the terms of the Apache License Version 2.0.