问题描述
我正在尝试使用 Symfit 解决并拟合这些耦合的 Differential equations。我试过按照 ODE fitting 的文档,查看 this answer 到一个类似的问题并搜索 the error in other contexts,但没有任何运气。我的代码应该用于实验数据,但作为开始,我刚刚创建了一些类似的数据。相关代码为:
#Imports
import numpy as np
import matplotlib.pyplot as plt
from symfit import variables,parameters,Parameter,Fit,D,ODEModel,cos,sin,atan,asin,sqrt
#Constants
m = 0.01 #Unit: kg - Mass of object
mt = 0.5 * m #Unit: kg - Related to moment of inertia (\tilde{m] in the equations)
R = 0.1 #Unit: m - Radius of object
g = 9.82 #Unit: m/s^2 - gravitational acceleration
#Variables and parameters
x_1,x_2,x_3,x_4,t = variables("x_1,t") #x-position,y-position,x-veLocity,y-veLocity,time
p = Parameter("p",min=1,max=4)
c,d = parameters("c,d") #d is \gamma in the equations
#Generate proxy data
t_vals = np.linspace(0,10,100) #Choosen time interval
x1_data = 8*np.sin(0.15*t_vals)+7*t_vals #x-data proxy
x2_data = 10*np.sin(0.15* t_vals)-2*t_vals #y-data proxy
x3_data = 0.15*8*np.cos(0.15*t_vals)-7 #vx-data proxy
x4_data = 0.15*10*np.cos(t_vals)-2 #vy-data proxy
#Model
model_dict = {
D(x_1,t) : x_3,D(x_2,t) : x_4,D(x_3,t) : -d/mt * (x_3**2 + x_4**2)**((p-1)/2) * x_3 + c*(x_3**2+x_4**2)/(mt*R) * cos(atan(x_1/x_2) + asin(R/sqrt(x_1**2+x_2**2))),D(x_4,t) : -d/mt * (x_3**2 + x_4**2)**((p-1)/2) * x_4 + c*(x_3**2+x_4**2)/(mt*R) * sin(atan(x_1/x_2) + asin(R/sqrt(x_1**2+x_2**2)))+m/mt * g,}
ode_model = ODEModel(model_dict,initial={t: 0.0,x_1: 0,x_2: 3,x_3: 0,x_4: 1})
fit = Fit(ode_model,t=tdata,x_1=x1_data,x_2=x2_data,x_3=x3_data,x_4=x4_data)
fit_result = fit.execute()
但我收到以下错误:
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-8-ef065dc7aeb9> in <module>
27
28 fit = Fit(ode_model,x_4=x4_data)
---> 29 fit_result = fit.execute()
30
31 t_vec = np.linspace(0,1000)
[...]
ValueError: operands Could not be broadcast together with remapped shapes [original->remapped]: (5,) and requested shape (100,)
我在使用 Symfit 方面相当缺乏经验,但我认为错误所指的两个数组是五个变量和每个变量的 100 个数据点。我不知道是什么导致了这个错误,除了可能缺乏初始猜测和错误生成的数据。 问题只是我使用了错误的时间数组
编辑:正如用户 mikuszefski 指出的那样,我在应该使用 tdata
的地方使用了 t_vals
。 这解决了我最初的错误。同一用户还正确预测了数值问题,因为一些拟合值的数量级为 10^300。一旦找到解决方案,我将编辑我的帖子。
我还编写了我的导入并定义了 g
。
解决方法
暂无找到可以解决该程序问题的有效方法,小编努力寻找整理中!
如果你已经找到好的解决方法,欢迎将解决方案带上本链接一起发送给小编。
小编邮箱:dio#foxmail.com (将#修改为@)