Instructions: Create a new notebook called IntegrationExercises in your Integration directory and solve the following problems inside it. Be sure to include the problem statements in a markdown cell above your solution. You don't need to put the "helper" code in the markdown cell, just implement the helper code in your code cell with your solution.
Preliminaries: At the top of your notebook, include a "Heading 1" cell with the title Integration Exercises. Then include the inline functions and libraries by adding a code cell that invokes the %pylab inline magic and imports the needed packages.
(a) Write two functions trapz(func,a,b,N) and simps(func,a,b,N) to compute the integral of the function func over the variable x using the trapezoidal rule and Simpson's rule to a file called Integrators.py. Do not use the scipy.integrate built-in functions. Include docstrings with each function that describe what they do. Then import the module and use the functions to answer the following questions. To avoid namespace conflicts, import your module as myint and then call the functions from that namespace. i.e.
import Integrators as myint
#...your code...here
I = myint.trapz(func,a,b,N)
Also, be sure that you import any needed modules inside your Integrators.py file.
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
(b) Use your myint.trapz and myint.simps functions from part (a) to calculate the integral of $x^4 - 2x + 1$ from $x$ = 0 to $x$ = 2 with $N$ = 10, $N$ = 100, and $N$ = 1000. Then compare your result to the known correct value of 4.4. What is the percent error in each case?
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
Consider the integral
$$ E(x) = \int_0^x e^{-t^2} dt $$This is the error function, commonly seen in probability and statistics. There is no known way to perform this particular integral analytically, although the integrand can be expanded in a Taylor series and terms computed to arbitrary order. For most applications, numerical approaches are the only way forward.
(a) Use scipy.integrate.cumtrapz to calculate $E(x)$ for values of $x$ from 0 to 3 in steps of 0.1. Print the result.
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
(b) Plot the integrand as a function of $t$ (which is just $x$) and $E(x)$ as a function of $x$ ranging from 0 to 3 on the same graph.
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
Consider a particle in a one-dimensional box of width $L$. The probability of finding the particle between $a$ and $b$ is given by
$$P(a,b) = \displaystyle \int _a^b|\psi(x)|^2dx $$where
$$\psi(x) = \sqrt{ \frac{2}{L} } \sin\left(\frac{n \pi x}{L}\right) $$is the wavefunction.
(a) What is the probability of finding the particle between $L/3$ and $L/2$ for the ground state ($n$ = 1) and for the first excited state ($n$ = 2)? Let $L$ = 1. Perform the integral using both scipy.integrate.trapz and scipy.integrate.quad with an accuracy of 6 sig figs. How many slices did you need to use for scipy.integrate.trapz?
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
(b) Make a plot of $|\psi(x)|^2$ vs. $x$ for the first two excited states with $L$=1.
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
Integrals of mass density lead to three interesting quantities:
$$ M = \int \rho dV \\ \vec{r}_{cm} = (x_{cm},y_{cm},z_{cm}) = \frac{1}{M}\int \vec{r} \rho dV \\ I_{cm} = \int (r-r_{cm})^2 \rho dV $$where $M$ is the total mass (a scalar), $\vec{r}_{cm}$ is the center of mass position (a vector of 3 components), $I_{cm}$ are the moments of inertia about the center of mass (diagonals of a 3x3 matrix), and the mass density, $\rho$, may be a function of the spatial variables.
Consider a rectangular box: length (in $x$) = 0.2 m, width (in $y$) = 0.2 m, and height (in $z$) = 1.0 m centered on the origin, (0,0,0) and with a mass density, $\rho(x, y, z)$ = (100 kg/m$^4$) ($y$ + 0.1) + (100 kg/m$^5$) $z^2$ for $x$, $y$, and $z$ in meters.
(a) Use numerical integration to find $M$, $\vec{r}_{cm}$, and $I_{cm}$. Note that you need $M$ to compute $\vec{r}_{cm}$ and $\vec{r}_{cm}$ to compute $I_{cm}$. The moments of inertia $I_{cm}$ are $I_{xx}$, $I_{yy}$, and $I_{zz}$, where e.g. $I_{xx} = \int ((y-y_{cm})^2 + (z-z_{cm})^2 )\rho dV$, etc.
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
(b) The mass, center of mass, and moments of inertia can be easily computed analytically for a rectangular box of uniform density ($\rho$ = constant). (Compute them yourself or look them up). Verify your algorithm from part (a) works by having it compute $M$, $\vec{r}_{cm}$ and $I_{cm}$ for this test case. How accurate (how many sig figs?) are the numerical results?
In [ ]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell
(c) Finally, change the density function to something of your choosing (different from that used in parts (a) and (b)) and recompute $M$, $\vec{r}_{cm}$ and $I_{cm}$ for that case.
In [1]:
#Copy the exercise statement to a markdown cell in your notebook and then implement a solution in a code cell