MATLAB中求解多个非线性独立方程的最快方法？

Question

MATLAB有两种求解非线性方程的方法：

• fzero：求解一个单一的非线性方程
• fsolve：解一个非线性方程组

因此，可以使用以下方法求解n非线性独立方程组：

1. 用一个循环用fzero分别求解方程
2. 用一个循环用fsolve分别求解方程
3. 使用fsolve一起解决它们

我的直觉是：

• 当复杂度(梯度计算)为0(n^2)时，循环方法比单个大n系统更快。
• 对于小型n，循环可能比较慢，因为一个循环在MATLAB中有很高的开销，并且可能有一些恒定的启动时间。
• fzero比fsolve更快，因为它是专为一个非线性方程而做的。

问题：，对于一个n非线性独立方程组，最快的方法是什么？还有其他方法吗？应该使用哪些选项来加速这一进程？

示例问题：(可用于对不同答案进行基准测试)

x_i^2 = 1

i在1和n之间

相关线程

• 解多个独立的非线性方程组
• 求解大量不带环的非线性方程组

Answer 1

我加入这个答案是为了详细说明我上面的评论。根据我的经验，最快的方法是使用文件交换：链接到向量化二分码上可用的fzero矢量化版本。

下面是几个基准，比较它的性能与(i)环fzero和(ii)独立的fsolve。

f = @(x) x.^2-1; %the set of non-linear equations
ns = 1e5; %size of the problem

% method 1: looped fzero
t = timeit(@() fzero(f, rand(1)));
loopFZero = t*ns

% method 2: independent fsolve
options=optimset('Display','off'); % disable displaying
options.Algorithm = 'trust-region-reflective';
options.JacobPattern = speye(ns);
options.PrecondBandWidth = 0;
indepFSolve = timeit(@() fsolve(f, rand(ns,1), options))

% method 3: vectorized bisection, available here:
% https://www.mathworks.com/matlabcentral/fileexchange/28150-bisection-method-root-finding
vecBisection = timeit(@() bisection(f, zeros(ns,1), 2))

结果

%---------------
% ns = 10

loopFZero =

    0.0027


indepFSolve =

    0.0049


vecBisection =

   5.0978e-05

%---------------
% ns = 1e5

loopFZero =

   28.7574


indepFSolve =

    7.7601


vecBisection =

    0.0013

Answer 2

评估某些方法性能的最佳方法是编写一个基准。审议了四起案件：

1. 循环fzero:使用一个循环分别用fzero求解方程
2. 循环解:使用fsolve分别用循环来求解方程
3. 默认解:将方程作为一个方程组一起求解。
4. 独立的解决方案:与默认的解决方案相同，但指定方程是独立的。

f = @(x) x.^2-1; % the set of non-linear equations
ns = 1:1:100; % the sizes for which the benchmark is performed
options=optimset('Display','off'); % disable displaying

figure
hold on
plot(ns, loopFSolve(f, ns, options), 'DisplayName', 'loop fsolve')
plot(ns, loopFZero(f, ns, options), 'DisplayName', 'loop fzero')
plot(ns, defaultFSsolve(f, ns, options), 'DisplayName', 'default fsolve')
plot(ns, independentFSolve(f, ns, options), 'DisplayName', 'independent fsolve')

legend ('Location', 'northwest')

function t = loopFZero(f, ns, options)
  t1 = timeit(@() fzero(f, rand(1), options));
  t = ns * t1;
end

function t = loopFSolve(f, ns, options)
  t1 = timeit(@() fsolve(f, rand(1), options));
  t = ns * t1;
end

function t = defaultFSsolve(f, ns, options)
  t = zeros(size(ns));
  for i=1:length(ns)
    n = ns(i);
    un = rand(n, 1);
    t(i) = timeit(@() fsolve(f, un, options));
  end
end

function t = independentFSolve(f, ns, options)
  t = zeros(size(ns));
  for i=1:length(ns)
    n = ns(i);
    un = rand(n, 1);
    options.Algorithm = 'trust-region-reflective';
    options.JacobPattern = speye(n);
    options.PrecondBandWidth = 0;

    t(i) = timeit(@() fsolve(f, un, options));
  end
end

结果

所有的图形都显示了整个系统在n函数中的计算时间，方程的个数。

前两个数字绘制的n高达1000，间隔为100。最后两幅图绘制的n高达100，间隔为1。对于每一幅图，第二幅图与第一幅图相同，但没有循环fzero，因为它比其他图慢得多。

​

结论

1. 循环f:不要使用它，启动时间太长。
2. 循环fzero:您可以将其用于小型n (n < ~20的最快方法)。
3. 默认n：您可以将其用于相对较小的~20 < n < ~50 (~20 < n < ~50最快的方法，但与2和3的差异相对较小)。
4. 独立解决方案:您应该将它用于大型n (~50 < n的最快方法)

在一般情况下，您应该使用独立的，只有对于小的n循环fzero才可以使用，即使用fsolve具有以下选项：

options.Algorithm = 'trust-region-reflective';
options.JacobPattern = speye(n);
options.PrecondBandWidth = 0;

Lazy people可能只是使用默认的解决方案，因为它对于中等数量的方程(n < ~200)具有合理的性能。

重要备注

鲍勃指出，fzero可能要快几个数量级。根据Bob的测试结果，他的方法对于问题的每一个大小，即大小n.，都是最快的方法。

备注

• 注意，默认情况下的时间复杂度是O(n^2)，而其他的是O(n)。
• 请注意，在某些边界情况下，fzero和fsolve的行为可能有所不同，f.e fzero在搜索符号更改的位置时不会为x^2找到解决方案。