<jupyter_start><jupyter_text>Kepler problem<jupyter_code>import sympy
sympy.init_printing()

import numpy

import matplotlib.pyplot as plt

import sys
sys.path.insert(0, './code')

from gauss_legendre import gauss_legendre
from evaluate_functional import evaluate_functional
# state vector (symbolic)
q_x = sympy.Symbol('q_x')
q_y = sympy.Symbol('q_y')
p_x = sympy.Symbol('p_x')
p_y = sympy.Symbol('p_y')
x = [q_x, q_y, p_x, p_y]
# physical parameters of Hamiltonian function
from sympy.abc import m, G
params = {
    # mass of the planet (kg)
    m: 4.0,
    
    # gravitational constant (m**3/s**2)
    G: 10.0
}
# Hamiltonian function
# (kinetic energy T + potential energy V)
T = sympy.Rational(1,2) * (p_x**2 / m + p_y**2 / m)
V = - G * m / sympy.sqrt(q_x**2 + q_y**2)
H = T + V
H
# efforts
e = sympy.Matrix([H.diff(co) for co in x])
e
# structure matrix
J = sympy.Matrix([[0, 0, 1, 0],
                  [0, 0, 0, 1],
                  [-1, 0, 0, 0],
                  [0, -1, 0, 0]])
# time evolution
xdot = J @ e
xdot<jupyter_output><empty_output><jupyter_text>compute stationary orbit$F_G \overset{!}{=} F_{\omega}$$\frac{G \, m}{r^2} \overset{!}{=} m \, \omega^2 \, r$$\omega = \frac{G}{r^3}$<jupyter_code>r = 5.0
ω = numpy.sqrt(params[G]/r**3)
# initial condition
x_0 = numpy.array([0.0, r, params[m] * r * ω, 0.0])

# duration of simulation
t_f = 260.0

# time step
dt = 0.1

# order of collocation polynomial
s = 2
%time time, solution = gauss_legendre(x, xdot, x_0, t_f, dt, s, params=params)
# plot trajectory (x-y-plane)
fig, ax = plt.subplots(dpi=200)
ax.set_aspect('equal')
ax.plot(solution[:, 0], solution[:, 1]);
# relative error of energy conservation
energy = evaluate_functional(x, H, solution, params=params)
abs(energy[0] - energy[-1]) / energy[0]<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.16916250693288962, 'fraction_numerical': 0.03161397670549085, 'mean_word_length': 3.7724867724867726, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '544365', 'license': 'MIT', 'n_tokens_mistral': 764, 'n_tokens_neox': 725, 'n_words': 188, 'path': 'kepler.ipynb', 'repo_name': 'MarkusLohmayer/master-thesis-code', 'size': '1803'}