Does Mathematica have any AD functionality or does it only support symbolic differentiation? If not, are there any packages or other third party implementations available?

Edit (J. M.)

Here is an example that might clarify how AD works. Consider the Horner method for evaluating a polynomial from its coefficients:

h[coeffs_?VectorQ, x_] := Module[{np1 = Length[coeffs], p, k},
   p = 0;
   Do[p = x p + coeffs[[k]], {k, np1, 1, -1}];

An AD method should be able to "differentiate" this procedure and automatically produce this:

hp[coeffs_?VectorQ, x_] := Module[{np1 = Length[coeffs], p, p$, k},
       p$ = 0; p = 0;
   Do[p$ = x p$ + p;
      p = x p + coeffs[[k]], {k, np1, 1, -1}];

Note the use of the product rule to differentiate the line p = x p + coeffs[[k]], as well as the ordering of the variable assignments in the "differentiated procedure".

There are more complicated examples, but traditionally, AD methods have dealt with routines that have been written procedurally. A challenge, then, is if AD methods can cope with more "idiomatic" Mathematica programs; that is, functional/rule-based algorithms. I believe that will be more difficult.

  • $\begingroup$ I don't believe there's any AD functionality built-in, but it looks to me that Mathematica is currently powerful enough that one should be able to write an AD procedure. $\endgroup$ – J. M.'s discontentment Jun 11 '15 at 18:34
  • $\begingroup$ I made a few additions; I hope you don't mind. BTW: have you seen this? $\endgroup$ – J. M.'s discontentment Jun 11 '15 at 18:57
  • 1
    $\begingroup$ Is there a reason why you want source code transformation instead of operator overloading (which is straightforward to do in MMA, I think)? $\endgroup$ – Niki Estner Jun 13 '15 at 9:23
  • $\begingroup$ Duplicate: mathematica.stackexchange.com/questions/1873/… (haven't voted to close yet since presumably we could do with a more complete answer). $\endgroup$ – Szabolcs Jul 13 '15 at 11:55

This is not a package, just an idea how such a package might work. You could define a custom value type autoDiffValue that stores a value and it's first derivative (with respect to some variable), then define arithmetic operations for that type:

autoDiffValue /: autoDiffValue[a_, adx_] + autoDiffValue[b_, bdx_] := 
 autoDiffValue[a + b, adx + bdx]
autoDiffValue /: autoDiffValue[a_, adx_] * autoDiffValue[b_, bdx_] := 
 autoDiffValue[a*b, adx*b + a*bdx]
autoDiffValue /: autoDiffValue[a_, adx_] + b_ := 
 autoDiffValue[a + b, adx]
autoDiffValue /: autoDiffValue[a_, adx_]*b_ := 
 autoDiffValue[a*b, adx*b]
autoDiffValue /: autoDiffValue[a_, adx_]^b_ := 
 autoDiffValue[a^b, a^(b - 1)*adx*b]

Using that and your function h, you can then call:

h[{1, 2, 3}, autoDiffValue[5, 1]]

and get

autoDiffValue[86, 32]

i.e. the value and the first derivative of h[{1, 2, 3}, x] at x=5

You could automate generating derivative definitions:

defineFunc[fn_] := 
 autoDiffValue /: fn[autoDiffValue[a_, adx_]] = 
  autoDiffValue[fn[a], D[fn[a + x*adx], x] /. x -> 0]
defineFunc /@ {Sin, Cos, Exp}; 

(it's probably possible to extend defineFunc so it can create definitions for Plus, Times, ... too)

This is more or less the same technique you would use to define AD in e.g. C++ using operator overloading. It should work in principle, although a lot of work would still be needed to implement multivariable and higher derivatives and to support vector and matrix arguments (so calculations can be done with packed arrays).

| improve this answer | |

I am not sure I understood the advantages of AD from your example. It would be nice to have a function that you can't differentiate properly or fast with Derivative.

This is not an answer, but too long for a comment. I suppose the automatic accuracy tracking is some sort of very limited built-in automatic differentiation, right? At least for combinations of real monotonically increasing functions.

makeDual[n_] := SetAccuracy[n, 1]
getDual[n_] := {N@n, 10^(1 - Accuracy[n])}

f[x_] := If[x > 0, Log[Sin[x]] + Tan[x]^3 - Sin[Cos[Sin[x]]]^17, x - 3]
r = 0.234;
{f[r], f'[r]} == getDual@f[makeDual[r]]
(* True *)
| improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.