Можете ли вы решить оду с переменными константами? - PullRequest
Новая
       73

Можете ли вы решить оду с переменными константами?

0 голосов
/ 07 сентября 2018

У меня есть универ, где я должен создать модель для лунного спуска. Модель начинается с высоты и опускается с двух двигателей, которые расположены под углом. Это часть кода, с которой у меня проблемы: 

F1 = f*np.array([1, 1, 0.5, 0, 0, 1, 1]) #f-thrust in newtons
F2 = f*np.array([1, 1, 0.5, 0, 0, 1, 1]) #F1 and F2 regime of engines for each time step(percentage of thrust)

t = np.arange(0, np.size(F1), dt) #time domain

#a,b - angles of engines according to global coordinate system
#s0x, v0x, s0y, v0y, fi0, omg0 - initial values
#m - mass
def model_engine(u, t):
    x, vx, y, vy, phi, w = u   
    d1 = vx
    d2 = (1/(m))*(F1*np.cos((a)) - F2*np.cos((b)))
    d3 = vy
    d4 = (1/(m))*(F1*np.sin((a)) + F2*np.sin((b))) + g 
    d5 = w
    d6 = (F1*np.cos(a)*(h/2) - F1*np.sin(a)*(w/2) + F2*np.sin(b)*(w/2) - 
             F2*np.cos(b)*(h/2))/((1/12)*(m)*h**2)
    return np.array([d1, d2, d3, d4, d5, d6])

U = odeint(model_engine, [s0x, v0x, s0y, v0y, fi0, omg0], t)

  1. Я не знаю, как реализовать разные F1 для каждого временного шага в определении? Я получаю ValueError: setting an array element with a sequence. То же самое с переменной массой, это зависит от режима работы двигателей.

  2. Решение дает угол космического корабля (d6). Уравнения, написанные таким образом, всегда дают a и b как 45 градусов, что, очевидно, не правильно, потому что космический корабль вращается. Как я могу взять решение (угол) из предыдущего шага и поместить его в текущий шаг? Смотрите фото для уточнения

photo

1 Ответ

0 голосов
/ 07 сентября 2018

Вот метод для первого вопроса: Как реализовать переменный коэффициент? Предполагается, что тяга является ступенчатой ​​функцией только t. Он получается путем интерполяции заданных команд тяги в соответствующее время. Функция интерполяции передается функции model_engine в качестве дополнительного аргумента.

Модель упрощена для вертикальной платформы, имеющей постоянную массу. 

import matplotlib.pylab as plt

import numpy as np
from scipy.interpolate import interp1d
from scipy.integrate import odeint

from scipy.constants import g

def model_engine(u, t, thrust, mass):
    y, vy = u 

    dydt = vy
    dvy_dt = 1/mass * thrust(t) - g 

    return np.array([dydt, dvy_dt])


# Thrust definition:
F_values = g*np.array([1, 1, 0.5, 0, 0, 2, 2, 1.5, 1]) # f-thrust in newtons
timesteps = np.linspace(0, 10, np.size(F_values)) # timesteps used to defined the thrust

thrust = interp1d(timesteps, F_values, kind='previous', fill_value='extrapolate')

# Solve
s0y, v0y = 102, 0
m0 = 1
t_span = np.linspace(0, 1.2*timesteps[-1], 111)
U = odeint(model_engine, [s0y, v0y], t_span, args=(thrust, m0))

# Graphs
plt.figure();
plt.plot(t_span, U[:, 0], label='altitude')
plt.axhline(y=0, color='black', label='moon surface');
plt.xlabel('time'); plt.legend();

plt.figure();
plt.plot(t_span, thrust(t_span), label='thrust');
plt.xlabel('time'); plt.legend();
...