From 1cfc5217717767e745bf4a742c3ea5fd0dd0ca52 Mon Sep 17 00:00:00 2001 From: "W. Trevor King" Date: Wed, 27 Feb 2013 07:51:53 -0500 Subject: [PATCH] posts:sympy: Add post demonstrating usage --- posts/SymPy.mdwn_itex | 121 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 121 insertions(+) create mode 100644 posts/SymPy.mdwn_itex diff --git a/posts/SymPy.mdwn_itex b/posts/SymPy.mdwn_itex new file mode 100644 index 0000000..e8a7cde --- /dev/null +++ b/posts/SymPy.mdwn_itex @@ -0,0 +1,121 @@ +[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]] -- 2.26.2