Подходящие дифференциальные уравнения: Curve_fit сходится к локальным минимумам - PullRequest
0 голосов
/ 07 февраля 2012

Я пытаюсь подогнать дифференциальное уравнение ay '+ by' '= 0 к кривой, изменяя a и b. Следующий код не работает. Проблема с curve_fit, похоже, заключается в том, что отсутствие первоначальной догадки приводит к неудаче при подгонке. Я также попробовал наименьшее количество. Может ли кто-нибудь предложить мне другие способы подбора такого дифференциального уравнения? Если у меня нет правильного предположения, то получится, что Curve_fit не удастся!

from scipy.integrate import odeint
from scipy.optimize import curve_fit
from numpy import linspace, random, array

time = linspace(0.0,10.0,100)
def deriv(time,a,b): 
    dy=lambda y,t : array([ y[1], a*y[0]+b*y[1] ])
    yinit = array([0.0005,0.2]) # initial values
    Y=odeint(dy,yinit,time)
    return Y[:,0]

z = deriv(time, 2, 0.1)
zn = z + 0.1*random.normal(size=len(time))

popt, pcov = curve_fit(deriv, time, zn)
print popt  # it only outputs the initial values of a, b!

Ответы [ 2 ]

2 голосов
/ 07 февраля 2012

Давайте перепишем уравнение:

ay' + by''=0
y'' = -a/b*y'

Таким образом, это уравнение может быть представлено следующим образом

dy/dt = y'
d(y')/dt = -a/b*y'

Код в Python 2.7:

from scipy.integrate import odeint
from pylab import *

a = -2
b = -0.1

def deriv(Y,t):
    '''Get the derivatives of Y at the time moment t
Y = [y, y' ]'''
    return array([ Y[1], -a/b*Y[1] ])

time = linspace(0.0,1.0,1000)
yinit = array([0.0005,0.2]) # initial values
y = odeint(deriv,yinit,time)
figure()
plot(time,y[:,0])
xlabel('t')
ylabel('y')
show()

Вы можете сравнить полученные графики с графиками в WolframAlpha

0 голосов
/ 09 февраля 2012

Если ваша проблема заключается в том, что начальные догадки по умолчанию считываются, прочитайте документацию curve_fit , чтобы узнать, как задать их вручную, указав параметр p0.Например, curve_fit(deriv, time, zn, p0=(12, 0.23)), если вы хотите, чтобы a=12 и b=0.23 были начальным предположением.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...