我正在寻找一种有效的方法来寻找多个（〉1百万）但独立的四阶多项式在给定/指定范围/边界的局部min。
我有两个要求：

R1：即使对于一百万个不同的多项式方程也是有效的
R2：本地最小值精确到0.01（即2dp）

这是我用scipy创建的一些代码。没问题，但我想知道**在我开始并行编程之前，**是否有其他更好的包来执行这样的任务。
为了说明我的问题，让我们先从一个多项式开始：
下面，我试图在区间（-5，5）内找到4x^4 + 6x^3 + 3x^2 + x + 5的局部极小值。
在我的笔记本电脑上，找到本地最小值（大约为-0.72770502）需要大约2毫秒。
时间是一个多项式的权利，但我会希望更快的东西，因为我需要执行这个操作超过100万次定期。

from scipy import optimize
import numpy as np

# Define a objective and gradient function for 4th order polynomial

# x is the value to be evaluated

# par is a numpy array of len 5 that specifies the polynomial coefficients.

def obj_grad_fun_custom(x,par):
    obj = (np.array([x**4,x**3,x**2,x**1,1]) * par).sum()
    grad = (np.array([4*x**3,3*x**2,2*x,1]) * par[:-1]).sum()
    return obj, grad

# Try minimise an example polynomial of 4x^4 + 6x^3 + 3x^2 + x + 5

# with contrainted bound

res = optimize.minimize(
    fun = obj_grad_fun_custom,
    x0 = 0,
    args=(np.array([4,6,3,1,5])), # polynomial coefficients
    jac=True ,
    bounds=[(-2, 10)],
    tol=1e-10)
print(res.x)

# Timing (this takes about 2 ms for me)

%timeit optimize.minimize(fun = obj_grad_fun_custom, x0 = 0, args=(np.array([4,6,3,1,5])), jac=True, bounds=[(-5, 5)], tol=1e-10)

下面是我计划对100万个不同的4阶多项式进行正则化的方法，我希望局部最小化。希望有人能给我一个比scipy更合适的软件包。或者有其他方法吗？谢谢！


# Multiple polynomials

result = [] # saving the local minima
poly_sim_no = 1000000 #ideally 1 million or even more
np.random.seed(0)
par_set = np.random.choice(np.arange(10), size=(poly_sim_no, 5), replace=True) #generate some order 4 polynomial coefficients 

for a in par_set:
    res = optimize.minimize(obj_grad_fun_custom, 0,args=(a), jac=True ,bounds=[(-5, 5)], tol=1e-10)
    result.append(res.x)

print(result)

import numpy as np from numpy.polynomial import Polynomial def find_extrema(poly, bounds): deriv = poly.deriv() extrema = deriv.roots() # Filter out complex roots extrema = extrema[np.isreal(extrema)] # Get real part of root extrema = np.real(extrema) # Apply bounds check lb, ub = bounds extrema = extrema[(lb <= extrema) & (extrema <= ub)] return extrema def find_minimum(poly, bounds): extrema = find_extrema(poly, bounds) # Note: initially I tried taking the 2nd derivative to filter out local maxima. # This ended up being slower than just evaluating the function. # Either bound could end up being the minimum. Check those too. extrema = np.concatenate((extrema, bounds)) # Check every candidate by evaluating the polynomial at each possible minimum, # and picking the minimum. value_at_extrema = poly(extrema) minimum_index = np.argmin(value_at_extrema) return extrema[minimum_index] # Warning: polynomial expects coeffients in the opposite order that you use. poly = Polynomial([5,1,3,6,4]) print(find_minimum(poly, (-5, 5)))

p = np.poly1d([4,6,3,1,5]) # Note: polynomials are opposite order of before def find_minimum2(poly): roots = np.real(np.roots(poly.deriv())) return roots[np.argmin(poly(roots))] print(find_minimum2(p))

1条答案

按热度按时间

lymgl2op1#

因为你要求多项式的最小值，所以你可以利用多项式求导很容易，而且有很多很好的算法可以求多项式的根。
以下是它的工作原理：
1.首先，求导数。所有最小值的点的导数都是零。
1.寻找那些零，也就是找到导数的根。
1.一旦我们有了候选人的名单，检查解决方案是真实的。
1.检查解决方案是否在您设置的范围内。（我不知道您添加范围是因为您实际上需要范围，还是为了使它运行得更快。如果是后者，请随意删除此步骤。）
1.用多项式实际评估候选项，并找到最小的一个。
下面是代码：

这在我的计算机上花费了162微秒，比scipy.optimize解决方案快了大约6倍。（问题中显示的解决方案在我的计算机上花费了1.12毫秒。）

Edit：更快的替代方法

这里有一个更快的方法。但是，它放弃了边界检查，使用了一个过时的API，并且通常更难阅读。

它的时钟频率为110微秒，比原来的快了大约10倍。

赞(0）回复(0）举报 2022-11-10

scipy Python高效求多个多项式的局部极大值/极小值

1条答案

Edit：更快的替代方法

相关问题

热门标签

最新问答