Symbolic Differentiation Tutorial

In the following, we like to program a simple machinery that allows us to differentiate univariate functions.

To understand the programs here, you should first have read and understood ProgrammingSPAD.

The aim of this page is rather to give a short account on how one can quickly achieve quite some sophisticated program with a few lines of SPAD code. The code deliberately does only build on a very small subset of the FriCAS library. In fact, we only use the category SetCategory and the domains Integer, Fraction, Complex, String, Symbol, List, Product, Kernel, OutputForm, Expression. Although the Expression domain already comes with differentiation, we are not going to use it. We could have done without Expression at the cost of implementing our own expression domain rule for identifying e1+2 with e2+e1, etc. However, since the purpose of this page is to show programming in SPAD, rather than entering the complicated area of low-level expression manipulation, we decided against it.

The first file mex.spad was actually used as the basis for a little project to teach some children (age 16-18) how Symbolic Differentiation basically works. They were supposed to "invent" the dex.spad file.

The file mex.spad is not really dealing with symbolic differentiation. It rather deals with very basic things like defining what a monoid, a ring, or a field is.

The most interesting thing is that it defines expressions that can be decomposed into a top-level symbol and its corresponding arguments (which might be complicated expressions themselves).

With the help of this MExpression domain, it is realatively easy to program symbolic differentiation by simply looking at the top-level symbol, and applying the corresponding rule for it and differentiating the arguments according to the chain rule.

(1) -> <spad>
FILE ==> "mex.spad"
-------------------------------------------------------------------
--
-- MEX
-- Copyright (C) 2012,  Ralf Hemmecke <ralf@hemmecke.de>
--
-------------------------------------------------------------------
-- This program is free software: you can redistribute it and/or modify
-- it under the terms of the GNU General Public License as published by
-- the Free Software Foundation, either version 3 of the License, or
-- (at your option) any later version.
--
-- This program is distributed in the hope that it will be useful,
-- but WITHOUT ANY WARRANTY; without even the implied warranty of
-- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-- GNU General Public License for more details.
--
-- You should have received a copy of the GNU General Public License
-- along with this program.  If not, see <http://www.gnu.org/licenses/>.
-------------------------------------------------------------------
-- This file contains basic definitions for a simple expression domain.
-------------------------------------------------------------------

-- These two macros are necessary to distinguish between Rep and %.
rep x ==> (x@%) pretend Rep
per x ==> (x@Rep) pretend %
Z ==> Integer
Q ==> Fraction Z
C ==> Complex Q
X ==> Expression C
D ==> Product(String, List %)

-------------------------------------------------------------------

)abbrev category MABMON MAbelianAdditiveMonoid
++ MAbelianAdditiveMonoid is the category of abelian monoids that
++ are written in an additive way.
MAbelianAdditiveMonoid: Category == SetCategory with
    0: %
        ++ 0 is the neutral element w.r.t. +.
    _+: (%, %) -> %
        ++ + is a commutative and associative operation.

-------------------------------------------------------------------

)abbrev category MABAGRP MAbelianAdditiveGroup
++ MAbelianAdditiveGroup is the category of abelian groups that
++ are written in an additive way.
MAbelianAdditiveGroup: Category == MAbelianAdditiveMonoid with
    _-: % -> %
        ++ x + (-x) == 0
    _-: (%, %) -> %
        ++ x - y = z if and only if x = y + z.
  add
    ((x: %) - (y: %)): % == x + (- y)

-------------------------------------------------------------------

)abbrev category MMMON MMultiplicativeMonoid
++ MMultiplicativeMonoid is the category of monoids that
++ are written in a multiplicative way.
MMultiplicativeMonoid: Category == SetCategory with
    1: %
        ++ 1 is the neutral element w.r.t. *.
    _*: (%, %) -> %
        ++ * is an associative operation.

-------------------------------------------------------------------

)abbrev category MABMMON MAbelianMultiplicativeMonoid
++ MAbelianMultiplicativeMonoid is the category of abelian monoids that
++ are written in a multiplicative way.
MAbelianMultiplicativeMonoid: Category == MMultiplicativeMonoid

-------------------------------------------------------------------

)abbrev category MRING MRing
++ MRing is the category of rings with unity.
MRing: Category == Join(MAbelianAdditiveGroup, MMultiplicativeMonoid)

