0

假设我正在求解一个非线性方程组。一个简单的例子是:

function example
    x0 = [15; -2];
    options = optimoptions('fsolve','Display','iter','TolFun',eps,'TolX',eps);
    [x,fval,exitflag,output] = fsolve(@P1a,x0,options);
end

function f1 = P1a(x)
    f1 = [x(1)+x(2)*(x(2)*(5-x(2))-2)- 13; x(1)+x(2)*(x(2)*(1+x(2))-14)-29];
end

如何确定收敛速度?'Display''iter'向我展示了每一步的规范,但我找不到提取这些值的方法。(对于这个特定的例子,我相信fsolve不会收敛到正确的解决方案,而是收敛到局部最小值。然而,这不是问题。我只是想找到一种方法来估计收敛速度。)

4

1 回答 1

0

你可以从中得到很多fsolve。但是,您需要做一些工作。阅读 Matlab 优化方法的'OutputFcn'选项和编写输出函数。这与 Matlab 的 ODE 求解器使用的同名选项有所不同。'Display','iter'这是一个复制选项值的示例fsolve(专门针对默认'trust-region-dogleg'算法):

function stop = outfun(x,optimValues,state)
% See private/trustnleqn
stop = false;

switch state
    case 'init'
        header = sprintf(['\n                                         Norm of      First-order   Trust-region\n',...
                          ' Iteration  Func-count     f(x)          step         optimality    radius']);
        disp(header);
    case 'iter'
        iter = optimValues.iteration;               % Iteration
        numFevals = optimValues.funccount;          % Func-count
        F = optimValues.fval;                       % f(x)
        normd = optimValues.stepsize;               % Norm of step
        normgradinf = optimValues.firstorderopt;    % First-order optimality
        Delta = optimValues.trustregionradius;      % Trust-region radius

        if iter > 0
            formatstr = ' %5.0f      %5.0f   %13.6g  %13.6g   %12.3g    %12.3g';
            iterOutput = sprintf(formatstr,iter,numFevals,F'*F,normd,normgradinf,Delta);
        else
            formatstr0 = ' %5.0f      %5.0f   %13.6g                  %12.3g    %12.3g';
            iterOutput = sprintf(formatstr0,iter,numFevals,F'*F,normgradinf,Delta);
        end
        disp(iterOutput);
    case 'done'

    otherwise

end

然后,您可以通过以下方式调用它:

function example
P1a=@(x)[x(1)+x(2)*(x(2)*(5-x(2))-2)- 13; x(1)+x(2)*(x(2)*(1+x(2))-14)-29];
x0 = [15; -2];
opts = optimoptions('fsolve','Display','off','OutputFcn',@outfun,'TolFun',eps,'TolX',eps);
[x,fval,exitflag,output] = fsolve(P1a,x0,opts);

这仍然只是打印到命令窗口。从这里开始创建一个可以将数据写入数组、文件或其他数据结构的输出函数。以下是使用全局变量执行此操作的方法(通常,这不是一个好主意):

function stop = outfun2(x,optimValues,state)
stop = false;
global out;   % Global variable, define in main function too

switch state
    case 'init'
        out = [];
    case 'iter'
        iter = optimValues.iteration;               % Iteration
        numFevals = optimValues.funccount;          % Func-count
        F = optimValues.fval;                       % f(x)
        normd = optimValues.stepsize;               % Norm of step
        normgradinf = optimValues.firstorderopt;    % First-order optimality
        Delta = optimValues.trustregionradius;      % Trust-region radius

        out = [out;iter numFevals F'*F normd normgradinf Delta];
    case 'done'

    otherwise

end

然后global out;在调用之前在你的 main 函数中声明fsolve。您还可以通过将输出函数设为嵌套函数来完成此操作,在这种情况下,out数组将与外部 main 函数共享。

第二个输出函数示例执行动态内存分配,而不是重新分配整个out数组。没有办法解决这个问题,因为我们和算法都不知道需要多少次迭代才能收敛。但是,对于几百次迭代,动态内存分配将非常快。

既然您已经掌握了工具,我将把“确定收敛速度”留给您...

于 2014-11-07T06:40:54.887 回答