【发布时间】:2019-02-14 17:15:56
【问题描述】:
我已经记录了固体表面五个位置的实验温度。在每个时间步,我想将这些读数拟合到由我的函数定义的理论曲线:Temp_Function_JLT(X,h)。
X 是一个多维数组,包括 x_coordinates 以及时间、初始温度和材料属性(所有独立变量)。 “h”是传热系数,出于本练习的目的,我试图对其进行优化(暂时搁置物理学。)
这是我的温度函数的定义:
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import pickle
import scipy.optimize as opt
from scipy.special import erfc
def Temp_Function_JLT(X,ht):
# Work around the fact that only one independent variable can be passed to optimize.curve_fit
x,t,T0,q,alpha,rho,c,k = X
term_a = q/ht
term_b = erfc(x/np.sqrt(4*alpha*t))
term_c = np.exp(((ht*x)/(np.sqrt(alpha)*np.sqrt(k*rho*c)))+((ht**2)/(k*rho*c)))
term_d = erfc((ht*np.sqrt(t))/(np.sqrt(k*rho*c)) + (x/np.sqrt(4*alpha*t)))
Temperature = (term_a * (term_b - term_c * term_d)) + T0 - 273
return Temperature
该功能有效。我可以使用一些初始参数运行它并获得合理的值。对于这个问题更重要的是,如果我用以下数据调用它:
t = 1
x_test = np.linspace(0.004,0.02,5) # TC locations
time_test = range(1,180,30)
T0_test = 25 + 273
q_test = 20000
h_test = 10
我将获得一个 numpy 数组作为形状 (1,) 的解,它给出 np.ndim 1 的答案(这已在以下先前问题中提到:
Fitting a vector function with curve_fit in Scipy
Fitting a 2D Gaussian function using scipy.optimize.curve_fit - ValueError and minpack.error
当我调用 opt.curve_fit() 时出现问题。 indepth_temperatures 是一个列表,其中包含每个测试作为一个数组。我迭代它(迭代每个测试),然后根据以下代码对每一行(每个时间步)执行拟合:
for i,test in enumerate(indepth_temperatures):
# Iterate over every row
for j,row in enumerate(test):
# Define tuple that contains all independent variables
X = (TC_depth,
times[i][j],
T0_temperatures[i] + 273,
20000,
pmma_alpha,
pmma_rho,
pmma_c,
pmma_k)
print(Temp_Function_JLT(X,h0))
print(row)
print('---')
# Call function to optimize curve fit on h
popt, pcov = opt.curve_fit(Temp_Function_JLT,X,row,h0)
print(popt)
对于第一次迭代,我得到以下结果:
[23.2034 23.2034 23.2034 23.2034 23.2034] # comes from print(Temp_Function_JLT(X,h0))
[23.937 22.619 22.59 24.884 21.987000000000002] # comes from print(row)
随后出现此错误:
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'safe'
---------------------------------------------------------------------------
error Traceback (most recent call last)
<ipython-input-67-9c4545fd257b> in <module>()
22 print('---')
23 # Call function to optimize curve fit on h
---> 24 popt, pcov = opt.curve_fit(Temp_Function_JLT,X,row,h0)
25 print(popt)
~\AppData\Local\Continuum\anaconda2\envs\py36\lib\site-packages\scipy\optimize\minpack.py in curve_fit(f, xdata, ydata, p0, sigma, absolute_sigma, check_finite, bounds, method, jac, **kwargs)
749 # Remove full_output from kwargs, otherwise we're passing it in twice.
750 return_full = kwargs.pop('full_output', False)
--> 751 res = leastsq(func, p0, Dfun=jac, full_output=1, **kwargs)
752 popt, pcov, infodict, errmsg, ier = res
753 cost = np.sum(infodict['fvec'] ** 2)
~\AppData\Local\Continuum\anaconda2\envs\py36\lib\site-packages\scipy\optimize\minpack.py in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
392 with _MINPACK_LOCK:
393 retval = _minpack._lmdif(func, x0, args, full_output, ftol, xtol,
--> 394 gtol, maxfev, epsfcn, factor, diag)
395 else:
396 if col_deriv:
error: Result from function call is not a proper array of floats.
我尝试从我的函数 np.ravel(Temperature) 或 Temperature.flatten() 返回,但没有成功。错误仍然存在,我无法弄清楚它为什么在那里。正如我所提到的,我已经检查了我的函数返回的维度,它是一个一维数组。
任何帮助将不胜感激!
更新:我意识到很难复制这段代码,所以这是一个简化版本:
Temp_Function_JLT(X,h0):保持不变。
pmma_rho = 1200 # kg/m3
pmma_c = 1500 # J/kgK
pmma_k = 0.16 # W/mK
pmma_alpha = pmma_k/(pmma_rho*pmma_c)
x_test = np.linspace(0.004,0.02,5) # TC locations
t = 1
T0_test = 25 + 273
q_test = 20000
h_test = 10
X = (x_test,t,T0_test,q_test,pmma_alpha,pmma_rho,pmma_c,pmma_k)
y_data = [23.937 22.619 22.59 24.884 21.987000000000002]
opt.curve_fit(Temp_Function_JLT, X, y_data, h_test)
【问题讨论】:
-
我无法运行您的代码,因为您没有指定
indepth_temperatures是什么。此外,还有很多其他未知变量:h0,以及所有其他以pmma___开头的变量 -
感谢您的回答! pmma_rho = 1200 pmma_c = 1500 pmma_k = 0.16 pmma_alpha = pmma_k/(pmma_rho * pmma_c) h0 = 10 indepth_temperatures 是一个列表,其中包含每次测试的温度读数(总共 20 次测试)。每个测试包含每个时间步长(大约 1200 个时间步长)的五个温度读数(我的 x_data)。这就是为什么我在我的结果中打印第一行,这是第一个时间步中曲线拟合的“y_data”。换句话说,我的曲线拟合使用:Temp_Function_JLT 作为函数,行(在我的结果中打印)作为 y_data,X 作为 x_data。我会简化我的帖子。
-
嗯嗯对。所以,出于某种原因,我的行 (y_data) 是一个 dtype = object 的 numpy 数组。我认为这是问题的根源。如果我打印行,这就是我得到的:array([23.937, 22.619, 22.59, 24.884, 21.987000000000002], dtype=object)
-
抱歉,
indepth_temperatures仍然没有为我定义,尽管您解释了它包含的内容。为了帮助您的人们节省时间,最好定义所有变量并将上面的代码复制粘贴到一个新的 python 文件中并执行它。 仅如果您没有收到任何未定义的变量错误,请发布代码。否则我们会一而再再而三地问你这个变量是缺失的还是那个变量是缺失的
标签: python numpy scipy curve-fitting