- Details
- Parent Category: Mathematics Assignments' Solutions
We Helped With This Mathematics Assignment: Have A Similar One?
| Category | Math |
|---|---|
| Subject | Other |
| Difficulty | Graduate |
| Status | Solved |
| More Info | Pay Someone To Do Your Math Homework |
Short Assignment Requirements
Assignment Description
FIRST Currently, the cheb function defined in cheb_tools.py only produces the correct answer for the interval [-1,1]. If you run pytest you'll see that the function test_cheb_big_2 fails. That's because it's trying to test the cheb function on the interval [-3.2, 2.38].
To fix the code, you need to do two things.
1. The x values need to be transformed (by a linear change of variables) that maps [-1,1] to [x_left,x_right].
- The D matrix needs to be rescaled so that it correctly approximates a first derivative on the requested interval.
Hint: The second part is a calculus problem. How does a linear change of variables affect derivatives?
Cheb_tools.py
""" | |
Tools for Chebyshev differentiation. | |
""" | |
| |
import numpy as np | |
| |
def cheb(n, x_left=-1, x_right=1): | |
""" | |
Computes the Chebyshev grid and differentiation matrix on n+1 points. | |
Parameters | |
---------- | |
n : int | |
n+1 is the number of grid points. | |
x_left : float | |
Location of the left grid point. | |
x_right : float | |
Location of the right grid point. | |
Returns | |
------- | |
D : array_like | |
The (n+1)-by-(n+1) Chebyshev differentiation matrix. | |
x : array_like | |
The Chebyshev nodes on the interval [-x_left,x_right]. | |
Notes | |
----- | |
This function is based on cheb.m from Trefethen's Spectral | |
Methods in MATLAB. Unlike Trefethen's, however, this function | |
orders the nodes from left to right. That is, x[0] = x_left | |
and x[-1] = x_right | |
""" | |
| |
# Trefethen's algorithm | |
n_range = np.arange(n+1) | |
x = np.cos(np.pi*n_range/n) | |
c = np.array([2] + (n-1)*[1] + [2]) * (-1)**n_range | |
| |
X = np.tile(x, (n+1, 1)).T | |
dX = X - X.T | |
| |
D = np.outer(c, 1/c)/(dX + np.eye(n+1)) | |
D -= np.diag(np.sum(D, axis=1)) | |
| |
# swap orders of nodes so x=-1 is first | |
x = x[::-1] | |
D = D[::-1, ::-1] | |
| |
# change variables for given interval | |
# TODO: FIX THIS | |
| |
| |
return D, x |
SECOND
Next,
take a look at the function solve_poisson in elliptic_tools.py. Notice that there
is a missing piece to this function (# fill in the top and bottom rows). Just like with
the finite difference solver, you need to fill in the top and bottom rows
of D2 depending on
the type of boundary condition.
Hint: For Neumann conditions, keep in mind that the top and bottom rows of D contain spectral approximations of the first derivative at the left and right endpoints respectively.
elliptic_tools.py
""" | |
Functions for solving 1d elliptic problems. | |
""" | |
| |
import numpy as np | |
| |
def solve_poisson(D, x, f, left, right): | |
""" | |
Solves the 1D Poisson equation u'' = f with the specified boundary conditions. | |
Parameters | |
---------- | |
D : array_like | |
The Chebyshev differentiation matrix. | |
x : array_like | |
The Chebyshev nodes. | |
f : array_like | |
The array of right-hand side values (interior only). | |
left : tuple | |
A tuple of the form (type,val) where type is either 'D' for Dirichlet | |
data or 'N' for Neumann data and val is a float containing the | |
left boundary value data. | |
right : tuple | |
Similar to left, but for the right endpoint. | |
Returns | |
------- | |
u : array | |
An array of solution values on the grid. | |
Notes | |
----- | |
The array f should not include values at the endpoints. The returned array | |
will have two more elements than f. | |
""" | |
| |
# x has length n+1 | |
n = len(x) - 1 | |
| |
# compute the second derivative matrix | |
D2 = D.dot(D) | |
| |
# zero out the top and bottom rows | |
D2[0, :] = 0 | |
D2[-1, :] = 0 | |
| |
# initialize the right-hand side vector | |
b = np.zeros(n+1) | |
b[1:-1] = f | |
b[0] = left[1] | |
b[-1] = right[1] | |
| |
# fill in the the top and bottom rows | |
# TODO: FIX THIS | |
| |
| |
# solve the system | |
u = np.linalg.solve(D2, b) | |
return u |
