MVC_boat_problem.py

# my love for math is driven by its beauty.
# but ultimately, i'm here to learn how to solve problems.
# sometimes, the best path forward is not the beautiful one -- and we have to learn how to follow that one too.

# what follows is the ugly solution, but frankly, a solution nonetheless.

import numpy as np
import matplotlib.pyplot as plt
from tqdm import tqdm

TIME_LIM = 500
TIME_STEP = 0.0001
EPSILON = 0.01

def calc_thrust(t):
    if (t >= 30):
        return 25
    return 25+np.log(-t+30)
    # return 30

def calc_current(y):
    return 10/(y+1)
    # return 20

def calc_landing(theta, plot_path=False):
    x, y, t = [0]*3

    if plot_path:
        xs = []
        ys = []

    while x < 2:
        if t > TIME_LIM:
            print("TLE!")
            return -1

        thrust = calc_thrust(t)

        thrust_y = (theta / 90) * thrust
        thrust_x = np.sqrt(thrust**2 - thrust_y**2)

        current = calc_current(y)

        x += thrust_x * TIME_STEP
        y += (thrust_y - current) * TIME_STEP

        t += TIME_STEP

        if plot_path:
            xs.append(x)
            ys.append(y)

    if abs(y-2) < EPSILON:
        print(f"theta: {theta}, y: {y}, t (min) {t*60}")
        if plot_path:
            plt.plot(xs, ys)
            plt.show()

    return y

theta = np.linspace(20,89,500)

landing = [calc_landing(t, plot_path=False) for t in tqdm(theta)]

plt.plot(theta, landing)
plt.ylim([-1, 6])
plt.show()

calc_landing(72.0099) #theta: 72.0099, y: 2.0004, t (min) 7.050