Skip to content

Commit

Permalink
Bump Dependencies and Update Equation (#28)
Browse files Browse the repository at this point in the history
* Bump dependencies and fix resulting errors

* Improve equation formatting and add spelling words

* More bumps!

* Fix errors
  • Loading branch information
bryanwweber committed Jun 16, 2024
1 parent 3e51bd7 commit b2b2a74
Show file tree
Hide file tree
Showing 14 changed files with 860 additions and 730 deletions.
64 changes: 36 additions & 28 deletions .vscode/settings.json
Original file line number Diff line number Diff line change
@@ -1,31 +1,39 @@
{
"cSpell.words": [
"apoapsis",
"deorbit",
"Hohmann",
"Kármán",
"periapsis",
"semimajor",
"Storrs"
],
"cSpell.ignoreWords": ["devcontainer"],
"cSpell.languageSettings": [
{
// use with Markdown files
"languageId": "markdown",
"dictionaries": ["latex"],
// Exclude code
// "ignoreRegExpList": [
// "\\$.*\\$",
// "\\$\\$([\\n\\(\\)]|[^(\\$\\$)])*\\$\\$",
// "`.*`",
// "```([\\n\\(\\)]|[^(\\$\\$)])*```",
// "\\{\\w+\\}`.+`",
// ":::\\{math\\}([\\n\\(\\)]|[^(\\$\\$)])*:::",
// "---([\\n\\(\\)]|[^(\\$\\$)])*---"
// ]
}
"cSpell.words": [
"apoapse",
"apoapsis",
"deorbit",
"Hohmann",
"Kármán",
"Laguerre",
"periapsis",
"Prussing",
"semimajor",
"Semiminor",
"Storrs",
"Stumpff"
],
"cSpell.ignoreWords": [
"devcontainer"
],
"cSpell.languageSettings": [
{
// use with Markdown files
"languageId": "markdown",
"dictionaries": [
"latex"
],
"editor.formatOnSave": true

// Exclude code
// "ignoreRegExpList": [
// "\\$.*\\$",
// "\\$\\$([\\n\\(\\)]|[^(\\$\\$)])*\\$\\$",
// "`.*`",
// "```([\\n\\(\\)]|[^(\\$\\$)])*```",
// "\\{\\w+\\}`.+`",
// ":::\\{math\\}([\\n\\(\\)]|[^(\\$\\$)])*:::",
// "---([\\n\\(\\)]|[^(\\$\\$)])*---"
// ]
}
],
"editor.formatOnSave": true
}
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,7 @@
warnings.simplefilter("ignore")
from myst_nb import glue

np.set_printoptions(legacy="1.25")
glue("orbital-elements-radius", r, display=False)
glue("orbital-elements-velocity", v, display=False)
glue("orbital-elements-v_r", v_r, display=False)
Expand Down
1 change: 1 addition & 0 deletions interplanetary-maneuvers/planetary-arrival-flyby.md
Original file line number Diff line number Diff line change
Expand Up @@ -250,6 +250,7 @@ h_1 = np.sqrt(mu_sun * R_Neptune * (1 - e_1))
from functools import partial
from myst_nb import glue as myst_glue
glue = partial(myst_glue, display=False)
np.set_printoptions(legacy="1.25")
glue("interplanetary-flyby-e_1", e_1)
glue("interplanetary-flyby-h_1", h_1)
```
Expand Down
1 change: 1 addition & 0 deletions orbital-maneuvers/comparison-of-bielliptic-and-hohmann.md
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,7 @@ v_3 = np.sqrt(mu / r_3) # km/s

```{code-cell} ipython3
:tags: [remove-cell]
np.set_printoptions(legacy="1.25")
glue("moon_circle_velocity", v_1)
glue("leo_circle_velocity", v_3)
```
Expand Down
1 change: 1 addition & 0 deletions orbital-maneuvers/plane-change-example.md
Original file line number Diff line number Diff line change
Expand Up @@ -66,6 +66,7 @@ Delta_v_1 = Delta_v_1_LEO + Delta_v_1_GEO

```{code-cell} ipython3
:tags: [remove-cell]
np.set_printoptions(legacy="1.25")
from functools import partial
from myst_nb import glue as myst_glue
glue = partial(myst_glue, display=False)
Expand Down
2 changes: 2 additions & 0 deletions orbital-maneuvers/single-impulse-example.md
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,8 @@ from matplotlib.patches import Circle, Arc
import numpy as np
from myst_nb import glue
np.set_printoptions(legacy="1.25")
R_E = 6378 # km
orbit_radius = 1000 # km
mu = 398_600 # km**2/s**3
Expand Down
1,479 changes: 795 additions & 684 deletions pdm.lock

Large diffs are not rendered by default.

4 changes: 3 additions & 1 deletion pyproject.toml
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ dependencies = [
"astropy>=6.0.0",
"docutils>=0.20.1",
"ipympl>=0.9.3",
"jupyter-book>=1.0.0",
"jupyter-book>=1.0",
"matplotlib>=3.8.3",
"myst-nb-bokeh>=2024.1.0",
"myst-nb>=1.0.0",
Expand All @@ -22,6 +22,8 @@ requires-python = ">=3.10,<3.13"
name = "orbital-mechanics-notes"
authors = [{ name = "Bryan Weber", email = "[email protected]" }]
license = { text = "CC-BY-SA-4.0" }
description = "A book for Orbital Mechanics"
version = "0.0.0"

[tool.pdm]
distribution = false
Expand Down
14 changes: 7 additions & 7 deletions the-n-body-problem/lagrange-points.md
Original file line number Diff line number Diff line change
Expand Up @@ -294,8 +294,8 @@ def init():
def animate(pi_2):
m_1 = -pi_2
m_2 = 1 - pi_2
m1_line.set_data(m_1, 0)
m2_line.set_data(m_2, 0)
m1_line.set_data([m_1], [0])
m2_line.set_data([m_2], [0])
x_2 = m_2 * cos
y_2 = m_2 * sin
x_1 = m_1 * cos
Expand All @@ -306,11 +306,11 @@ def animate(pi_2):
L_1 = newton(func=collinear_lagrange, x0=0, args=(pi_2,))
L_3 = newton(func=collinear_lagrange, x0=-1, args=(pi_2,))
L_4 = L_5 = 0.5 - pi_2
L1_line.set_data(L_1, 0)
L2_line.set_data(L_2, 0)
L3_line.set_data(L_3, 0)
L4_line.set_data(L_4, np.sqrt(3) / 2)
L5_line.set_data(L_5, -np.sqrt(3) / 2)
L1_line.set_data([L_1], [0])
L2_line.set_data([L_2], [0])
L3_line.set_data([L_3], [0])
L4_line.set_data([L_4], [np.sqrt(3) / 2])
L5_line.set_data([L_5], [-np.sqrt(3) / 2])
equil.set_data([m_1, L_4, m_2, L_5, m_1], [0, np.sqrt(3) / 2, 0, -np.sqrt(3) / 2, 0])
ann.set_text(fr"$\pi_2$ = {pi_2:.4G}")
Expand Down
2 changes: 1 addition & 1 deletion the-orbit-equation/hyperbolic_trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -59,7 +59,7 @@ def init():

def animate(t):
r = a * (e**2 - 1) / (1 + e * np.cos(t))
point.set_data(-a - r_p + r * np.cos(t), r * np.sin(t))
point.set_data([-a - r_p + r * np.cos(t)], [r * np.sin(t)])
ann.set_text(rf"$\nu$ = {np.degrees(t):.2F}°")
return (point, ann)

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,6 @@ e = \frac{r_a - r_p}{r_a + r_p}
:::

```{code-cell} ipython3
# %matplotlib notebook
import numpy as np
from scipy.optimize import newton
Expand All @@ -51,6 +50,7 @@ e = (r_a - r_p)/(r_a + r_p)
:tags: [remove-cell]
from functools import partial
from myst_nb import glue as myst_glue
np.set_printoptions(legacy="1.25")
glue = partial(myst_glue, display=False)
glue("ellipse-time-since-periapsis-e", e)
```
Expand Down Expand Up @@ -149,7 +149,7 @@ def kepler(E, M_e, e):
def d_kepler_d_E(E, M_e, e):
"""The derivative of Kepler's equation, to be used in a Newton solver.
Note that the argument M_e is unused, but must be present so the function
arguments are consistent with the kepler function.
"""
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,7 @@ from matplotlib.patches import Ellipse, Circle, Arc, Rectangle
import numpy as np
glue = partial(myst_glue, display=False)
np.set_printoptions(legacy="1.25")
mu = 3.986004418E5 # km**3/s**2
R_E = 6378 # km
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,8 @@ e = h**2 / (r_p * mu) - 1
:tags: [remove-cell]
from functools import partial
from myst_nb import glue as myst_glue
np.set_printoptions(legacy="1.25")
glue = partial(myst_glue, display=False)
glue("hyperbolic-time-since-perigee-h", h)
Expand Down Expand Up @@ -117,7 +119,7 @@ def kepler(F, M_h, e):
def d_kepler_d_F(F, M_h, e):
"""The derivative of Kepler's equation, to be used in a Newton solver.
Note that the argument M_h is unused, but must be present so the function
arguments are consistent with the kepler function.
"""
Expand All @@ -129,7 +131,7 @@ F_2 = newton(func=kepler, fprime=d_kepler_d_F, x0=np.pi, args=(M_h2, e))
With this value for $F$, we can calculate the value for $\nu$. To avoid quadrant ambiguity problems, we will use Eq. {eq}`eq:eccentric-anomaly-true-anomaly-hyperbola`.

```{code-cell} ipython3
sqrt_e = np.sqrt((e + 1) / (e - 1))
sqrt_e = np.sqrt((e + 1) / (e - 1))
nu_2 = (2 * np.arctan(sqrt_e * np.tanh(F_2 / 2))) % (2 * np.pi)
```

Expand Down Expand Up @@ -181,12 +183,12 @@ function kepler
v_p = 15; % km/s
h = r_p * v_p;
e = h^2 / (mu * r_p) - 1;
nu_1 = deg2rad(100);
F_1 = 2 * atanh(sqrt((e - 1) / (e + 1)) * tan(nu_1 / 2));
M_h1 = e * sinh(F_1) - F_1;
t_1 = h^3 / mu^2 * 1 / (e^2 - 1)^(3 / 2) * M_h1;
t_2 = t_1 + 3 * 3600;
M_h2 = mu^2 / h^3 * (e^2 - 1)^(3 / 2) * t_2;
Expand All @@ -198,7 +200,7 @@ function kepler
t2 = 2 * atan(sqrt((e + 1) / (e - 1)) * tanh(F_2 / 2));
nu_2 = mod(t2, 2 * pi);
disp(rad2deg(nu_2))
r_2 = h^2 / mu * 1 / (1 + e * cos(nu_2));
v_perp = h / r_2;
v_r = mu / h * e * sin(nu_2);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -306,13 +306,13 @@ Alternatively, a more refined estimate can be determined using a secant estimate
where $f(\chi^{+})$ is the solution of Kepler's equation with the $\chi^{+}$ value.

::::{note}
Prussing and Conway {cite}`Prussing2013`, citing Conway {cite}`Conway1986` suggest that faster convergence in the solution of Kepler's equation can be achieved by using the **Laguerre algorithm**, rather than Newton's algorithm. Another advantage of the Laguerre algorithm is that it is relatively insensitive to the value of the initial guess.
Prussing and Conway {cite}`Prussing2013`, citing Conway {cite}`Conway1986` suggest that faster convergence in the solution of Kepler's equation can be achieved by using the [**Laguerre algorithm**](https://en.wikipedia.org/wiki/Laguerre%27s_method), rather than Newton's algorithm. Another advantage of the Laguerre algorithm is that it is relatively insensitive to the value of the initial guess.

The Laguerre algorithm can be implemented as:

:::{math}
:label:
\chi_{i + 1} = \chi_{i} - \frac{n f(\chi_i)}{f'(\chi_i) \pm \left[\left(n - 1\right)^2 \left(f'(\chi_i)\right)^2 - n\left(n - 1\right) f(\chi_i)f''(\chi_i)\right]^{1/2}}
\chi_{i + 1} = \chi_{i} - \frac{n f(\chi_i)}{f'(\chi_i) \pm \sqrt{\left(n - 1\right)^2 \left[f'(\chi_i)\right]^2 - n\left(n - 1\right) f(\chi_i)f''(\chi_i)}}
:::

The sign ambiguity in the denominator is determined by taking the sign of the numerical value of $f'(\chi_i)$. In addition, the solution is relatively insensitive to the choice of the value of $n$, which is an integer constant. It seems as though $n = 5$ is a reasonable value. Choosing $n = 1$ gives the standard Newton's algorithm.
Expand Down

0 comments on commit b2b2a74

Please sign in to comment.