Нелогичные параметры, возвращаемые scipy.curve_fit - PullRequest
0 голосов
/ 30 апреля 2018

Я моделирую шар, падающий через жидкость в Python, и подгоняю функцию модели к набору точек данных, используя коэффициенты демпфирования (a и b) и плотность жидкости, но подходящее значение для плотность жидкости возвращается отрицательной, и я понятия не имею, что не так в коде. Мой код ниже:

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
from scipy.optimize import curve_fit


##%%Parameters and constants 
m       = 0.1                               #mass of object in kg
g       = 9.81                              #acceleration due to gravity in m/s**2
rho     = 700                               #density of object in kg/m**3
v0      = 0                                 #velocity at t=0
y0      = 0                                 #position at t=0
V       = m / rho                           #volume in cubic meters
r       = ((3/4)*(V/np.pi))**(1/3)          #radius of the sphere
asample = 0.0001                            #sample value for a 
bsample = 0.0001                            #sample value for b

#%%Defining integrating function 
##fuction => y'' = g*(1-(rhof/rho))-((a/m)y'+(b/m)y'**2)
##   y' = v
##   v' = g*(1-rhof/rho)-((a/m)v+(b/m)v**2)

def sinkingball(Y, time, a, b, rhof):
    return [Y[1],(1/m)*(V*g*(rho-rhof)-a*Y[1]-b*(Y[1]**2))]

def balldepth(time, a, b, rhof):
    solutions = odeint(sinkingball, [y0,v0], time, args=(a, b, rhof))
    return solutions[:,0] 

time    = np.linspace(0,15,151)

# imported some experimental values and named the array data

a, b, rhof = curve_fit(balldepth, time, data, p0=(asample, bsample, 100))[0]
print(a,b,rhof) 

1 Ответ

0 голосов
/ 30 апреля 2018

Предоставление результата, который вы действительно получите, было бы полезно, и комментарий о том, что time не используется sinkingball(), заслуживает внимания.

Вам может пригодиться lmfit (https://lmfit.github.io/lmfit-py). Это обеспечивает более высокий уровень интерфейса для подбора кривой, который позволяет, помимо прочего, устанавливать границы для параметров, чтобы они могли оставаться физически чувствительными. Я думаю, что ваша проблема будет переводить с curve_fit на lmfit как:

from lmfit import Model
def balldepth(time, a, b, rhof):
    solutions = odeint(sinkingball, [y0,v0], time, args=(a, b, rhof))
    return solutions[:,0] 

# create a model based on the model function "balldepth"
ballmodel = Model(balldepth)

# create parameters, which will be named using the names of the
# function arguments, and provide initial values
params = bollmodel.make_params(a=0.001, b=0.001, rhof=100)

# you wanted rhof to **not** vary in the fit:
params['rhof'].vary = False

# set min/max values on `a` and `b`:
params['a'].min = 0
params['b'].min = 0

# run the fit
result = ballmodel.fit(data, params, time=time)

# print out full report of results
print(result.fit_report())

# get / print out best-fit parameters:
for parname, param in result.params.items():
    print("%s = %f +/- %f" % (parname, param.value, param.stderr))
...