Demo of Python in Jupyter Notebook¶
In [1]:
%%html
<iframe width="560" height="315" src="https://www.youtube.com/embed/HW29067qVWk" frameborder="0" allow="autoplay; encrypted-media" allowfullscreen></iframe>
In [2]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
N = 50
x = np.random.rand(N)
y = np.random.rand(N)
colors = np.random.rand(N)
area = np.pi * (15 * np.random.rand(N))**2
plt.scatter(x, y, s=area, c=colors, alpha=0.5)
plt.show()
In [3]:
import pandas as pd
import numpy as np
df = pd.DataFrame(np.random.randn(10,5))
df
Out[3]:
| 0 | 1 | 2 | 3 | 4 | |
|---|---|---|---|---|---|
| 0 | -0.553045 | 1.058916 | -0.963079 | -1.092180 | -0.486151 |
| 1 | 0.592896 | 0.007607 | -0.798579 | 0.338088 | 0.371761 |
| 2 | -1.606119 | -1.177547 | -1.651392 | 1.448293 | 1.116540 |
| 3 | 0.794082 | 1.083283 | 0.078592 | -1.417528 | -1.014568 |
| 4 | -0.576371 | -0.305854 | -1.712976 | 0.956424 | 0.086864 |
| 5 | -0.181922 | -0.243525 | 0.619816 | -0.040076 | 1.660956 |
| 6 | -0.331784 | 0.463404 | -1.211695 | -0.150694 | 0.332713 |
| 7 | -0.083598 | 1.024621 | -1.420418 | -0.245272 | 0.094245 |
| 8 | 1.064590 | -0.435657 | 1.348334 | 0.214171 | 0.447584 |
| 9 | -0.470874 | 0.197894 | 1.014463 | -1.507990 | -1.676122 |
In [4]:
from __future__ import division
from IPython.display import display
from sympy.interactive import printing
printing.init_printing(use_latex='mathjax')
import sympy as sym
from sympy import *
x, y, z = symbols("x y z")
k, m, n = symbols("k m n", integer=True)
f, g, h = map(Function, 'fgh')
Elementary operations¶
In [5]:
Rational(3,2)*pi + exp(I*x) / (x**2 + y)
Out[5]:
$\displaystyle \frac{3 \pi}{2} + \frac{e^{i x}}{x^{2} + y}$
In [6]:
exp(I*x).subs(x,pi).evalf()
Out[6]:
$\displaystyle -1.0$
In [7]:
e = x + 2*y
In [8]:
srepr(e)
Out[8]:
"Add(Symbol('x'), Mul(Integer(2), Symbol('y')))"
In [9]:
exp(pi * sqrt(163)).evalf(50)
Out[9]:
$\displaystyle 262537412640768743.99999999999925007259719818568888$
Algebra¶
In [10]:
eq = ((x+y)**2 * (x+1))
eq
Out[10]:
$\displaystyle \left(x + 1\right) \left(x + y\right)^{2}$
In [11]:
expand(eq)
Out[11]:
$\displaystyle x^{3} + 2 x^{2} y + x^{2} + x y^{2} + 2 x y + y^{2}$
In [12]:
a = 1/x + (x*sin(x) - 1)/x
a
Out[12]:
$\displaystyle \frac{x \sin{\left(x \right)} - 1}{x} + \frac{1}{x}$
In [13]:
simplify(a)
Out[13]:
$\displaystyle \sin{\left(x \right)}$
In [14]:
eq = Eq(x**3 + 2*x**2 + 4*x + 8, 0)
eq
Out[14]:
$\displaystyle x^{3} + 2 x^{2} + 4 x + 8 = 0$
In [15]:
solve(eq, x)
Out[15]:
$\displaystyle \left[ -2, \ - 2 i, \ 2 i\right]$
In [16]:
a, b = symbols('a b')
Sum(6*n**2 + 2**n, (n, a, b))
Out[16]:
$\displaystyle \sum_{n=a}^{b} \left(2^{n} + 6 n^{2}\right)$
Calculus¶
In [17]:
limit((sin(x)-x)/x**3, x, 0)
Out[17]:
$\displaystyle - \frac{1}{6}$
In [18]:
(1/cos(x)).series(x, 0, 6)
Out[18]:
$\displaystyle 1 + \frac{x^{2}}{2} + \frac{5 x^{4}}{24} + O\left(x^{6}\right)$
In [19]:
diff(cos(x**2)**2 / (1+x), x)
Out[19]:
$\displaystyle - \frac{4 x \sin{\left(x^{2} \right)} \cos{\left(x^{2} \right)}}{x + 1} - \frac{\cos^{2}{\left(x^{2} \right)}}{\left(x + 1\right)^{2}}$
In [20]:
integrate(x**2 * cos(x), (x, 0, pi/2))
Out[20]:
$\displaystyle -2 + \frac{\pi^{2}}{4}$
In [21]:
eqn = Eq(Derivative(f(x),x,x) + 9*f(x), 1)
display(eqn)
dsolve(eqn, f(x))
$\displaystyle 9 f{\left(x \right)} + \frac{d^{2}}{d x^{2}} f{\left(x \right)} = 1$
Out[21]:
$\displaystyle f{\left(x \right)} = C_{1} \sin{\left(3 x \right)} + C_{2} \cos{\left(3 x \right)} + \frac{1}{9}$
Illustrating Taylor series¶
We will define a function to compute the Taylor series expansions of a symbolically defined expression at various orders and visualize all the approximations together with the original function
In [22]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
In [23]:
# You can change the default figure size to be a bit larger if you want,
# uncomment the next line for that:
#plt.rc('figure', figsize=(10, 6))
In [24]:
def plot_taylor_approximations(func, x0=None, orders=(2, 4), xrange=(0,1), yrange=None, npts=200):
"""Plot the Taylor series approximations to a function at various orders.
Parameters
----------
func : a sympy function
x0 : float
Origin of the Taylor series expansion. If not given, x0=xrange[0].
orders : list
List of integers with the orders of Taylor series to show. Default is (2, 4).
xrange : 2-tuple or array.
Either an (xmin, xmax) tuple indicating the x range for the plot (default is (0, 1)),
or the actual array of values to use.
yrange : 2-tuple
(ymin, ymax) tuple indicating the y range for the plot. If not given,
the full range of values will be automatically used.
npts : int
Number of points to sample the x range with. Default is 200.
"""
if not callable(func):
raise ValueError('func must be callable')
if isinstance(xrange, (list, tuple)):
x = np.linspace(float(xrange[0]), float(xrange[1]), npts)
else:
x = xrange
if x0 is None: x0 = x[0]
xs = sym.Symbol('x')
# Make a numpy-callable form of the original function for plotting
fx = func(xs)
f = sym.lambdify(xs, fx, modules=['numpy'])
# We could use latex(fx) instead of str(), but matploblib gets confused
# with some of the (valid) latex constructs sympy emits. So we play it safe.
plt.plot(x, f(x), label=str(fx), lw=2)
# Build the Taylor approximations, plotting as we go
apps = {}
for order in orders:
app = fx.series(xs, x0, n=order).removeO()
apps[order] = app
# Must be careful here: if the approximation is a constant, we can't
# blindly use lambdify as it won't do the right thing. In that case,
# evaluate the number as a float and fill the y array with that value.
if isinstance(app, sym.core.numbers.Number):
y = np.zeros_like(x)
y.fill(app.evalf())
else:
fa = sym.lambdify(xs, app, modules=['numpy'])
y = fa(x)
tex = sym.latex(app).replace('$', '')
plt.plot(x, y, label=r'$n=%s:\, %s$' % (order, tex) )
# Plot refinements
if yrange is not None:
plt.ylim(*yrange)
plt.grid()
plt.legend(loc='best').get_frame().set_alpha(0.8)
With this function defined, we can now use it for any sympy function or expression
In [25]:
plot_taylor_approximations(sin, 0, [2, 4, 6], (0, 2*pi), (-2,2))
In [26]:
plot_taylor_approximations(cos, 0, [2, 4, 6], (0, 2*pi), (-2,2))
This shows easily how a Taylor series is useless beyond its convergence radius, illustrated by a simple function that has singularities on the real axis:
In [27]:
# For an expression made from elementary functions, we must first make it into
# a callable function, the simplest way is to use the Python lambda construct.
plot_taylor_approximations(lambda x: 1/cos(x), 0, [2,4,6], (0, 2*pi), (-5,5))
Basic Numerical Integration: the Trapezoid Rule¶
A simple illustration of the trapezoid rule for definite integration:
$$
\int_{a}^{b} f(x)\, dx \approx \frac{1}{2} \sum_{k=1}^{N} \left( x_{k} - x_{k-1} \right) \left( f(x_{k}) + f(x_{k-1}) \right).
$$
First, we define a simple function and sample it between 0 and 10 at 200 points
In [28]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
In [29]:
def f(x):
return (x-3)*(x-5)*(x-7)+85
x = np.linspace(0, 10, 200)
y = f(x)
Choose a region to integrate over and take only a few points in that region
Plot both the function and the area below it in the trapezoid approximation
In [31]:
plt.plot(x, y, lw=2)
plt.axis([0, 9, 0, 140])
plt.fill_between(xint, 0, yint, facecolor='gray', alpha=0.4)
plt.text(0.5 * (a + b), 30,r"$\int_a^b f(x)dx$", horizontalalignment='center', fontsize=20);
Compute the integral both at high accuracy and with the trapezoid approximation
In [32]:
from __future__ import print_function
from scipy.integrate import quad
integral, error = quad(f, a, b)
integral_trapezoid = sum( (xint[1:] - xint[:-1]) * (yint[1:] + yint[:-1]) ) / 2
print("The integral is:", integral, "+/-", error)
print("The trapezoid approximation with", len(xint), "points is:", integral_trapezoid)
The integral is: 565.2499999999999 +/- 6.275535646693696e-12 The trapezoid approximation with 5 points is: 559.890625