Homework 1¶
INSERT YOUR NAME HERE
Directions: Add work to this notebook to solve the problems below.
Be careful to name the objects you create as described in the problem. Some problems may include some tests you can run to test your code.
Some problems are derived from problems in the book:
- CM: Computational Mathematics with SageMath by Paul Zimmermann, Alexandre Casamayou, Nathann Cohen, Guillaume Connan, Thierry Dumont, Laurent Fousse, François Maltey, Matthias Meulien, Marc Mezzarobba, Clément Pernet, Nicolas M. Thiéry, Erik Bray, John Cremona, Marcelo Forets, Alexandru Ghitza, and Hugh Thomas. (freely available for download)
1. Numerical approximation¶
Compute the quantity
$$\sin(\sqrt[3]{2})$$
correct to $20$ digits of precision. Store the result in the variable value
.
Tests for your code:
If nothing is printed by the following commands then two tests are passed. This doesn't mean that your answer is necessarily correct!
try:
value_as_string = str(value) # Convert value to a string
assert value_as_string[1] == '.', 'The second character should be a decimal'
assert len(value_as_string) == 22, 'There should be 20 digits after the decimal.'
except Exception as e:
print('When trying to work with `value` an error occurred: {e}')
2. CM-Exercise 2¶
(This is a modification of CM Exercise 2)
Compute, without using the sum
command, the sum of $p$-powers of integers from $0$ to $n$, for $p = 0, \dots, 4$:
The following recurrence can be useful to compute this sum:
$$ S_p(n) = \frac{1}{p+1} \left( (n+1)^{p+1} - \sum_{j=0}^{p-1} \binom{p+1}{j} S_j(n) \right). $$This recurrence is easily obtained when computing by two different methods the telescopic sum $$\sum_{0 \leq k \leq n} \left( (k+1)^{p+1} - k^{p+1} \right).$$
To be concrete, define Sage
symbolic functions S_0(n)
, …, S_4(n)
. The first function S_0(n)
has been done for you. (Note that $0^0=1$.)
S_0(n) = n + 1
Tests for your code:
try:
failed = 0
if S_1(2) != 3:
print('S_1(2) is incorrect')
failed += 1
if S_2(2) != 5:
print('S_2(2) is incorrect')
failed += 1
if S_3(2) != 9:
print('S_3(2) is incorrect')
failed += 1
if S_4(2) != 17:
print('S_4(2) is incorrect')
failed += 1
if failed == 0:
print('Passed tests of the value S_p(2)')
except Exception as e:
print(f'When trying to check the values S_p(2), an error occurred: {e}')
try:
# These commands nicely print your answers:
show(LatexExpr('S_0(n)='), S_0(var('n')).factor())
show(LatexExpr('S_1(n)='), S_1(var('n')).factor())
show(LatexExpr('S_2(n)='), S_2(var('n')).factor())
show(LatexExpr('S_3(n)='), S_3(var('n')).factor())
show(LatexExpr('S_4(n)='), S_4(var('n')).factor())
except Exception as e:
print(f'When trying to print the values S_p(2), an error occurred: {e}')
3. Trigonometry to polynomial¶
There is a polynomial $p_{23}(x)$ such that $$p_{23}(\cos \theta) = \cos (23 \theta).$$
Find $p_{23}(x)$ and store it in the Python variable p23
. (There is a polynomial $p_n$ for every integer $n$, but I'm just asking for $p_{23}$.)
Hint: I was able to do this by manipulating expressions using trig_expand
, subs
, and expand
.
Tests for your code:
When $\cos \theta=1$, we have $\cos 23 \theta = 1$.
try:
assert bool(p23(x=1) == 1), 'Error: p_23(1) is not 1.'
except Exception as e:
print(f'When trying to evaluate `p23(x=1)`, an exception occurred: {e}')
4. Line from two circles¶
Let $x_1$, $y_1$, and $r_1$ be real numbers with $r_1 > 0$. This information determines a circle $C_1$, namely the solution set to the equation $$(x-x_1)^2 + (y-y_1)^2 = r_1^2.$$ The data $(x_2, y_2, r_2)$ determines another circle $C_2$.
Assuming that $C_1 \cap C_2$ consists of two points $P$ and $Q$, we may form a line $L=\overline{P Q}$. Further assuming that $L$ is not vertical, the line $L$ has a real slope $m$ and a $y$-intercept $b$. Write a function named circles_to_line
which given (x1, y1, r1, x2, y2, r2)
returns the pair (m, b)
.
You can use the following to format your answer:
def circles_to_line(x1, y1, r1, x2, y2, r2):
# Do some calculations to produce m and b
m = ???
b = ???
return (m, b)
(Remark: It makes sense to use Sage to compute the answer, but in the function you probably want simple formulas. Actually, if you use some geometry it is fairly easy to get formulas for $m$ and $b$ by hand.)
Tests for your code:
The circles with center $(0,0)$ and $(0,1)$ of radius $1$ intersect in the points $\left(\pm \frac{\sqrt{3}}{2}, \frac{1}{2}\right)$, so the following test should pass:
try:
m,b = circles_to_line(0, 0, 1, 0, 1, 1)
if m != 0:
print (f'The returned slope should be zero but is {m}')
elif b != 1/2:
print(f'The returned y-intercept should be 1/2 but is {b}')
else:
print('The return seems to be correct.')
except Exception as e:
print(f'When trying to work with `circles_to_line(0, 0, 1, 0, 1, 1)`, an exception occurred: {e}')
Now we'll plot a picture. If this looks correct, then your code probably works.
try:
m, b = circles_to_line(1, 2, 3, 4, 6, 5)
y(x) = m*x + b
show(circle((1,2), 3) + circle((4,6), 5) + plot(y, (x, -3, 8), color='red'))
except Exception as e:
print(f'When trying to work with `circles_to_line(1, 2, 3, 4, 5, 7)`, an exception occurred: {e}')
5. Tangent planes¶
Let $f(x,y)$ be a differentiable function and $(x_0, y_0)$ be a point in the domain of $f$. Write a function that returns z
where $z(x,y)$ is a function whose graph is the tangent plane to the graph of $f$ at $(x_0, y_0)$. Your solution should have the following form:
def tangent_plane(f, x0, y0):
x,y = var('x y')
# Here we do some calculations
z(x, y) = # some function of x and y and constants defined above.
return z
(Hint: In my solution, I used f.derivative()
.)
Tests for your code:
f(x,y) = x^2 + y^2
try:
z = tangent_plane(f, 0, 0) # This should be the zero function.
if bool(z==0):
print('Your code seems to handle f(x,y)=x^2+y^2 case correctly.')
else:
print(f'Your code seems to handle f(x,y)=x^2+y^2 incorrectly. The tangent plane is $z(x,y)=0$ but your code returned {z}.')
except Exception as e:
print(f'When calling and working with `z = tangent_plane(f, 0, 0)` an exception occurred: {e}')
Below we do a plot:
f(x, y) = cos(x)*sin(y)
x0 = pi/4
y0 = -pi/6
try:
z = tangent_plane(f, x0, y0) # This should be the zero function.
plt = plot3d(f, (x, x0-1/2, x0+2), (y, y0-1/2, y0+2))
plt += plot3d(z, (x, x0-1/2, x0+2), (y, y0-1/2, y0+2), color='yellow')
plt += point3d([(x0, y0, f(x0, y0))], color='red', size=10)
show(plt)
except Exception as e:
pass