**LL**:*Programming for Computations - Python*by Svein Linge and Hans Petter Langtangen, 2nd edition.**L**:*A Primer on Scientific Programming with Python*by Hans Petter Langtangen, 2nd edition.**TAK**:*Applied Scientific Computing With Python*by Peter R. Turner, Thomas Arildsen, and Kathleen Kavanagh.

The final exam will focus on material not tested on the midterm. However, this class is cumulative, and it is necessary to understand earlier material as well. In particular, it is expected that students have developed a working knowledge of Python.

The following is a list of topics that may be explicitly tested on the final exam.

Numerical Differentiation following TAK § 3.2. Notebook: Numerical Differentiation [download] [view] [azure] |

Numerical Integration following TAK § 3.3. |

Numerical Integration following TAK § 3.4. |

Matrices and vectors in Numpy following LL § 2.3. |

Gaussian Elimination following TAK § 4.2. |

LU Factorization and Differential Equations (TAK § 4.3). |

Least Squares (TAK § 4.5). |

Eigenvalues (TAK § 4.6). Notebook: Eigenvalues and Eigenvectors [download] [view] [azure] |

Dictionaries (L § 6.2). |

Begin discussing classes following L Chapter 7. |

Random Numbers following L Chapter 8. |

Monte Carlo Integration following L Chapter 8. |

- The homework assignments and solutions.
- The textbooks listed above and exercises on the topics above.
- The notebooks available on the Class' Calendar page
- Class videos, available from Blackboard Collaborate.

Below some problems are included that are aimed at reviewing the topics above. They do not necessarily represent problems that would be assigned on the final exam.

In [ ]:

```
# Standard imports
import numpy as np
import matplotlib.pyplot as plt
import math as m
from mpmath import mp, iv
from scipy import linalg
from numpy import random
import random as random_number
```

Write a funtion `derivative(f, h)`

which takes as input:

- A function
`f(x)`

, which takes as input a floating point real number and produces a floating point real number as its output. - A float $h>0$.

The `derivative`

function should return a new function `df(x)`

. The function `df(x)`

should take as input a floating point real number `x`

and return the symmetric two-point approximation of the derivative of $f$ at $x$ using points at distance $h$ from $x$.

The `derivative`

function should make use of currying.

In [ ]:

```
# Test
f = lambda x: np.sin(x)
df = derivative(f, 10**-8)
for x in np.linspace(0, 3, 201):
assert abs( df(x) - np.cos(x) ) < 10**-8, \
"The value df({:0.2f}) is incorrect.".format(x)
```

A malfunctioning drone crashes destroying itself. The dictionary `altitude`

below maps time (in seconds) to measurements of a drone's altitude (in meters). The drone crashes shortly after the $5$ second mark.

- Approximately fast was the drone falling when it crashed into the ground? Store the answer (which would be in meters per second) in the variable
`crash_speed`

. - Before the drone crashed, it accellerated upward for some unknown reason. What was the fastest speed it was moving upward over the provided time interval? At approximately what time was it moving at that speed? Store the respective answers in
`up_speed`

and`up_time`

.

In [ ]:

```
altitude = {0.00: 2.0438, 0.10: 2.1013, 0.20: 2.1577, 0.30: 2.2136, 0.40: 2.2699,
0.50: 2.3271, 0.60: 2.3859, 0.70: 2.4472, 0.80: 2.5117, 0.90: 2.5799,
1.00: 2.6528, 1.10: 2.7311, 1.20: 2.8154, 1.30: 2.9066, 1.40: 3.0053,
1.50: 3.1123, 1.60: 3.2281, 1.70: 3.3535, 1.80: 3.4888, 1.90: 3.6347,
2.00: 3.7915, 2.10: 3.9593, 2.20: 4.1384, 2.30: 4.3285, 2.40: 4.5294,
2.50: 4.7406, 2.60: 4.9613, 2.70: 5.1904, 2.80: 5.4263, 2.90: 5.6672,
3.00: 5.9109, 3.10: 6.1543, 3.20: 6.3941, 3.30: 6.6263, 3.40: 6.8462,
3.50: 7.0480, 3.60: 7.2256, 3.70: 7.3716, 3.80: 7.4775, 3.90: 7.5341,
4.00: 7.5307, 4.10: 7.4553, 4.20: 7.2945, 4.30: 7.0336, 4.40: 6.6560,
4.50: 6.1436, 4.60: 5.4763, 4.70: 4.6323, 4.80: 3.5874, 4.90: 2.3156,
5.00: 0.7886, }
```

