# Presentation: 228676

## mxyzptlk An automatic Differentiation Library:Recent Developments 

- Jean-Francois Ostiguy
- AD/AP
- May 16 2007

## Acknowledgments

- mxyzpltk is one among of a set of libraries that Leo Michelotti started developing in the early 1990s.
- Almost everything I know about automatic differentiation, I learned from Leo, with whom I still collaborate closely
- Thanks to J. Amunson (CD) for some useful bug reports.

## What is mxyzptlk ?

**mxyzptlk** is an Automatic Differentiation engine. CHEF, ( a framework for accelerator physics related computations) is layered on top of mxyzptlk ...
**To reiterate**: mxyzptlk does not know **anything** about beamlines, lattice functions, etc. It is purely used for the numerical computation of derivatives.
- This brief report is specifically about mxyzptlk, not about CHEF ...

## CHEF Libraries Hierarchy

## What is Automatic Differentiation ?

- In a nutshell: systematic application of rules such as Leibnitz’s rule to compute derivatives to machine precision without having to perform finite difference approximations
- Basic operations are performed on mathematical abstractions that carry information about the value of all derivatives.  In mxyzptlk, the concrete type for this abstraction is a **Jet.**
- A first order Jet carries the same information as a matrix (in fact, the Jet is more general because it includes a translation )

## Using Leibnitz’z Rule

- Problem: let h(x) = f(x) * g(x)
- Compute:
-                     d/dx [ h(x) ] at x = x0
- What is conventionally done (finite difference):
-   d/dx  [ h(x)]x0  =  h(x0+Δx) – h(x0) / Δx
-     roundoff error introduced after subtracting two nearly equal numbers
- Alternate approach (Leibnitz’s rule):
-     d/dx [h(x)]x0 = f’(x0)g(x0) + f(x0)g’(x0)
-      Assuming that f’(x0) and g’(x0) are known, h’(x0) is computed to machine precision

## Example: sqrt(x+y+z+px+py+pz)

-      Expand the above function in a Taylor Series about
-      (x=1,y=1,z=1  px=0,py=0,pz=0)
-      at order 12:
-      result has 12564
-      non-zero terms

## Performance

- This is a rather specialized issue and there is not enough time to go into details. I am planning a Division level Seminar.
- From a user standpoint, _**performance**_ is one thing (the only thing ?) that matters:
- 1 st order performance must be competitive with “standard” matrix-based codes ( e.g. beam optics)
- For high order computations: performance must be competitive with other existing specialized codes
- a few years ago, performance was not acceptable
- Incremental improvements and algorithmic changes were made;  performance steadily improved
- Late 2006- early 2007: quantum leap in performance due in part to a new internal monomial indexing scheme
- mxyzptlk performance is now _**very competitive**_ with other AD engines; **it is basically 2-3 orders of magnitudes better than > 2004**.
- Further performance improvements are still possible !

## Sample Benchmark

| Order | No of monomials
(6 variables) | CPU time (s)
P4 3 GHz |
| --- | --- | --- |
| 1 | 7 | 0.08 |
| 2 | 28 | 0.52 |
| 3 | 84 | 1.42 |
| 4 | 210 | 4.37 |
| 5 | 462 | 16.08 |
| 6 | 924 | 51.35 |
| 7 | 1716 | 146.76 |
| 8 | 3003 | 367.58 |

-    The benchmark code computes a 6D phase space map for a complete Booster synchrotron. While this is not a rigorous benchmark, it provides a useful indication of real-world performance

## CPU Time is Dominated by Jet Multiplications

## What is this good for ?

- Derivatives are the backbone of perturbation theory. Useful tool to perform perturbation analysis. Simply perform computations by declaring variables as Jet to automatically get all derivatives without adding additional code.
- Understanding/suppressing optical  aberrations
- Coupled dispersion computations (no special matrix formalism required, just propagate a “JetParticle” an look at the coefficients). Dispersion can be computed to high orders.
- Expansion of the momentum compaction to arbitrary order (e.g. in isochronous rings).
- Propagate small phase space “patches” (small phase space volume) instead of single particles. Can be Used to speed up emittance dilution computations

## Remarks and Conclusions

- mxyzpltk performance now makes higher order computations possible on a std desktop computer
- mxyzpltk performance now compares **more than favorably** with other existing alternatives.
-  mxyzptlk is native C++. This may or may not be considered to be an advantage by some.
- mxyzpltlk does not involve a special language an a special compiler
- mxyzptlk has a very reasonable memory footprint