-------------------------------------------------------------------

)abbrev category MCRING MCommutativeRing
++ MCommutativeRing is the category of commutative rings with unity.
MCommutativeRing: Category == Join(MRing, MAbelianMultiplicativeMonoid)

-------------------------------------------------------------------

)abbrev category MINTDOM MIntegralDomain
++ MIntegralDomain is the category of commutative rings with unity
++ that have no zero-divisors.
MIntegralDomain: Category == MCommutativeRing

-------------------------------------------------------------------

)abbrev category MFIELD MField
++ MField is the category of fields.
MField: Category == MIntegralDomain with
    inv: % -> %
        ++ x * inv(x) = 1
    _/: (%, %) -> %
        ++ Only defined if the second argument is non-zero.
        ++ x/y = z if and only if x = z*y
    _^: (%, Integer) -> %
        ++ x^n is an abbreviation for x*x*...*x (n-fold product of x) if n>0.
        ++ x^0 := 1 for any x (even for x=0),
        ++ x^n := 1/(x^(-n)) if n<0.
  add
    inv(x: %): % == 1/x

-------------------------------------------------------------------

)abbrev category MEXCAT MExpressionCategory
++ MExpressionCategory is the category of formal expressions domains
++ that are closed under field operations and taking sin, cos, exp,
++ log, and sqrt.
MExpressionCategory: Category == MField with
    coerce: % -> X
        ++ coerce(z) transforms z into an equivalent value in
        ++ Expression(Complex Fraction Integer).
    coerce: Z -> %
        ++ coerce(z) transforms the integer z into an element of this type.
    coerce: Q -> %
        ++ coerce(z) transforms the rational number z into an element
        ++ of this type.
    coerce: Complex Z -> %
        ++ coerce(z) transforms the complex number z into an element
        ++ of this type.
    coerce: Complex Q -> %
        ++ coerce(z) transforms the complex number z into an element
        ++ of this type.
    coerce: Symbol -> %
        ++ coerce(z) returns z as a variable of this type.
    variable?: % -> Boolean
        ++ variable?(z) returns true if z is a variable and false otherwise.
    complex?: % -> Boolean
        ++ complex?(z) returns true if z can be represented as a
        ++ complex number and false otherwise
    _=: (%, %) -> Boolean
        ++ y=z tests whether two expressions are (after mild simplification)
        ++ syntactically equal.
    decompose: % -> Product(String, List %)
        ++ decompose(z) return ["complex", [z]] if complex?(z) is true.
        ++ It returns ["variable", [z]] if variable?(z) is true.
        ++ In all other cases z is an expression of the form foo(a1,...,an)
        ++ and ["foo", [a1, ..., an]] will be returned.
        ++ Note that for infix expression like a1+a2+...+an, this functions
        ++ returns ["+",[a1,...,an]].
    _^: (%, %) -> %
        ++ y^z is a formal exponentiation expression.
    sqrt: % -> %
        ++ sqrt(z) is an expresssion that denotes square root of z.
    sin: % -> %
        ++ sin(z) is an expression that denotes the sine of z.
    cos: % -> %
        ++ cos(z) is an expression that denotes the cosine of z.
    exp: % -> %
        ++ exp(z) is an expression that denotes the e^z where e is
        ++ Euler's constant.
    log: % -> %
        ++ log(z) is an expression that denotes the natural logarithm of z.
    tan: % -> %
        ++ tan(z) is an expression that denotes the tangent of z.
    cot: % -> %
        ++ cot(z) is an expression that denotes the cotangent of z.
    asin: % -> %
        ++ asin(z) is an expression that denotes the arc-sine of z.
    acos: % -> %
        ++ acos(z) is an expression that denotes the arc-cosine of z.
    atan: % -> %
        ++ atan(z) is an expression that denotes the arc-tangent of z.
    acot: % -> %
        ++ acot(z) is an expression that denotes the arc-cotangent of z.

-------------------------------------------------------------------

)abbrev domain MEX MExpression

