Skip to content

Commit f3cb32a

Browse files
committed
update square root lecture
1 parent 5edbb02 commit f3cb32a

File tree

1 file changed

+38
-35
lines changed

1 file changed

+38
-35
lines changed

lectures/greek_square.md

Lines changed: 38 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -30,6 +30,8 @@ Invariant subspace methods are used throughout applied economic dynamics, for ex
3030

3131
Our approach here is to illustrate the method with an ancient example, one that ancient Greek mathematicians used to compute square roots of positive integers.
3232

33+
In this lecture we assume that we have yet
34+
3335
## Perfect Squares and Irrational Numbers
3436

3537
An integer is called a **perfect square** if its square root is also an integer.
@@ -528,20 +530,22 @@ In equation {eq}`eq:invariantsub101`, the $i$th $\lambda_i$ equals the $V_{i, 1
528530
The following graph verifies this for our example.
529531
530532
```{code-cell} ipython3
533+
:tags: [hide-input]
534+
531535
# Plotting the eigenvectors
532536
plt.figure(figsize=(8, 8))
533537
534-
plt.quiver(0, 0, V[0, 0], V[0, 1], angles='xy', scale_units='xy',
538+
plt.quiver(0, 0, V[0, 0], V[1, 0], angles='xy', scale_units='xy',
535539
scale=1, color='C0', label=fr'$\lambda_1={np.round(Λ[0], 4)}$')
536-
plt.quiver(0, 0, V[1, 0], V[1, 1], angles='xy', scale_units='xy',
540+
plt.quiver(0, 0, V[0, 1], V[1, 1], angles='xy', scale_units='xy',
537541
scale=1, color='C1', label=fr'$\lambda_2={np.round(Λ[1], 4)}$')
538542
539543
# Annotating the slopes
540-
plt.text(V[0, 0]-0.5, V[0, 1]*1.2,
541-
r'slope=$\frac{V_{1,1}}{V_{1,2}}=$'+f'{np.round(V[0, 0] / V[0, 1], 4)}',
544+
plt.text(V[0, 0]-0.5, V[1, 0]*1.2,
545+
r'slope=$\frac{V_{1,1}}{V_{1,2}}=$'+f'{np.round(V[0, 0] / V[1, 0], 4)}',
542546
fontsize=12, color='C0')
543-
plt.text(V[1, 0]-0.5, V[1, 1]*1.2,
544-
r'slope=$\frac{V_{2,1}}{V_{2,2}}=$'+f'{np.round(V[1, 0] / V[1, 1], 4)}',
547+
plt.text(V[0, 1]-0.5, V[1, 1]*1.2,
548+
r'slope=$\frac{V_{2,1}}{V_{2,2}}=$'+f'{np.round(V[0, 1] / V[1, 1], 4)}',
545549
fontsize=12, color='C1')
546550
547551
# Adding labels
@@ -590,8 +594,8 @@ Notice that it follows from
590594
591595
$$
592596
\begin{bmatrix} V^{1,1} & V^{1,2} \cr
593-
V^{2,2} & V^{2,2} \end{bmatrix} \begin{bmatrix} V_{1,1} & V_{1,2} \cr
594-
V_{2,2} & V_{2,2} \end{bmatrix} = \begin{bmatrix} 1 & 0 \cr 0 & 1 \end{bmatrix}
597+
V^{2,1} & V^{2,2} \end{bmatrix} \begin{bmatrix} V_{1,1} & V_{1,2} \cr
598+
V_{2,1} & V_{2,2} \end{bmatrix} = \begin{bmatrix} 1 & 0 \cr 0 & 1 \end{bmatrix}
595599
$$
596600
597601
that
@@ -650,8 +654,6 @@ Let's verify {eq}`eq:deactivate1` and {eq}`eq:deactivate2` below
650654
To deactivate $\lambda_1$ we use {eq}`eq:deactivate1`
651655
652656
```{code-cell} ipython3
653-
Λ, V = np.linalg.eig(M)
654-
655657
xd_1 = np.array((x_0[0],
656658
V[1,1]/V[0,1] * x_0[0]),
657659
dtype=np.float64)
@@ -678,47 +680,48 @@ np.round(V_inv @ xd_2, 8)
678680
We find $x_{2,0}^* = 0$.
679681
680682
```{code-cell} ipython3
681-
# Simulate the difference equation with muted λ1
683+
# Simulate with muted λ1 λ2.
682684
num_steps = 10
683-
xs_λ1, Λ, _, _ = iterate_M(xd_1, M, num_steps)
685+
xs_λ1 = iterate_M(xd_1, M, num_steps)[0]
686+
xs_λ2 = iterate_M(xd_2, M, num_steps)[0]
684687
685-
# Simulate the difference equation with muted λ2
686-
xs_λ2, _, _, _ = iterate_M(xd_2, M, num_steps)
688+
# Compute ratios y_t / y_{t-1}
689+
ratios_λ1 = xs_λ1[1, 1:] / xs_λ1[1, :-1]
690+
ratios_λ2 = xs_λ2[1, 1:] / xs_λ2[1, :-1]
691+
```
687692
688-
# Compute ratios y_t / y_{t-1} with higher precision
689-
ratios_λ1 = xs_λ1[1, 1:] / xs_λ1[1, :-1] # Adjusted to second component
690-
ratios_λ2 = xs_λ2[1, 1:] / xs_λ2[1, :-1] # Adjusted to second component
693+
```{code-cell} ipython3
694+
:tags: [hide-input]
691695
692-
# Plot the ratios for y_t with higher precision
696+
# Plot the ratios for y_t
693697
plt.figure(figsize=(14, 6))
694698
695699
plt.subplot(1, 2, 1)
696-
plt.plot(np.round(ratios_λ1, 6), label='Ratios $y_t / y_{t-1}$ after muting $\lambda_1$', color='blue')
697-
plt.axhline(y=Λ[1], color='red', linestyle='--', label='$\lambda_2$')
698-
plt.xlabel('Time')
699-
plt.ylabel('Ratio $y_t / y_{t-1}$')
700-
plt.title('Ratios after Muting $\lambda_1$')
700+
plt.plot(np.round(ratios_λ1, 6),
701+
label=r'$\frac{y_t}{y_{t-1}}$', linewidth=3)
702+
plt.axhline(y=Λ[1], color='red',
703+
linestyle='--', label='$\lambda_2$', alpha=0.5)
704+
plt.xlabel('t', size=18)
705+
plt.ylabel(r'$\frac{y_t}{y_{t-1}}$', size=18)
706+
plt.title(r'$\frac{y_t}{y_{t-1}}$ after Muting $\lambda_1$',
707+
size=13)
701708
plt.legend()
702709
703710
plt.subplot(1, 2, 2)
704-
plt.plot(ratios_λ2, label='Ratios $y_t / y_{t-1}$ after muting $\lambda_2$', color='orange')
705-
plt.axhline(y=Λ[0], color='green', linestyle='--', label='$\lambda_1$')
706-
plt.xlabel('Time')
707-
plt.ylabel('Ratio $y_t / y_{t-1}$')
708-
plt.title('Ratios after Muting $\lambda_2$')
711+
plt.plot(ratios_λ2,
712+
label=r'$\frac{y_t}{y_{t-1}}$', linewidth=3)
713+
plt.axhline(y=Λ[0], color='green',
714+
linestyle='--', label='$\lambda_1$', alpha=0.5)
715+
plt.xlabel('t', size=18)
716+
plt.ylabel(r'$\frac{y_t}{y_{t-1}}$', size=18)
717+
plt.title(r'$\frac{y_t}{y_{t-1}}$ after Muting $\lambda_2$',
718+
size=13)
709719
plt.legend()
710720
711721
plt.tight_layout()
712722
plt.show()
713723
```
714724
715-
Here we compare $V_{2,2} V_{1,2}^{-1}$ and $V_{2,1} V_{1,1}^{-1}$ with the roots we computed above
716-
717-
```{code-cell} ipython3
718-
np.round((V[1,1]/V[0,1],
719-
V[1,0]/V[0,0]), 8)
720-
```
721-
722725
## Concluding Remarks
723726
724727
This lecture sets the stage for many other applications of the **invariant subspace** methods.

0 commit comments

Comments
 (0)