This problem concerns the following quadrature rule, defined for any interval $[a, b]$: $$Q_{[a,b]}(f) = \frac{1}{4} f\left(\frac{2a+b}{3}\right) + \frac{1}{2} f\left(\frac{a+b}{2}\right) + \frac{1}{4} f\left(\frac{a+2b}{3}\right).$$

Define a function

`q(a , b, f)`

which takes as input floating point real numbers $a$ and $b$ with $a<b$ and a function`f`

representing a real-valued function $f:{\mathbb R} \to {\mathbb R}$. The function should return $Q_{[a,b]}(f)$.What is the degree of of precision of $Q_{[a,b]}$? Store the answer in the variable

`q_degree`

.

Write a function `q_compound(n, f)`

which takes as input a positive integer $n$ and a function $f$ as above. The function should divide the interval $[0,1]$ into $n$ equal sized subintervals. It should return the sum of $Q_{\ast}(f)$ where $\ast$ varies over the $n$ subintervals.

Write a function `manipulate1(m)`

which takes as input a $2$-dimensional numpy array `m`

representing a matrix $M$. The function should make a copy of `m`

(so that `m`

is not changed). Then the following operations should be performed, in order:

- All diagonal entries should be doubled.
- The top right and bottom left entries should be swapped.
- The entries which are strictly closer to the bottom right courner than any other corner should be negated.

For example, the if the following matrix is input for `manipulate(m)`

$$\left(\begin{matrix}
2 & 3 & 4 \\
1 & -1 & 3 \\
8 & 5 & 11
\end{matrix}\right)$$
then the following matrix should be output:
$$\left(\begin{matrix}
4 & 3 & 8 \\
1 & -2 & 3 \\
4 & 5 & -22
\end{matrix}\right).$$

In [ ]:

```
# Test
m = np.array([[2, 3, 4],
[1, -1, 3],
[8, 5, 11]])
ret = manipulate1(m)
assert m!=ret, "m should not change!"
assert ret == np.array([[4, 3, 8],
[1, -2, 3],
[4, 5, -22]])
```

Write a function `manipulate2(m)`

which takes as input a $2$-dimensional numpy array `m`

representing a matrix $M$. The function should return the matrix
$$\left(\begin{array}{r|r}
4I & 3 M \\
\hline
M^T & M^T M
\end{array}\right).$$
Here, $M^T$ represents the transpose of $M$. The matrix $I$ represents the identity matrix whose sizes are chosen to make the above block matrix form well-defined. The matrix should be returned as a $2$-dimensional numpy array.

**Minor changes to this problem's wording**

Write a function `elimination(m)`

which takes as input a $2$-dimensional numpy array `m`

of floats representing a matrix $M$ with at least as many columns as rows. The function should make a copy of `m`

so that this matrix is unchanged. The function should then perform operations of adding multiples of one column to another until the matrix has an upper triangular form. Other operations are not permissible. The upper triangular matrix should be returned.

Consider the differential equation $$y''(x) + \sin(x) y(x) = x^2$$ on the interval $[0, \pi]$ subject to the boundary conditions $y(0)=1$ and $y(\pi)=0$.

Given an integer $n \geq 3$, consider cutting the interval into $n$ equal-sized subintervals. You are interested in approximating the solution $y$ at the cut points $$x_0 = \frac{1}{n}, x_1 = \frac{2}{n}, \ldots, x_{n-1}=\frac{n-1}{n}.$$ Set $y_i = y(x_i)$. Using numerical derivatives, find an $n-1 \times n-1$ matrix $M$ and a vector $v \in {\mathbb R}^{n-1}$ such that solving $M y= v$ leads to an approximate solution $y=(y_0, \ldots, y_{n-1})$ for the differential equation.

Write a function `linearization(n)`

which returns a pair of numpy arrays `(m,v)`

. Here `m`

should represent the matrix $M$ from the previous paragraph and `v`

should reprensent the vector $v$.

Consider the $L^2$ inner product on $[0, 1]$ given by
$$\langle f, g \rangle = \int_0^1 f(x)g(x)~dx.$$
This inner product has been implemented using the `quad`

function from `scipy.integrate`

below.

In [ ]:

```
from scipy.integrate import quad
def ip(f, g):
return quad(lambda x:f(x)*g(x), 0, 1)[0]
```

