Chapter 5

Fill lists with function values. A function with many applications in science is defined as $$ h(x)=\frac{1}{\sqrt{2 \pi}} e^{-\frac{1}{2} x^{2}} $$ Fig. 5.13 Plot of the dimensionless temperature \(T(z, t)\) in the ground for two different \(t\) values and \(b=2\). Fill lists xlist and hlist with \(x\) and \(h(x)\) values for uniformly spaced \(x\) coordinates in \([-4,4]\). You may adapt the example in Chapter 5.2.1. Name of program file: fill_lists.py.

Apply a function to a vector. Given a vector \(v=(2,3,-1)\) and a function \(f(x)=x^{3}+x e^{x}+1\), apply \(f\) to each element in \(v .\) Then calculate \(f(v)\) as \(v^{3}+v * e^{v}+1\) using vector computing rules. Show that the two results are equal.

Simulate by hand a vectorized expression. Suppose \(\mathrm{x}\) and \(\mathrm{t}\) are two arrays of the same length, entering a vectorized expression, $$ y=\cos (\sin (x))+\exp (1 / t) $$ If \(\mathrm{x}\) holds two elements, 0 and 2, and \(\mathrm{t}\) holds the elements 1 and 1.5, calculate by hand (using a calculator) the y array. Thereafter, write a program that mimics the series of computations you did by hand (typically a sequence of operations of the kind we listed on page 182 - use explicit loops, but at the end you can use Numerical Python functionality to check the results). Name of program file: simulate_vector_computing.py.

Plot a formula. Make a plot of the function \(y(t)=v_{0} t-0.5 g t^{2}\) for \(v_{0}=10, g=9.81\) and \(t \in\left[0,2 v_{0} / g\right] .\) The label on the \(x\) axis should be 'time (s)' and the label on the \(y\) axis should be 'height \((\mathrm{m}) '\). Name of program file: plot_ball1.py.

Plot a formula for several parameters. Make a program that reads a set of \(v_{0}\) values from the command line and plots the corresponding curves \(y(t)=v_{0} t-0.5 g t^{2}\) in the same figure (set \(g=9.81\) ). Let \(t \in\left[0,2 v_{0} / g\right]\) for each curve, which implies that you need a different vector of \(t\) coordinates for each curve. Name of program file: plot_bal12.py.

Plot exact and inexact Fahrenheit-Celsius formulas. Exercise \(2.2\) introduces a simple rule to quickly compute the Celsius temperature from the Fahreheit degrees: \(C=(F-30) / 2\). Compare this curve against the exact curve \(C=(F-32) 5 / 9\) in a plot. Let \(F\) vary between \(-20\) and \(120 .\) Name of program file: f2c_shortcut_plot.py. \(\diamond\)

Implement Lagrange's interpolation formula. Imagine we have \(n+1\) measurements of some quantity \(y\) that depends on \(x\) : \(\left(x_{0}, y_{0}\right),\left(x_{1}, y_{1}\right), \ldots,\left(x_{n}, y_{n}\right)\). We may think of \(y\) as a function of \(x\) and ask what \(y\) is at some arbitrary point \(x\) not coinciding with any of the \(x_{0}, \ldots, x_{n}\). This problem is known as interpolation. One way to solve this problem is to fit a continuous function that goes through all the \(n+1\) points and then evaluate this function for any desired \(x\). A candidate for such a function is the polynomial of degree \(n\) that goes through all the points. This polynomial can be written $$ p_{L}(x)=\sum_{k=0}^{n} y_{k} L_{k}(x) $$ where $$ L_{k}(x)=\prod_{i=0, i \neq k}^{n} \frac{x-x_{i}}{x_{k}-x_{j}} $$ The \(\Pi\) notation corresponds to \(\sum\), but the terms are multiplied, e.g., $$ \prod_{i=0, i \neq k}^{n} x_{i}=x_{0} x_{1} \cdots x_{k-1} x_{k+1} \cdots x_{n} $$ The polynomial \(p_{L}(x)\) is known as Lagrange's interpolation formula, and the points \(\left(x_{0}, y_{0}\right), \ldots,\left(x_{n}, y_{n}\right)\) are called interpolation points. Make a function Lagrange \(\left(x\right.\), points) that evaluates \(p_{L}\) at the point x, given \(n+1\) interpolation points as a two-dimensional array points, such that points \([i, 0]\) is the \(x\) coordinate of point number \(i\) and points \([i, 1]\) is the corresponding \(y\) coordinate. To verify the program, we observe that \(L_{k}\left(x_{k}\right)=1\) and that \(L_{k}\left(x_{i}\right)=0\) for \(i \neq k\), implying that \(p_{L}\left(x_{k}\right)=y_{k}\). Write a function verify (points) that computes \(\left|p_{L}\left(x_{k}\right)-y_{k}\right|\) at all the interpolation points and checks that the value is approximately zero. Call verify with 5 equally spaced points along the curve \(y=\sin (x)\) for \(x \in[0, \pi]\). Then evaluate \(p_{L}(x)\) for an \(x\) in the middle of two interpolation points and compare the value of \(p_{L}(x)\) with the exact one: \(\sin (x)\). Name of program file: Lagrange_polynomial1.py.

Plot a wave packet. The function. $$ f(x, t)=e^{-(x-3 t)^{2}} \sin (3 \pi(x-t)) $$ describes for a fixed value of \(t\) a wave localized in space. Make a program that visualizes this function as a function of \(x\) on the interval \([-4,4]\) when \(t=0\). Name of program file: plot_wavepacket.py.

Plot the viscosity of water. The viscosity of water, \(\mu\), varies with the temperature \(T\) (in Kelvin) according to $$ \mu(T)=A \cdot 10^{B /(T-C)} $$ where \(A=2.414 \cdot 10^{-5} \mathrm{~Pa} \mathrm{~s}, B=247.8 \mathrm{~K}\), and \(C=140 \mathrm{~K}\). Plot \(\mu(T)\) for \(T\) between 0 and 100 degrees Celsius. Label the \(x\) axis with 'temperature (C)' and the \(y\) axis with 'viscosity (Pa s)'. Note that \(T\) in the formula for \(\mu\) must be in Kelvin. Name of program file: water_viscosity.py.

Explore a function graphically. The wave speed \(c\) of water surface waves depends on the length \(\lambda\) of the waves. The following formula relates \(c\) to \(\lambda\) : $$ c(\lambda)=\sqrt{\frac{g \lambda}{2 \pi}\left(1+s \frac{4 \pi^{2}}{\rho g \lambda^{2}}\right) \tanh \left(\frac{2 \pi h}{\lambda}\right)} $$ where \(g\) is the acceleration of gravity, \(s\) is the air-water surface tension \(\left(7.9 \cdot 10^{-4} \mathrm{~N} / \mathrm{cm}\right), \rho\) is the density of water (can be taken as \(\left.1 \mathrm{~kg} / \mathrm{cm}^{3}\right)\), and \(h\) is the water depth. Let us fix \(h\) at \(50 \mathrm{~m}\). First make a plot of \(c(\lambda)\) for small \(\lambda(1 \mathrm{~mm}\) to \(10 \mathrm{~cm})\). Then make a plot \(c(\lambda)\) for larger \(\lambda\) ) (1 m to \(2 \mathrm{~km}\) ). Name of program file: water_wave_velocity.py.