查找用户定义函数的局部最大值和最小值

Question

我想要的

我想找到一个包含固定点、它们的值和位置以及它们是最小值还是最大值的列表。

我的函数如下所示：

import numpy as np

def func(x,y):
  return (np.cos(x*10))**2 + (np.sin(y*10))**2

方法

以下是我正在考虑使用的方法：

实际上我已经在 Mathematica 上做了类似的事情。我一次又一次地区分函数。我查看一阶导数为 0 的点，计算它们的值和位置。然后我在这些位置取二阶导数并检查它们是最小值还是最大值。
我也想知道是否只制作 x 和 y 中函数值的二维数组，并找到该数组的最大值和最小值。但这需要我知道如何精细地定义 x 和 y 网格以可靠地捕获函数的行为

对于后一种情况，我已经找到了一些方法，例如 this one。

我只是想知道，在 Python 中，哪种方法在效率、速度、准确性甚至优雅方面更有意义？

Answer 1

find a list of the stationary points, of their values and locations, and of whether they are minima or maxima.

这通常是一个无法解决的问题。方法 1（符号）适用于此，但对于复杂的函数，没有固定点的符号解（没有符号解两个方程的一般系统的方法）。

使用 SymPy 的符号解决方案

对于像您的示例这样的简单函数，SymPy 可以正常工作。这是一个找到固定点并根据 Hessian 的特征值对它们进行分类的完整示例。

import sympy as sym
x, y = sym.symbols("x y")
f = sym.cos(x*10)**2 + sym.sin(y*10)**2
gradient = sym.derive_by_array(f, (x, y))
hessian = sym.Matrix(2, 2, sym.derive_by_array(gradient, (x, y)))

至此，Hessian 矩阵为 2×2 的符号矩阵：[[200*sin(10*x)**2 - 200*cos(10*x)**2, 0], [0, -200*sin(10*y)**2 + 200*cos(10*y)**2]]。接下来，我们通过使gradient等于零来找到驻点，并将它们一一插入到Hessian中。

stationary_points = sym.solve(gradient, (x, y))
for p in stationary_points:
    value = f.subs({x: p[0], y: p[1]})
    hess = hessian.subs({x: p[0], y: p[1]})
    eigenvals = hess.eigenvals()
    if all(ev > 0 for ev in eigenvals):
        print("Local minimum at {} with value {}".format(p, value))
    elif all(ev < 0 for ev in eigenvals):
        print("Local maximum at {} with value {}".format(p, value))
    elif any(ev > 0 for ev in eigenvals) and any(ev < 0 for ev in eigenvals):
        print("Saddle point at {} with value {}".format(p, value))
    else:
        print("Could not classify the stationary point at {} with value {}".format(p, value))

最后一个子句是必要的，因为当Hessian只是semidefinite时，我们无法判断那个(x**2 + y**4和x**2 - y**4 在 (0, 0) 处具有相同的 Hessian 但行为不同）。输出：

Saddle point at (0, 0) with value 1
Local maximum at (0, pi/20) with value 2
Saddle point at (0, pi/10) with value 1
Local maximum at (0, 3*pi/20) with value 2
Local minimum at (pi/20, 0) with value 0
Saddle point at (pi/20, pi/20) with value 1
Local minimum at (pi/20, pi/10) with value 0
Saddle point at (pi/20, 3*pi/20) with value 1
Saddle point at (pi/10, 0) with value 1
Local maximum at (pi/10, pi/20) with value 2
Saddle point at (pi/10, pi/10) with value 1
Local maximum at (pi/10, 3*pi/20) with value 2
Local minimum at (3*pi/20, 0) with value 0
Saddle point at (3*pi/20, pi/20) with value 1
Local minimum at (3*pi/20, pi/10) with value 0
Saddle point at (3*pi/20, 3*pi/20) with value 1

显然，solve 没有找到所有解（其中有无穷多个）。考虑 solve vs solveset 但无论如何，处理无限多的解决方案是困难的。

使用 SciPy

进行数值优化

SciPy 提供了很多 numerical minimization routines, including brute force（这是您的方法 2；通常它非常非常慢）。这些都是强大的方法，但请考虑以下几点。

每个运行只会找到一个最小值。
将 f 替换为 -f 也可以找到最大值。
更改搜索的起点（minimize 的参数 x0）可能会产生另一个最大值或最小值。不过，您永远不会知道是否还有其他您还没有看到的极值。
None 这将找到鞍点。

混合策略

使用 lambdify 可以将符号表达式转换为 Python 可以传递给 SciPy 数值求解器的函数。

from scipy.optimize import fsolve
grad = sym.lambdify((x, y), gradient)
fsolve(lambda v: grad(v[0], v[1]), (1, 2))

这个returns一些固定点，在这个例子中是[0.9424778 , 2.04203522]。它是哪一点取决于最初的猜测，即 (1, 2)。通常（但不总是）您会得到接近初始猜测的解决方案。

这比直接最小化方法的优势在于也可以检测到鞍点。尽管如此，还是很难找到所有的解决方案，因为 fsolve 中的每个运行只能提出一个。

查找用户定义函数的局部最大值和最小值

Finding local maxima and minima of user defined functions

python

numpy

minimization

我想要的

方法

使用 SymPy 的符号解决方案

使用 SciPy

混合策略