2023年11月25日

Python代写 | AMATH 301 – Autumn 2020 Homework 8

本次Python代写是编程完成一些数学应用问题

AMATH 301 – Autumn 2020
Homework 8

1. Recall the pendulum equations we encountered in Homework 7. The nonlinear
pendulum equation with damping was
¨θ = −
g
L
sin(θ) − σ ˙θ. (1)
The linear pendulum equation with damping was
¨θ = −
g
L
θ − σ ˙θ. (2)
Here g = 9.8 is the acceleration due to gravity, L = 21 is the length of the pendulum,
and σ = 0.08 is a damping coefficient. The unknown that we must solve for is θ
(theta), the angle of deflection of the pendulum from the vertical. Throughout we
will use the initial conditions θ(0) = π/6 and v(0) = 0.1.
In this problem we are going to find the largest ∆t for which forward Euler is stable
on the linear pendulum equation (2). We will then use this timestep to calculate
the solution of (2) and (1) using forward Euler.
(a) Recall that we can convert (2) to a system of linear equations and write it in
the matrix form x˙ = Ax. Define A in your code. (You do not need to submit
anything here, but you will need to use it in what follows).
1
(b) Recall that the forward Euler formula may be rewritten as
xk+1 = (I + ∆tA) xk,
where I is the 2 × 2 identity matrix. Forward Euler is stable for this problem
if |λ| < 1 for all eigenvalues λ of I + ∆tA. Find the eigenvalue of I + ∆tA that
is largest in magnitude (or absolute value) when ∆t = 0.02. Save the absolute
value of this largest eigenvalue to the variable A1. You may want to recall how
the infinity norm works for calculating this.
(c) We now want to find the largest ∆t for which forward Euler is stable. Using
20 logarithmically spaced points (logspace) between 100 and 10−4
for the
possible values of ∆t, find the largest ∆t for which forward Euler is stable.
Save this value of ∆t to the variable A2. In what follows, we will call this
∆topt, the optimal ∆t.
(d) We are now going to solve (2) using ∆topt. Create an array/vector that holds
the values of t for 0 ≤ t ≤ 100 with a step size equal to ∆topt. Warning for
python users. In python you usually have to use T + dt when creating this
array. But doing so here will result in an array whose final value is greater
than 100. Remove the + dt term to make sure that this array does not go
over 100.
You should find that the last value in this array/vector is not actually equal
to 100, but it should be close. Save this array/vector as an N × 1 vector to
the variable A3, where N is the size of the array/vector holding the t values at
which we find the solution.
(e) Solve (2) for 0 ≤ t ≤ 100 using forward Euler with ∆topt. Save the resulting
2 × N matrix to the variable A4.
(f) We now want to solve the nonlinear problem using ∆topt. Remember that we
are not guaranteed that forward Euler will be stable for the nonlinear problem
since the stability results are about the linear problem, but it is our best hope.
Using ∆topt, solve (1) using forward Euler for 0 ≤ t ≤ 100. Save the resulting
2 × N matrix to the variable A5.
(g) We now want to solve the two problems for ∆t just larger than ∆topt. In
problem 1(c) we used 20 logarithmically spaced points between 100 and 10−4
.
If ∆topt was the nth value in that list, we want to use the n−1st value in that
list. Call this ∆t, ∆tbad. You should have ∆tbad > ∆topt. (Hint: you may find
the functions find (MATLAB) and np.argwhere() (python) helpful here).
Solve (2) using forward Euler and ∆tbad for 0 ≤ t ≤ 100. Save the resulting
2 × M matrix to the variable A6, where M is the length of the array/vector
that holds the t values for 0 ≤ t ≤ 100 with spacing ∆tbad.
(h) Solve (1) using forward Euler for and ∆tbad 0 ≤ t ≤ 100. Save the resulting
2 × M matrix to the variable A7.
2. We are now going to solve the linear and the nonlinear equations using RK4 (not
ode45 or scipy.integrate.solve_ivp). First, note that RK4, applied to a linear
2
ODE of the form x˙ = Ax can be written as
xk+1 =

I + ∆tA +
(∆t)
2
2
A
2 +
(∆t)
3
6
A
3 +
(∆t)
4
24
A
4

xk.
Therefore, RK4 is stable if all of the eigenvalues of the matrix
B = I + ∆tA +
(∆t)
2
2
A
2 +
(∆t)
3
6
A
3 +
(∆t)
4
24
A
4
(3)
are less than 1 in magnitude.
(a) Calculate the magnitude of the largest eigenvalue of B using ∆t = 0.2. Save
this value to the variable A8.
(b) Solve the linear system (2) for 0 ≤ t ≤ 100 using RK4 with ∆t = 0.2. Use the
same initial conditons defined in Problem 1. Save the result at t = 100 as a
2 × 1 column vector to the variable A9.
(c) Solve the nonlinear system (1) for 0 ≤ t ≤ 100 using RK4 with ∆t = 0.2. Use
the same initial conditions defined in Problem 1. Save the result at t = 100 as
a 2 × 1 column vector to the variable A10