fmin_slsqp returns 初始猜测找到三次样条的最小值
fmin_slsqp returns initial guess finding the minimum of cubic spline
我正在寻找自然三次样条的最小值。我编写了以下代码来查找自然三次样条。 (我已经得到了测试数据,并确认这个方法是正确的。)现在我不知道如何找到这个函数的最小值。
这是数据
xdata = np.linspace(0.25, 2, 8)
ydata = 10**(-12) * np.array([1,2,1,2,3,1,1,2])
这是函数
import scipy as sp
import numpy as np
import math
from numpy.linalg import inv
from scipy.optimize import fmin_slsqp
from scipy.optimize import minimize, rosen, rosen_der
def phi(x, xd,yd):
n = len(xd)
h = np.array(xd[1:n] - xd[0:n-1])
f = np.divide(yd[1:n] - yd[0:(n-1)],h)
q = [0]*(n-2)
for i in range(n-2):
q[i] = 3*(f[i+1] - f[i])
A = np.zeros(((n-2),(n-2)))
#define A for j=0
A[0,0] = 2*(h[0] + h[1])
A[0,1] = h[1]
#define A for j = n-2
A[-1,-2] = h[-2]
A[-1,-1] = 2*(h[-2] + h[-1])
#define A for in the middle
for j in range(1,(n-3)):
A[j,j-1] = h[j]
A[j,j] = 2*(h[j] + h[j+1])
A[j,j+1] = h[j+1]
Ainv = inv(A)
B = Ainv.dot(q)
b = (n)*[0]
b[1:(n-1)] = B
# now we find a, b, c and d
a = [0]*(n-1)
c = [0]*(n-1)
d = [0]*(n-1)
s = [0]*(n-1)
for r in range(n-1):
a[r] = 1/(3*h[r]) * (b[r + 1] - b[r])
c[r] = f[r] - h[r]*((2*b[r] + b[r+1])/3)
d[r] = yd[r]
#solution 1 start
for m in range(n-1):
if xd[m] <= x <= xd[m+1]:
s = a[m]*(x - xd[m])**3 + b[m]*(x-xd[m])**2 + c[m]*(x-xd[m]) + d[m]
return(s)
#solution 1 end
我想在我的 xdata 域中找到最小值,所以 fmin 不起作用,因为你不能在那里定义边界。我尝试了 fmin_slsqp 和最小化。它们与我编写的 phi
函数不兼容,因此我重写了 phi(x, xd,yd)
并添加了一个额外的变量,使 phi 为 phi(x, xd,yd, m)
。 M 指示我们正在计算解决方案的样条函数的哪个子函数(从 x_m 到 x_m+1)。在代码中,我们将 #solution 1
替换为以下
# solution 2 start
return(a[m]*(x - xd[m])**3 + b[m]*(x-xd[m])**2 + c[m]*(x-xd[m]) + d[m])
# solution 2 end
要找到域中的最小值 x_m 到 x_(m+1) 我们使用以下代码:(我们使用一个实例,其中 m=0,所以 x 从 0.25 到 0.5。初始猜测是 0.3)
fmin_slsqp(phi, x0 = 0.3, bounds=([(0.25,0.5)]), args=(xdata, ydata, 0))
然后我会做的(我知道这很粗糙)是用 for 循环迭代它以找到所有子域的最小值,然后取整体最小值。但是,函数fmin_slsqp
不断地returns将初始猜测作为最小值。所以出了点问题,我不知道如何解决。如果你能帮助我,我将不胜感激。感谢您阅读到这里。
当我绘制你的函数 phi
和你输入的数据时,我看到它的范围是 1e-12 的数量级。但是,fmin_slsqp
无法处理该级别的精度,并且无法在您的 objective.
中找到任何更改
我建议的解决方案是按相同的精度顺序缩放 objective 的 return,如下所示:
return(s*1e12)
然后你会得到很好的结果。
>>> sol = fmin_slsqp(phi, x0=0.3, bounds=([(0.25, 0.5)]), args=(xdata, ydata))
>>> print(sol)
Optimization terminated successfully. (Exit mode 0)
Current function value: 1.0
Iterations: 2
Function evaluations: 6
Gradient evaluations: 2
[ 0.25]
我正在寻找自然三次样条的最小值。我编写了以下代码来查找自然三次样条。 (我已经得到了测试数据,并确认这个方法是正确的。)现在我不知道如何找到这个函数的最小值。
这是数据
xdata = np.linspace(0.25, 2, 8)
ydata = 10**(-12) * np.array([1,2,1,2,3,1,1,2])
这是函数
import scipy as sp
import numpy as np
import math
from numpy.linalg import inv
from scipy.optimize import fmin_slsqp
from scipy.optimize import minimize, rosen, rosen_der
def phi(x, xd,yd):
n = len(xd)
h = np.array(xd[1:n] - xd[0:n-1])
f = np.divide(yd[1:n] - yd[0:(n-1)],h)
q = [0]*(n-2)
for i in range(n-2):
q[i] = 3*(f[i+1] - f[i])
A = np.zeros(((n-2),(n-2)))
#define A for j=0
A[0,0] = 2*(h[0] + h[1])
A[0,1] = h[1]
#define A for j = n-2
A[-1,-2] = h[-2]
A[-1,-1] = 2*(h[-2] + h[-1])
#define A for in the middle
for j in range(1,(n-3)):
A[j,j-1] = h[j]
A[j,j] = 2*(h[j] + h[j+1])
A[j,j+1] = h[j+1]
Ainv = inv(A)
B = Ainv.dot(q)
b = (n)*[0]
b[1:(n-1)] = B
# now we find a, b, c and d
a = [0]*(n-1)
c = [0]*(n-1)
d = [0]*(n-1)
s = [0]*(n-1)
for r in range(n-1):
a[r] = 1/(3*h[r]) * (b[r + 1] - b[r])
c[r] = f[r] - h[r]*((2*b[r] + b[r+1])/3)
d[r] = yd[r]
#solution 1 start
for m in range(n-1):
if xd[m] <= x <= xd[m+1]:
s = a[m]*(x - xd[m])**3 + b[m]*(x-xd[m])**2 + c[m]*(x-xd[m]) + d[m]
return(s)
#solution 1 end
我想在我的 xdata 域中找到最小值,所以 fmin 不起作用,因为你不能在那里定义边界。我尝试了 fmin_slsqp 和最小化。它们与我编写的 phi
函数不兼容,因此我重写了 phi(x, xd,yd)
并添加了一个额外的变量,使 phi 为 phi(x, xd,yd, m)
。 M 指示我们正在计算解决方案的样条函数的哪个子函数(从 x_m 到 x_m+1)。在代码中,我们将 #solution 1
替换为以下
# solution 2 start
return(a[m]*(x - xd[m])**3 + b[m]*(x-xd[m])**2 + c[m]*(x-xd[m]) + d[m])
# solution 2 end
要找到域中的最小值 x_m 到 x_(m+1) 我们使用以下代码:(我们使用一个实例,其中 m=0,所以 x 从 0.25 到 0.5。初始猜测是 0.3)
fmin_slsqp(phi, x0 = 0.3, bounds=([(0.25,0.5)]), args=(xdata, ydata, 0))
然后我会做的(我知道这很粗糙)是用 for 循环迭代它以找到所有子域的最小值,然后取整体最小值。但是,函数fmin_slsqp
不断地returns将初始猜测作为最小值。所以出了点问题,我不知道如何解决。如果你能帮助我,我将不胜感激。感谢您阅读到这里。
当我绘制你的函数 phi
和你输入的数据时,我看到它的范围是 1e-12 的数量级。但是,fmin_slsqp
无法处理该级别的精度,并且无法在您的 objective.
我建议的解决方案是按相同的精度顺序缩放 objective 的 return,如下所示:
return(s*1e12)
然后你会得到很好的结果。
>>> sol = fmin_slsqp(phi, x0=0.3, bounds=([(0.25, 0.5)]), args=(xdata, ydata))
>>> print(sol)
Optimization terminated successfully. (Exit mode 0)
Current function value: 1.0
Iterations: 2
Function evaluations: 6
Gradient evaluations: 2
[ 0.25]