MExpression: MExpressionCategory == add
    Rep ==> X

    -- exported functions
    coerce(a: %): X == a pretend X

    coerce(a: %): OutputForm == coerce rep a
    coerce(z: Z): % == per(z::X)
    coerce(q: Q): % == per(q::X)
    coerce(c: C): % == per(c::X)
    coerce(c: Complex Z): % == complex(real(c)::Q, imag(c)::Q)$C::%
    coerce(s: Symbol): % == per(s::X)
    variable?(x: %): Boolean == empty? second decompose x
    complex?(a: %): Boolean == number? rep a
    0: % == (0$Z)::%
    1: % == (1$Z)::%
    ((a: %) = (b: %)): Boolean == rep a = rep b

    _-(a: %): % == per(- rep a)
    ((a: %) + (b: %)): % == per(rep a + rep b)
    ((a: %) * (b: %)): % == per(rep a * rep b)
    ((a: %) / (b: %)): % == per(rep a / rep b)
    ((a: %) ^ (z: Z)): % == per(rep a ^ z)
    ((a: %) ^ (b: %)): % == per(rep a ^ rep b)
    sqrt(a: %): % == per sqrt rep a
    sin(a: %): % == per sin rep a
    cos(a: %): % == per cos rep a
    exp(a: %): % == per exp rep a
    log(a: %): % == per log rep a
    tan(a: %): % == per tan rep a
    cot(a: %): % == per cot rep a
    asin(a: %): % == per asin rep a
    acos(a: %): % == per acos rep a
    atan(a: %): % == per atan rep a
    acot(a: %): % == per acot rep a

    decompose(a: %): D ==
        complex? a => ["complex", [a]]$D
        ra := rep a
        il := isPlus ra
        il case List(X) => ["+", (il::List(X)) pretend List(%)]$D
        il := isTimes ra
        il case List(X) => ["*", (il::List(X)) pretend List(%)]$D
        ip := isPower ra
        ip case Record(val: X, exponent: Z) and ip.exponent ~= 1 =>
            ["^", [ip.val, (ip.exponent)::X] pretend List(%)]$D
        ks: List(Kernel X) := kernels ra
        #(ks)~=1 =>
            print("decompose: more than one kernel detected"::OutputForm)
            print(ra::OutputForm)
            print(ks::OutputForm)
            error "decompose error"
        k: Kernel(X) := first ks
        [string name k, [per z for z in argument k]]$D</spad>

Compiling FriCAS source code from file 
      /var/lib/zope2.10/instance/axiom-wiki/var/LatexWiki/8483626945223991797-25px001.spad
      using old system compiler.
   MABMON abbreviates category MAbelianAdditiveMonoid 
------------------------------------------------------------------------
   initializing NRLIB MABMON for MAbelianAdditiveMonoid 
   compiling into NRLIB MABMON 

;;;     ***       |MAbelianAdditiveMonoid| REDEFINED
Time: 0.01 SEC.

   finalizing NRLIB MABMON 
   Processing MAbelianAdditiveMonoid for Browser database:
--------constructor---------
--------((Zero) (%) constant)---------
--->-->MAbelianAdditiveMonoid(((Zero) (%) constant)): Improper first word in comments: 
"0 is the neutral element \\spad{w}.\\spad{r}.\\spad{t}. +."
--------(+ (% % %))---------
; compiling file "/var/aw/var/LatexWiki/MABMON.NRLIB/MABMON.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MABMON.NRLIB/MABMON.fasl written
; compilation finished in 0:00:00.006
------------------------------------------------------------------------
   MAbelianAdditiveMonoid is now explicitly exposed in frame initial 
   MAbelianAdditiveMonoid will be automatically loaded when needed from
      /var/aw/var/LatexWiki/MABMON.NRLIB/MABMON

   MABAGRP abbreviates category MAbelianAdditiveGroup 
------------------------------------------------------------------------
   initializing NRLIB MABAGRP for MAbelianAdditiveGroup 
   compiling into NRLIB MABAGRP 

;;;     ***       |MAbelianAdditiveGroup| REDEFINED
Time: 0 SEC.

   MABAGRP- abbreviates domain MAbelianAdditiveGroup& 
------------------------------------------------------------------------
   initializing NRLIB MABAGRP- for MAbelianAdditiveGroup& 
   compiling into NRLIB MABAGRP- 
   compiling exported - : (S,S) -> S
Time: 0 SEC.

(time taken in buildFunctor:  0)

