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;次に、を呼び出す前にメイン関数で宣言するだけfsolveです。出力関数をネストされた関数にすることでこれを達成することもできます。この場合、out配列は外側のメイン関数と共有されます。

out2 番目の出力関数の例では、配列全体を再割り当てする代わりに、動的メモリ割り当てを実行します。これを回避する方法はありません。なぜなら、私たちもアルゴリズムも、収束するまでに何回の反復が必要かを知らないからです。ただし、数百回の反復では、動的メモリ割り当てはかなり高速になります。

ツールが手元にあるので、「収束率の決定」はあなたに任せます...

于 2014-11-07T06:40:54.887 に答える