Example: A Symbolic Differentiator
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk0">Introduction</Text-field> This section illustrates the various module concepts through a symbolic differentiator example. Since Maple provides a built-in differentiator diff, the example symbolic differentiator is named differentiate. Its (final) implementation is in the module DiffImpl (later in this chapter), which holds all the local state for the program. Much of the code for the differentiator is designed to implement either a standard rule (such as the rule that the derivative of a sum is the sum of the derivatives of the summands), or special case rules for mathematical functions such as sin and exp. The example differentiator handles only real valued functions of a single real variable. The following example shows several steps in the development of the module, from a very simple first try to the final, fully functional program. The final form of the differentiator is a good illustration of a very common Maple design pattern. This pattern arises when you have a single top-level routine that dispatches a number of subroutines to handle special cases using special purpose algorithms.
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk1">The First Attempt</Text-field> This initial example presents the differentiator as an ordinary procedure, not a module. differentiate := proc( expr, var ) local a, b; if type( expr, 'constant' ) then 0; elif expr = var then 1; elif type( expr, '`+`' ) then map( procname, args ); elif type( expr, '`^`' ) then a, b := op( expr ); if a = var and not has( b, var ) then b * a ^ ( b - 1 ); else 'procname( args )'; end if; elif type( expr, '`*`' ) then a, b := op( 1, expr ), subsop( 1 = 1, expr ); procname( a, var ) * b + a * procname( b, var ); else 'procname( args )'; end if; end proc: Trivial cases are handled first: The derivative of a constant expression is equal to 0, and the derivative of the variable with respect to which we are differentiating is equal to 1. The additivity of the derivative operator is expressed by mapping the procedure over sums, using the command map( procname, args ); LUkpcHJvY25hbWVHNiI2I0klYXJnc0dGJA== This is commonly used to map a procedure over its first argument, passing along all the remaining arguments. Only the simple case of powers of the differentiation variable is handled so far, provided also that the power is independent of the differentiation variable. The product rule for derivatives is expressed by splitting expressions of type product into two pieces: the first factor in the product, and the product of all the remaining factors. This is achieved by the double assignment of a, b := op( 1, expr ), subsop( 1 = 1, expr ); PjYkSSJhRzYiSSJiR0YlNiRJJWV4cHJHRiUiIiI= so the input expression expr is expressed as expr = a * b. The standard technique of returning unevaluated is used so that computation can proceed symbolically on expressions that the procedure is unable to differentiate. This first example is simple, but it is already able to handle polynomials with numeric coefficients. differentiate( 2 - x + x^2 + 3*x^9, x ); LCgqJEkieEc2IiIiKSIjRkYkIiIjISIiIiIi However, it fails on expressions containing calls to standard mathematical functions. differentiate( sin( x ), x ); LUkuZGlmZmVyZW50aWF0ZUc2IjYkLUkkc2luRzYkSSpwcm90ZWN0ZWRHRilJKF9zeXNsaWJHRiQ2I0kieEdGJEYs It is also unable to deal successfully with symbolic coefficients. differentiate( a*x^2 + b*x + c, x ); LCoqJi1JLmRpZmZlcmVudGlhdGVHNiI2JEklZXhwckdGJkkieEdGJiIiIkYpIiIjRioqJkYoRipGKUYqRistRiU2JEkiY0dGJkYpRipGKkYq
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk2">Adding Missing Functionality</Text-field> To add the missing functionality, add a case for expressions of type function. differentiate := proc( expr, var ) local a, b; if not has( expr, var ) then 0; elif expr = var then 1; elif type( expr, '`+`' ) then map( procname, args ); elif type( expr, '`^`' ) then a, b := op( expr ); if not has( b, var ) then b * a ^ ( b - 1 ) * procname( a, var ); else 'procname( args )'; end if elif type( expr, '`*`' ) then a, b := op( 1, expr ), subsop( 1 = 1, expr ); procname( a, var ) * b + a * procname( b, var ) elif type( expr, 'function' ) and nops( expr ) = 1 then # functions of a single variable; chain rule b := op( 0, expr ); # the name of the function a := op( 1, expr ); # the argument if b = 'sin' then cos( a ) * procname( a, var ); elif b = 'cos' then -sin( a ) * procname( a, var ); elif b = 'exp' then exp( a ) * procname( a, var ); elif b = 'ln' then ( 1 / a ) * procname( a, var ); else 'procname( args )'; end if; else 'procname( args )'; end if; end proc: This uses the chain rule to compute the derivatives of calls to known functions. differentiate( sin( x ) + cos( exp( x ) ), x ); LCYtSSRjb3NHNiRJKnByb3RlY3RlZEdGJkkoX3N5c2xpYkc2IjYjSSJ4R0YoIiIiKiYtSSRzaW5HNiRGJkYnNiMtSSRleHBHNiRGJkYnRilGK0YxRishIiI= differentiate( sin( x^2 ) + cos( x^2 ), x ); LCYqJi1JJGNvc0c2JEkqcHJvdGVjdGVkR0YnSShfc3lzbGliRzYiNiMqJEkieEdGKSIiIyIiIkYsRi5GLSomLUkkc2luRzYkRidGKEYqRi5GLEYuISIj differentiate( sin( x )^2 + cos( x )^3, x ); LCYqJi1JJHNpbkc2JEkqcHJvdGVjdGVkR0YnSShfc3lzbGliRzYiNiNJInhHRikiIiItSSRjb3NHNiRGJ0YoRipGLCIiIyomRi1GMEYkRiwhIiQ= At the same time, this has also improved the handling of expressions independent of the variable of differentiation. differentiate( a*x^2 + b*x + c, x ); LCYqJkklZXhwckc2IiIiIkkieEdGJUYmIiIjRiZGJg== This is effected by using the expression has( expr, var ) instead of the weaker test type( expr, 'constant' ). The power rule now handles more than just powers of var. differentiate( sin( x )^2, x ); LCQqJi1JJHNpbkc2JEkqcHJvdGVjdGVkR0YnSShfc3lzbGliRzYiNiNJInhHRikiIiItSSRjb3NHNiRGJ0YoRipGLCIiIw== However, adding new functions to the differentiator is tedious and error prone, and the job of handling the chain rule must be repeated for each function recognized by it.
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk3">Introducing a Function Table</Text-field> Many functions (that you need to add) and the rules used for their differentiation can be stored in a table as follows: differentiate := proc( expr, var ) local a, b, functab; functab := table(); functab[ 'sin' ] := 'cos'; functab[ 'cos' ] := x -> -sin( x ); functab[ 'exp' ] := exp; functab[ 'ln' ] := x -> 1 / x; if not has( expr, var ) then 0; elif expr = var then 1; elif type( expr, '`+`' ) then map( procname, args ); elif type( expr, '`^`' ) then a, b := op( expr ); if a = var and not has( b, var ) then b * a ^ ( b - 1 ) * procname( a, var ); else 'procname( args )'; end if elif type( expr, '`*`' ) then a, b := op( 1, expr ), subsop( 1 = 1, expr ); procname( a, var ) * b + a * procname( b, var ); elif type( expr, 'function' ) and nops( expr ) = 1 then # functions of a single variable; chain rule b := op( 0, expr ); # the name of the function a := op( 1, expr ); # the argument if assigned( functab[ b ] ) then # This is a \"known\" function functab[ b ]( a ) * procname( a, var ); else # This function is not known; return unevaluated 'procname( args )'; end if; else 'procname( args )'; end if; end proc: This not only simplifies the code used for the function case, but also makes it very easy to add new functions.
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk4">Drawbacks</Text-field> Unfortunately, this implementation has serious drawbacks. It is not extensible. The known functions are hard-coded as part of the procedure definition for differentiate. New functions cannot be added without editing this source code. A second problem relates to performance. A complete implementation would require a table of dozens or hundreds of functions. That large table would need to be created and initialized each time differentiate is invoked.
<Text-field style="Heading 2" layout="Heading 2" bookmark="bkmrk5">Encapsulation and Extensibility</Text-field> One way to fix both problems is to make the table of functions a global variable. However, using global variables can be dangerous, because they pollute the user namespace and are subject to unwanted inspection and tampering.