--- /dev/null
+[SymPy][] is a [[Python]] library for symbolic mathematics. To give
+you a feel for how it works, lets extrapolate the extremum location
+for $f(x)$ given a quadratic model:
+
+\[
+ f(x) = A x^2 + B x + C
+\]
+
+and three known values:
+
+\[
+\begin{aligned}
+ f(a) &= A a^2 + B a + C \\
+ f(b) &= A b^2 + B b + C \\
+ f(c) &= A c^2 + B c + C
+\end{aligned}
+\]
+
+Rephrase as a matrix equation:
+
+\[
+\begin{pmatrix}
+ f(a) \\
+ f(b) \\
+ f(c)
+\end{pmatrix}
+ =
+\begin{pmatrix}
+ a^2 & a & 1 \\
+ b^2 & b & 1 \\
+ c^2 & c & 1
+\end{pmatrix}
+ \cdot
+\begin{pmatrix}
+ A \\
+ B \\
+ C
+\end{pmatrix}
+\]
+
+So the solutions for $A$, $B$, and $C$ are:
+
+\[
+\begin{pmatrix}
+ A \\
+ B \\
+ C
+\end{pmatrix}
+ =
+\begin{pmatrix}
+ a^2 & a & 1 \\
+ b^2 & b & 1 \\
+ c^2 & c & 1
+\end{pmatrix}^{-1}
+ \cdot
+\begin{pmatrix}
+ f(a) \\
+ f(b) \\
+ f(c)
+\end{pmatrix}
+ =
+\begin{pmatrix}
+ \text{long} \\
+ \text{complicated} \\
+ \text{stuff}
+\end{pmatrix}
+\]
+
+Now that we've found the model parameters, we need to find the $x$
+coordinate of the extremum.
+
+\[
+ \frac{\mathrm{d}f}{\mathrm{d}x} = 2 A x + B \;,
+\]
+
+which is zero when
+
+\[
+\begin{aligned}
+ 2 A x &= -B \\
+ x &= \frac{-B}{2 A}
+\end{aligned}
+\]
+
+Here's the solution in SymPy:
+
+ >>> from sympy import Symbol, Matrix, factor, expand, pprint, preview
+ >>> a = Symbol('a')
+ >>> b = Symbol('b')
+ >>> c = Symbol('c')
+ >>> fa = Symbol('fa')
+ >>> fb = Symbol('fb')
+ >>> fc = Symbol('fc')
+ >>> M = Matrix([[a**2, a, 1], [b**2, b, 1], [c**2, c, 1]])
+ >>> F = Matrix([[fa],[fb],[fc]])
+ >>> ABC = M.inv() * F
+ >>> A = ABC[0,0]
+ >>> B = ABC[1,0]
+ >>> x = -B/(2*A)
+ >>> x = factor(expand(x))
+ >>> pprint(x)
+ 2 2 2 2 2 2
+ a *fb - a *fc - b *fa + b *fc + c *fa - c *fb
+ ---------------------------------------------
+ 2*(a*fb - a*fc - b*fa + b*fc + c*fa - c*fb)
+ >>> preview(x, viewer='pqiv')
+
+Where `pqiv` is the executable for [pqiv][], my preferred image
+viewer. With a bit of additional factoring, that is:
+
+\[
+ x = \frac{a^2 [f(b)-f(c)] + b^2 [f(c) - f(a)] + c^2 [f(a) - f(b)]}
+ {2\cdot\{a [f(b) - f(c)] + b [f(c) - f(a)] + c [f(a) - f(b)]\}}
+\]
+
+[SymPy]: http://sympy.org/
+[pqiv]: http://www.pberndt.com/Programme/Linux/pqiv/
+
+[[!tag tags/python]]
+[[!tag tags/teaching]]
+[[!tag tags/tools]]