RuntimeWarning: недопустимое значение, возникающее при делении

Мне нужно сделать программу, используя метод Эйлера для модели "ball in a spring

from pylab import*
from math import*
m=0.1
Lo=1
tt=30
k=200
t=20
g=9.81
dt=0.01
n=int((ceil(t/dt)))
km=k/m
r0=[-5,5*sqrt(3)]
v0=[-5,5*sqrt(3)]
a=zeros((n,2))
r=zeros((n,2))
v=zeros((n,2))
t=zeros((n,2))
r[1,:]=r0
v[1,:]=v0
for i in range(n-1):
    rr=dot(r[i,:],r[i,:])**0.5
    a=-g+km*cos(tt)*(rr-L0)*r[i,:]/rr
    v[i+1,:]=v[i,:]+a*dt
    r[i+1,:]=r[i,:]+v[i+1,:]*dt
    t[i+1]=t[i]+dt

    #print norm(r[i,:])

plot(r[:,0],r[:,1])
xlim(-100,100)
ylim(-100,100)
xlabel('x [m]')
ylabel('y [m]')

show()

Я продолжаю получать эту ошибку:

a=-g+km*cos(tt)*(rr-L0)*r[i,:]/rr
RuntimeWarning: invalid value encountered in divide

Я не могу понять, что не так с кодом?

Ответы

Ответ 1

Я думаю, что ваш код пытается "делить на ноль" или "делить на NaN". Если вы знаете об этом и не хотите, чтобы это вас беспокоило, вы можете попробовать:

import numpy as np
np.seterr(divide='ignore', invalid='ignore')

Подробнее см.:

http://docs.scipy.org/doc/numpy/reference/generated/numpy.seterr.html

Ответ 2

Индексирование Python начинается с 0 (вместо 1), поэтому ваше задание "r [1,:] = r0" определяет элемент second (т.е. индекс 1) из r и оставляет первый ( index 0) как пара нулей. Первое значение я в вашем цикле for равно 0, поэтому rr получает квадратный корень из точечного произведения первой записи в r с самим собой (что равно 0), а деление на rr в последующей строке вызывает ошибку.

Ответ 3

Вы делите на rr, который может быть 0.0. Убедитесь, что rr равен нулю и делает что-то разумное, кроме использования его в знаменателе.

Ответ 4

Чтобы предотвратить деление на ноль, вы можете предварительно инициализировать вывод 'out', где происходит ошибка div0, например, np.where не обрезает его, так как полная строка вычисляется независимо от условия.

пример с предварительной инициализацией:

a = np.arange(10).reshape(2,5)
a[1,3] = 0
print(a)    #[[0 1 2 3 4], [5 6 7 0 9]]
a[0]/a[1]   # errors at 3/0
out = np.ones( (5) )  #preinit
np.divide(a[0],a[1], out=out, where=a[1]!=0) #only divide nonzeros else 1

Ответ 5

Это можно исправить с помощью

a = -g + km * cos (tt) * (rr-Lo) * r [i,:] / (rr + (1)) #