;;;     ***       |MAbelianAdditiveGroup&| REDEFINED
Time: 0 SEC.

   Cumulative Statistics for Constructor MAbelianAdditiveGroup&
      Time: 0 seconds

   finalizing NRLIB MABAGRP- 
   Processing MAbelianAdditiveGroup& for Browser database:
--------constructor---------
--------(- (% %))---------
--->-->MAbelianAdditiveGroup&((- (% %))): Improper initial operator in comments: +
"\\spad{x} + (\\spad{-x}) \\spad{==} 0"
--------(- (% % %))---------
; compiling file "/var/aw/var/LatexWiki/MABAGRP-.NRLIB/MABAGRP-.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MABAGRP-.NRLIB/MABAGRP-.fasl written
; compilation finished in 0:00:00.011
------------------------------------------------------------------------
   MAbelianAdditiveGroup& is now explicitly exposed in frame initial 
   MAbelianAdditiveGroup& will be automatically loaded when needed from
      /var/aw/var/LatexWiki/MABAGRP-.NRLIB/MABAGRP-
   finalizing NRLIB MABAGRP 
   Processing MAbelianAdditiveGroup for Browser database:
--------constructor---------
--------(- (% %))---------
--->-->MAbelianAdditiveGroup((- (% %))): Improper initial operator in comments: +
"\\spad{x} + (\\spad{-x}) \\spad{==} 0"
--------(- (% % %))---------
; compiling file "/var/aw/var/LatexWiki/MABAGRP.NRLIB/MABAGRP.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MABAGRP.NRLIB/MABAGRP.fasl written
; compilation finished in 0:00:00.003
------------------------------------------------------------------------
   MAbelianAdditiveGroup is now explicitly exposed in frame initial 
   MAbelianAdditiveGroup will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MABAGRP.NRLIB/MABAGRP

   MMMON abbreviates category MMultiplicativeMonoid 
------------------------------------------------------------------------
   initializing NRLIB MMMON for MMultiplicativeMonoid 
   compiling into NRLIB MMMON 

;;;     ***       |MMultiplicativeMonoid| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MMMON 
   Processing MMultiplicativeMonoid for Browser database:
--------constructor---------
--------((One) (%) constant)---------
--->-->MMultiplicativeMonoid(((One) (%) constant)): Improper first word in comments: 
"1 is the neutral element \\spad{w}.\\spad{r}.\\spad{t}. *."
--------(* (% % %))---------
; compiling file "/var/aw/var/LatexWiki/MMMON.NRLIB/MMMON.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MMMON.NRLIB/MMMON.fasl written
; compilation finished in 0:00:00.002
------------------------------------------------------------------------
   MMultiplicativeMonoid is now explicitly exposed in frame initial 
   MMultiplicativeMonoid will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MMMON.NRLIB/MMMON

   MABMMON abbreviates category MAbelianMultiplicativeMonoid 
------------------------------------------------------------------------
   initializing NRLIB MABMMON for MAbelianMultiplicativeMonoid 
   compiling into NRLIB MABMMON 

;;;     ***       |MAbelianMultiplicativeMonoid| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MABMMON 
   Processing MAbelianMultiplicativeMonoid for Browser database:
--------constructor---------
; compiling file "/var/aw/var/LatexWiki/MABMMON.NRLIB/MABMMON.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MABMMON.NRLIB/MABMMON.fasl written
; compilation finished in 0:00:00.003
------------------------------------------------------------------------
   MAbelianMultiplicativeMonoid is now explicitly exposed in frame 
      initial 
   MAbelianMultiplicativeMonoid will be automatically loaded when 
      needed from /var/aw/var/LatexWiki/MABMMON.NRLIB/MABMMON

   MRING abbreviates category MRing 
------------------------------------------------------------------------
   initializing NRLIB MRING for MRing 
   compiling into NRLIB MRING 

;;;     ***       |MRing| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MRING 
   Processing MRing for Browser database:
--------constructor---------
; compiling file "/var/aw/var/LatexWiki/MRING.NRLIB/MRING.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MRING.NRLIB/MRING.fasl written
; compilation finished in 0:00:00.002
------------------------------------------------------------------------
   MRing is now explicitly exposed in frame initial 
   MRing will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MRING.NRLIB/MRING

   MCRING abbreviates category MCommutativeRing 
