This is the solution to Homework 8: Problems - Python array computing and plotting.

The following figure illustrates the grade distribution for this homework.

This homework aims at giving you some experience with Python’s array computing and plotting features.

1.  The while-loop implementation of a for-loop. Consider the following mathematical function resembling a Hat function,

A scalar implementation of this function would be,

def hatFunc(x):
if x < 0:
return 0.0
elif 0 <= x < 1:
return x
elif 1 <= x < 2:
return 2 - x
elif x >= 2:
return 0.0


Write a vectorized version of this function. (Hint: you may need numpy’s logical_and method for building the vectorized version of this function.)

def hatFunc(x):
condition1 = x < 0
condition2 = np.logical_and(0 <= x, x < 1)
condition3 = np.logical_and(1 <= x, x < 2)
condition4 = x >= 2
r = np.zeros(len(x))
r[condition1] = 0.0
r[condition2] = x[condition2]
r[condition3] = 2-x[condition3]
r[condition4] = 0.0
return r


2.  The vertical position $y(t)$ of a ball thrown upward is given by $y(t)=v_0t-\frac{1}{2}gt^2$, where $g$ is the acceleration of gravity and $v_0$ is the initial vertical velocity at $t=0$. Two important physical quantities in this context are the potential energy, obtained by doing work against gravity, and the kinetic energy, arising from motion. The potential energy is defined as $P=mgy$, where $m$ is the mass of the ball. The kinetic energy is defined as $K=\frac{1}{2}mv^2$, where $v$ is the velocity of the ball, related to $y$ by $v(t)=y’(t)$.

Write a program that can plot $P(t)$ and $K(t)$ in the same plot, along with their sum $E = P + K$. Let $t\in[0,2v_0/g]$. Write your program such that $m$ and $v_0$ are read from the command line. Run the program with various choices of $m$ and $v_0$ and observe that $P+K$ always remains constant in this motion, regardless of initial conditions. This is in fact, the fundamental principle of conservation of energy in Physics.

A sample code can be downloaded from here. Here is an example output figure of the code:

3.  Integration by midpoint rule: The idea of the Midpoint rule for integration is to divide the area under a curve $f(x)$ into $n$ equal-sized rectangles. The height of the rectangle is determined by the value of $f$ at the midpoint of the rectangle. The figure below illustrates the idea,

To implement the midpointrule, one has to compute the area of each rectangle, sum them up, just as in the formula for the Midpoint rule,

where $h=(b-a)/n$ is the width of each rectangle. Implement this formula as a Python function midpoint(f, a, b, n) and test the integrator with the following example input mathematical functions.

An example code can be downloaded from here. Here is the output of the code,

In [38]: run midpoint.py
The exact integral of exp(x) between 0.00000 and 1.09861 is 2.00000.            The approximate answer is 1.99899 giving an error of 0.00101
The exact integral of cos(x) between 0.00000 and 3.14159 is 0.00000.            The approximate answer is 0.00000 giving an error of 0.00000
The exact integral of sin(x) between 0.00000 and 3.14159 is 2.00000.            The approximate answer is 2.00825 giving an error of 0.00825
The exact integral of sin(x) between 0.00000 and 1.57080 is 1.00000.            The approximate answer is 1.00103 giving an error of 0.00103


4.  Visualize approximations in the Midpoint integration rule Now consider the following function,

which we wish to integrate using the midpoint integrator that you wrote in the previous example. Now write a new code that visualizes the midpoint rule, similar to in the following figure. (Hint: you will need to use the Matplotlib function fill_between and use this function to create the filled areas between f(x) and the approximating rectangles)