2016-01-08 8 views
1

Я пытаюсь написать код, который будет отображать симуляцию шара, выпавшего с высоты h, и составить график положения во времени с использованием кинематических уравнений y = y_0 Моего код заключается в следующем:Шаг падения шага по времени с помощью python

из matplotlib.pylab импортного шоу, xlabel, ylabel, разброс, сюжет от Numpy импорта пустого

def drop(): 

    """ 
    This function calculates and creates arrays for the velocity at eac time interval as well as the position and plots it. Assuming no drag. 
    """ 
    #Define the constants in the problem 
    h_0 = 10 
    g = -9.8 #gravitational constant N 
    dt = 0.1 #timestep 

    #Now need to create arrays to hold the positins, time, and velocities 
    vel = empty(1000,float) 
    time = empty(1000,float) 
    y = empty(1000,float) 

    time[0] = 0 
    vel[0] = 0 
    y[0] = h_0 

    #code for the kinematic equations for the calculation of time, velocity and position 
    for i in range of (1000-1): 
     time[i+1] = time[i] + dt 
     vel[i+1] = vel[i] + (g * dt) 
     y[i+1] = time[i] + (vel[i+1] * dt) 

     if y[i] > 0: 
     #ensures that the graph will not keep going when the ball hits the ground 
      break 


    plot(time,y, '.') 
    xlabel("Time(s)") 
    ylabel("Position") 
    show() 

Однако мой график участков три точки по одному в каждом углу графика когда это он должен выглядеть как кривая, и мой график меняется каждый раз, когда он не должен, так как нет ne переменных меняется

ответ

1

Хорошо, давайте попробуем синтаксическую ошибку. for i in range of (1000-1) на самом деле for i in range(1000-1), но я предполагаю, что это была опечатка с вашей стороны, так как вы могли запустить код.

Теперь ваши уравнения движения ошибочны.

y[i+1] = time[i] + (vel[i+1] * dt) 

# should be 
y[i+1] = y[i] + (vel[i] * dt) 

Ваше состояние для выхода из моделирования также является ошибочным.

if y[i] > 0: 

# you have to stop when the height becomes until negative 
if y[i+1] < 0: 

Ваши ошибки до сих пор означает, что вы будете выйти из цикла после одной итерации, имея фактически не измененную ваш y массив. Последняя проблема здесь. numpy.empty() создает массив без инициализации значений. Это означает, что исходные значения будут независимо от того, что находится в памяти в этой точке. Если вы напечатаете y после разрыва цикла, вы можете заметить, что большинство значений равно 0, в то время как некоторые из них очень малы, но не близки к 0, например. 3.18377034e-308. Поскольку они являются самыми высокими значениями в вашем массиве, они будут масштабировать ваш график до их диапазона. Но поскольку они являются произвольными значениями, каждый раз, когда вы запускаете код, он будет производить разные числа.

У вас есть два варианта исправить это. Используйте либо numpy.zeros(), либо нарисуйте только первые значения y[:i], которые являются точкой в ​​петле, которую вы ломаете для удара по земле.


Поскольку мы имеем аналитическое решение уравнений в ваших проблемах, вы можете избавиться от петель и vectorise все с массивами. Мы можем решить уравнение смещения относительно t (квадратичное), чтобы выяснить, когда мы ударимся по земле. Затем мы инициализируем массив времени и используем его для вычисления смещения (скорость необязательна).

import numpy as np 
import matplotlib.pyplot as plt 

def time_to_hit_ground(y0, v0, a): 
    discriminant = np.sqrt(v0**2 - 2*a*y0) 
    t1 = (-v0 - discriminant)/a 
    t2 = (-v0 + discriminant)/a 
    if t1 >=0: 
     return t1 
    return t2 

def drop(y0, v0=0.0, dt=0.1, g=-9.8): 
    if y0 < 0: 
     print('Object is underground.') 
     return 
    # if you also allow the user to change `dt` and `g` from the arguments, 
    # you want to check they have the correct sign. 

    t_ground = time_to_hit_ground(y0, v0, g) 

    t = np.arange(0, t_ground+dt, dt) 
    v = v0 + g * t 
    y = y0 + v0 * t + g * t**2/2. 

    plt.plot(t, y, '.') 
    plt.axvline(t_ground, color='g') 
    plt.axhline(0, color='g') 
    plt.xlabel('Time (s)') 
    plt.ylabel('Height (m)') 
    plt.show() 

 Смежные вопросы

  • Нет связанных вопросов^_^