------------------------------------------------------------------------
   initializing NRLIB MCRING for MCommutativeRing 
   compiling into NRLIB MCRING 

;;;     ***       |MCommutativeRing| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MCRING 
   Processing MCommutativeRing for Browser database:
--------constructor---------
; compiling file "/var/aw/var/LatexWiki/MCRING.NRLIB/MCRING.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MCRING.NRLIB/MCRING.fasl written
; compilation finished in 0:00:00.002
------------------------------------------------------------------------
   MCommutativeRing is now explicitly exposed in frame initial 
   MCommutativeRing will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MCRING.NRLIB/MCRING

   MINTDOM abbreviates category MIntegralDomain 
------------------------------------------------------------------------
   initializing NRLIB MINTDOM for MIntegralDomain 
   compiling into NRLIB MINTDOM 

;;;     ***       |MIntegralDomain| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MINTDOM 
   Processing MIntegralDomain for Browser database:
--------constructor---------
; compiling file "/var/aw/var/LatexWiki/MINTDOM.NRLIB/MINTDOM.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MINTDOM.NRLIB/MINTDOM.fasl written
; compilation finished in 0:00:00.003
------------------------------------------------------------------------
   MIntegralDomain is now explicitly exposed in frame initial 
   MIntegralDomain will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MINTDOM.NRLIB/MINTDOM

   MFIELD abbreviates category MField 
------------------------------------------------------------------------
   initializing NRLIB MFIELD for MField 
   compiling into NRLIB MFIELD 

;;;     ***       |MField| REDEFINED
Time: 0 SEC.

   MFIELD- abbreviates domain MField& 
------------------------------------------------------------------------
   initializing NRLIB MFIELD- for MField& 
   compiling into NRLIB MFIELD- 
   compiling exported inv : S -> S
Time: 0.01 SEC.

(time taken in buildFunctor:  0)

;;;     ***       |MField&| REDEFINED
Time: 0 SEC.

   Cumulative Statistics for Constructor MField&
      Time: 0.01 seconds

   finalizing NRLIB MFIELD- 
   Processing MField& for Browser database:
--------constructor---------
--------(inv (% %))---------
--->-->MField&((inv (% %))): Improper initial operator in comments: *
"\\spad{x} * inv(\\spad{x}) = 1"
--------(/ (% % %))---------
--->-->MField&((/ (% % %))): Improper first word in comments: Only
"Only defined if the second argument is non-zero. x/y = \\spad{z} if and only if \\spad{x} = z*y"
--------(^ (% % (Integer)))---------
; compiling file "/var/aw/var/LatexWiki/MFIELD-.NRLIB/MFIELD-.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MFIELD-.NRLIB/MFIELD-.fasl written
; compilation finished in 0:00:00.008
------------------------------------------------------------------------
   MField& is now explicitly exposed in frame initial 
   MField& will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MFIELD-.NRLIB/MFIELD-
   finalizing NRLIB MFIELD 
   Processing MField for Browser database:
--------constructor---------
--------(inv (% %))---------
--->-->MField((inv (% %))): Improper initial operator in comments: *
"\\spad{x} * inv(\\spad{x}) = 1"
--------(/ (% % %))---------
--->-->MField((/ (% % %))): Improper first word in comments: Only
"Only defined if the second argument is non-zero. x/y = \\spad{z} if and only if \\spad{x} = z*y"
--------(^ (% % (Integer)))---------
; compiling file "/var/aw/var/LatexWiki/MFIELD.NRLIB/MFIELD.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MFIELD.NRLIB/MFIELD.fasl written
; compilation finished in 0:00:00.002
------------------------------------------------------------------------
   MField is now explicitly exposed in frame initial 
   MField will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MFIELD.NRLIB/MFIELD

   MEXCAT abbreviates category MExpressionCategory 
------------------------------------------------------------------------
   initializing NRLIB MEXCAT for MExpressionCategory 
   compiling into NRLIB MEXCAT 

;;;     ***       |MExpressionCategory| REDEFINED
Time: 0 SEC.

   finalizing NRLIB MEXCAT 
   Processing MExpressionCategory for Browser database:
--------constructor---------
--------(coerce ((Expression (Complex (Fraction (Integer)))) %))---------
--------(coerce (% (Integer)))---------
--------(coerce (% (Fraction (Integer))))---------
--------(coerce (% (Complex (Integer))))---------
--------(coerce (% (Complex (Fraction (Integer)))))---------
--------(coerce (% (Symbol)))---------
--------(variable? ((Boolean) %))---------
--------(complex? ((Boolean) %))---------
--------(= ((Boolean) % %))---------
--------(decompose ((Product (String) (List %)) %))---------
--------(^ (% % %))---------
--------(sqrt (% %))---------
--------(sin (% %))---------
--------(cos (% %))---------
--------(exp (% %))---------
--------(log (% %))---------
--------(tan (% %))---------
--------(cot (% %))---------
--------(asin (% %))---------
--------(acos (% %))---------
--------(atan (% %))---------
--------(acot (% %))---------
; compiling file "/var/aw/var/LatexWiki/MEXCAT.NRLIB/MEXCAT.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MEXCAT.NRLIB/MEXCAT.fasl written
; compilation finished in 0:00:00.003
------------------------------------------------------------------------
   MExpressionCategory is now explicitly exposed in frame initial 
   MExpressionCategory will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MEXCAT.NRLIB/MEXCAT

   MEX abbreviates domain MExpression 
------------------------------------------------------------------------
   initializing NRLIB MEX for MExpression 
   compiling into NRLIB MEX 
   processing macro definition Rep ==> Expression Complex Fraction Integer 
   compiling exported coerce : $ -> Expression Complex Fraction Integer
      MEX;coerce;$E;1 is replaced by a 
Time: 0.06 SEC.

   compiling exported coerce : $ -> OutputForm
Time: 0 SEC.

   compiling exported coerce : Integer -> $
Time: 0 SEC.

   compiling exported coerce : Fraction Integer -> $
Time: 0 SEC.

   compiling exported coerce : Complex Fraction Integer -> $
Time: 0.01 SEC.

   compiling exported coerce : Complex Integer -> $
Time: 0 SEC.

   compiling exported coerce : Symbol -> $
Time: 0 SEC.

   compiling exported variable? : $ -> Boolean
Time: 0.02 SEC.

   compiling exported complex? : $ -> Boolean
Time: 0.01 SEC.

   compiling exported Zero : () -> $
Time: 0 SEC.

   compiling exported One : () -> $
Time: 0 SEC.

   compiling exported = : ($,$) -> Boolean
Time: 0 SEC.

   compiling exported - : $ -> $
Time: 0 SEC.

   compiling exported + : ($,$) -> $
Time: 0 SEC.

   compiling exported * : ($,$) -> $
Time: 0 SEC.

   compiling exported / : ($,$) -> $
Time: 0.01 SEC.

   compiling exported ^ : ($,Integer) -> $
Time: 0 SEC.

   compiling exported ^ : ($,$) -> $
Time: 0.01 SEC.

   compiling exported sqrt : $ -> $
Time: 0 SEC.

   compiling exported sin : $ -> $
Time: 0.01 SEC.

   compiling exported cos : $ -> $
Time: 0 SEC.

   compiling exported exp : $ -> $
Time: 0.01 SEC.

   compiling exported log : $ -> $
Time: 0 SEC.

   compiling exported tan : $ -> $
Time: 0 SEC.

   compiling exported cot : $ -> $
Time: 0 SEC.

   compiling exported asin : $ -> $
Time: 0 SEC.

   compiling exported acos : $ -> $
Time: 0.01 SEC.

   compiling exported atan : $ -> $
Time: 0 SEC.

   compiling exported acot : $ -> $
Time: 0 SEC.

   compiling exported decompose : $ -> Product(String,List $)
Time: 0.06 SEC.

(time taken in buildFunctor:  0)

;;;     ***       |MExpression| REDEFINED

;;;     ***       |MExpression| REDEFINED
Time: 0.01 SEC.

   Warnings: 
      [1] variable?:  $$ has no value
      [2] decompose:  exponent has no value
      [3] decompose:  val has no value

   Cumulative Statistics for Constructor MExpression
      Time: 0.22 seconds

   finalizing NRLIB MEX 
   Processing MExpression for Browser database:
--->-->MExpression(): Missing Description
; compiling file "/var/aw/var/LatexWiki/MEX.NRLIB/MEX.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/MEX.NRLIB/MEX.fasl written
; compilation finished in 0:00:00.075
------------------------------------------------------------------------
   MExpression is now explicitly exposed in frame initial 
   MExpression will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/MEX.NRLIB/MEX

FILE ==> "dex.spad"
-------------------------------------------------------------------
--
-- DEX
-- Copyright (C) 2012,  Ralf Hemmecke <ralf@hemmecke.de>
--
-------------------------------------------------------------------
-- This program is free software: you can redistribute it and/or modify
-- it under the terms of the GNU General Public License as published by
-- the Free Software Foundation, either version 3 of the License, or
-- (at your option) any later version.
--
-- This program is distributed in the hope that it will be useful,
-- but WITHOUT ANY WARRANTY; without even the implied warranty of
-- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-- GNU General Public License for more details.
--
-- You should have received a copy of the GNU General Public License
-- along with this program.  If not, see <http://www.gnu.org/licenses/>.
-------------------------------------------------------------------
-- This file extends a basic expression domain by differentation.
-------------------------------------------------------------------
)abbrev domain DEX DExpression

D ==> Product(String, List %)

DExpression: MExpressionCategory with
    differentiate: (%, %) -> %
  == MExpression add
    f a ==> a::OutputForm

    -- auxiliary functions
    plusD(l: List %, x: %): % ==
        r: % := 0
        for z in l repeat r := r +  differentiate(z, x)
        r

    timesD(l: List %, x: %): % ==
        s: % := 0
        len := #(l)
        for i in 1..len repeat
            p: % := 1
            for z in l for j in 1..len repeat
                p := p*(if i=j then differentiate(z, x) else z)
            s := s + p
        s

    powerD(l: List %, x: %): % ==
        #(l) ~= 2 => error "powerD: not exactly 2 arguments"
        z: % := first l
        e: % := first rest l
        not complex? e => error "powerD: second argument is not complex"
        e*(z^(e-1))*differentiate(z, x)

    powerXD(l: List %, x: %): % ==
        #(l) ~= 2 => error "powerXD: not exactly 2 arguments"
        z: % := first l
        e: % := first rest l
        complex? z => log z * z^e * differentiate(e, x)
        differentiate(exp(log z * e), x)

    chainD(l: List %, x: %, d: % -> %): % ==
        #(l) ~= 1 => error "chainD: not exactly one argument"
        z := first l
        d(z)*differentiate(z, x)

    nthRootD(l: List %, x: %): % ==
        #(l) ~= 2 => error "nthRootD: not exactly two arguments"
        z := first l
        e := inv first rest l
        not complex? e => error "nthRootD: second argument is not complex"
        e*(z^(e-1))*differentiate(z, x)

    -- exported functions
    differentiate(z: %, x: %): % ==
        not variable? x => error "second argument must be a variable"
        complex? z => 0
        z=x => 1
        variable? z => 0
        d: D := decompose z
        o: String := first d
        a: List % := second d
        o = "+" => plusD(a, x)
        o = "*" => timesD(a, x)
        o = "^" => powerD(a, x)
        o = "%power" => powerXD(a, x)
        o = "nthRoot" => nthRootD(a, x)
        o = "sin" => chainD(a, x, cos)
        o = "cos" => chainD(a, x, t +-> - sin t)
        o = "exp" => chainD(a, x, exp)
        o = "log" => chainD(a, x, t +-> 1/t)
        o = "tan" => chainD(a, x, t +-> tan(t)^2+1)
        o = "cot" => chainD(a, x, t +-> -cot(t)^2-1)
        o = "asin" => chainD(a, x, t +-> 1/sqrt(1-t^2))
        o = "acos" => chainD(a, x, t +-> -1/sqrt(1-t^2))
        o = "atan" => chainD(a, x, t +-> 1/(1+t^2))
        o = "acot" => chainD(a, x, t +-> -1/(1+t^2))
        print("differentiate expression: "::OutputForm)
        print(z::OutputForm)
        print(((z=x)@Boolean)::OutputForm)
        print(d::OutputForm)
        print(o::OutputForm)
        print(a::OutputForm)
        error "Cannot differentiate"

   Compiling FriCAS source code from file 
      /var/lib/zope2.10/instance/axiom-wiki/var/LatexWiki/6808451264496032819-25px002.spad
      using old system compiler.
   DEX abbreviates domain DExpression 