Given continuous functions $f$ and $g$ (which are not zero on $[0,1]$), there is always a unique constant $c$ such that $f$ and $g+cf$ are orthogonal with respect to $\langle, \rangle$. Here $g+cf$ denotes the function sending $x$ to $g(x)+c f(x)$.

Write a function `orth(f, g)`

which returns this constant `c`

.

Write a function `projection(f, g, h)`

which takes as input three functions representing functions ${\mathbb R} \to {\mathbb R}$. The function should return a pair of real numbers $(a,b)$ so that the function $a f + b g$ is the orthogonal projection of $h$ to the span of $f$ and $g$. Orthogonality should be measured with respect to the inner product from problem 9.

The matrix $A$ below has a real eigenvalue $\lambda$ which is very close to $-3$.

In [ ]:

```
A = np.array([[ 4, 1, 2, 4],
[-2, 4, -1, 2],
[ 1, -2, 3, 4],
[ 3, 4, 2, 1]
])
```

We define $B = (A+3I)^{-1}$. This matrix is approximately given below.

In [ ]:

```
B = np.array([[-0.36363636, -0.54545455, -0.27272727, 0.90909091],
[-0.81818182, -0.72727273, -0.36363636, 1.54545455],
[-1.40909091, -1.36363636, -0.43181818, 2.52272727],
[ 1.79545455, 1.81818182, 0.78409091, -3.23863636]])
```

Use powers of $B$ to compute the eigenvalue `lambda_A`

of $A$ which is closest to $-3$. Store the associated eigenvector in `eigenvector_A`

. Both quantities should be correct to six decimal places.

The matrix $C$ below has four distinct eigenvalues. Three of them are given below (correct to several digits of accuracy). Find the fourth and store the result in `eigvalue_3`

to approximately the same accuracy.

In [ ]:

```
C = np.array([[ 5, 1, 2, 4],
[-2, 4, -1, 2],
[ 1, -2, 3, 4],
[ 3, 4, 2, 10]
])
eigvalue_0 = 12.92391499
eigvalue_1 = 1.23951179
eigvalue_2 = 1.60202972
```

The deltoid is the region in the plane satisfying $$(x^2+y^2)^2 + 18(x^2+y^2)-27 < 8(x^3-3xy^2).$$ The deltoid is contained in the rectangle $$R=[-2, 3] \times [-3, 3] = \{(x,y):~-2 \leq x \leq 3 \text{ and } -3 \leq y \leq 3\}.$$

Write a function `deltoid_area(n)`

which takes as input a positive integer $n$. The function should choose $n$ uniformly random points from $R$ and check if they are contained in the deltoid. The number of points in the deltoid should be used to estimate its area using Monte Carlo Integration, and this area estimate should be returned.

**Classes will not be tested!**

Write a class `QField`

such that objects of the class represent numbers of the form $a+b\sqrt{2}$ where $a$ and $b$ are integers. For example, `Qfield(3, 5)`

should return a `QField`

object representing $3 + 5 \sqrt{2}$.

If `x`

and `y`

are two `QField`

objects representing $a_x+b_x \sqrt{2}$ and $a_y+b_y \sqrt{2}$ respectively, the following should be possible:

`x.float()`

should return a floating point approximation for $a_x+b_x \sqrt{2}$.`x.plus(y)`

should return the`QField`

object representing $a_x+a_y+(b_x+b_y) \sqrt{2}$.`x.times(y)`

should similarly return the`QField`

object representing the product (which also can be expressed as an integer plus an integer times the square root of two).

Recall that a polynomial
$$p(x) = \sum_{k=0}^d c_k x^k$$
can be represented by a dictionary with one entry for every non-zero coefficient. The dictionary representing $p$ would then map a degree $k$ with a non-zero $c_k$ to the floating point number $c_k$. For example `{9:3.0, 7:2.0}`

would represent $3 x^9 + 2x^7$.

Write a function `shift(p, t)`

which takes as input a dictionary `p`

representing a polynomial $p(x)$ and a floating point real number $t$. The function should return the dictionary representing the polynomial
$$q(x) = p(x+t).$$

*Help:* To do this, you need to know the coefficients of $(x+t)^k$. We have
$$(x+t)^k = \sum_{i=0}^k \binom{k}{i} t^{k-i} x^i.$$
The binomial coefficient $\binom{k}{i}$ can be computed using the `binom`

function from the `scipy.special`

package. For example, $\binom{4}{2}=6$ and this can be seen using the following python code:

```
from scipy.special import binom
binom(4, 2)
```