------------------------------------------------------------------------
   initializing NRLIB DEX for DExpression 
   compiling into NRLIB DEX 
   processing macro definition f a ==> ::(a,OutputForm) 
   compiling local plusD : (List $,$) -> $
Time: 0 SEC.

   compiling local timesD : (List $,$) -> $
Time: 0.01 SEC.

   compiling local powerD : (List $,$) -> $
Time: 0 SEC.

   compiling local powerXD : (List $,$) -> $
Time: 0.01 SEC.

   compiling local chainD : (List $,$,$ -> $) -> $
Time: 0 SEC.

   compiling local nthRootD : (List $,$) -> $
Time: 0 SEC.

   compiling exported differentiate : ($,$) -> $
Time: 0.01 SEC.

(time taken in buildFunctor:  0)

;;;     ***       |DExpression| REDEFINED

;;;     ***       |DExpression| REDEFINED
Time: 0.03 SEC.

   Cumulative Statistics for Constructor DExpression
      Time: 0.06 seconds

   finalizing NRLIB DEX 
   Processing DExpression for Browser database:
--->-->DExpression(constructor): Not documented!!!!
--->-->DExpression((differentiate (% % %))): Not documented!!!!
--->-->DExpression(): Missing Description
; compiling file "/var/aw/var/LatexWiki/DEX.NRLIB/DEX.lsp" (written 21 FEB 2023 12:47:03 PM):

; /var/aw/var/LatexWiki/DEX.NRLIB/DEX.fasl written
; compilation finished in 0:00:00.162
------------------------------------------------------------------------
   DExpression is now explicitly exposed in frame initial 
   DExpression will be automatically loaded when needed from 
      /var/aw/var/LatexWiki/DEX.NRLIB/DEX

After having defined an expression domain and having equipped it with a differentiate function, we can now use it.

Let's start with some abbreviations.

DEX ==> DExpression

x := 'x :: DEX

(1)

d(e) ==> differentiate(e, x)$DEX

Now we can differentiate any polynomial formed with the x from above.

d x

(2)

d(x^2)

(3)

d(5*(x^2-2)^2)

(4)

Of course, we can also differentiate all the functions that are included in the code of differentiate.

d sin(x)

(5)

d acos(x)

(6)

d(x^(-1/3))

(7)

Clearly, our code was designed with the chain rule, so we can also differentiate more complicated composed functions.

d sin acot exp(x^(1/3))

(8)

Comments on the above program

• Comments in SPAD are introduced with double minus -- and reach till the end of the line.
• SPAD knows docstrings. They are started with double plus ++ and reach until the end of line. Docstrings must be attached to actual program code. As can be seen from above, one can document whole categories and domains, as well as function signatures. These docstrings are shown in the HyperDoc system.
• Although MAbelianAdditiveGroup is a category, it has an implementation part (introduced with the keyword add). Function implementations that go in such a place are generic in the sense that they cannot refer to any specific representation of a corresponding domain.
• There is a convention is SPAD that all functions that decide something (like variable?), i.e., return a Boolean value, should end with a question mark. The question mark is part of the name.
• The form

        B => T
        E1
        E2
  

  can be read as

        if B then
                  T
             else
                  E1
                  E2
  

• A form like

        for z in l for j in 1..len repeat foo(z, j)
  

  runs in parallel through each element z of the list l and from j=1 until j=len. It will call

        f(l.1, 1)
        f(l.2, 3)
        f(l.3, 3)
        ...
  

  The function foo(z,j) is called at most len times, but if the list has fewer than len elements, then there will be only as many invocations as the list has elements.

• The function # is defined for lists in the FriCAS library. It returns the length of a list.
• The error function prints the string that it gets as an argument to the console and then aborts the program. There is not yet an exception mechanism is SPAD.
• In SPAD one can write f n or f.n for f(n). The difference is in how multiple compositions associate. We have:
  • x.y.z means (x.y)(z) which is the same as (x(y))(z), and
  • x y z means x(y(z)).