メインコンテンツ

lsqlin のコード生成

求解する線形最小二乗問題

範囲指定と線形不等式制約を課された C*x – d のノルムを最小化する問題のための疑似ランダム データを作成します。lb = –1ub = 1 の範囲指定と 150 の線形制約 A*x <= b を課された 15 の変数用の問題を作成します。

N = 15; % Number of variables
rng default % For reproducibility
A = randn([10*N,N]);
b = 5*ones(size(A,1),1);
Aeq = []; % No equality constraints
beq = [];
ub = ones(N,1);
lb = -ub;
C = 10*eye(N) + randn(N);
C = (C + C.')/2; % Symmetrize the matrix
d = 20*randn(N,1);

lsqlin による求解

コードの生成には "active-set" アルゴリズムが必要で、さらにこのアルゴリズムには初期点 x0 が必要です。コード生成に必要なアルゴリズムを使用して MATLAB® で問題を解くには、オプションと初期点を設定します。また、生成されるコードと同じコードを使用して MATLAB で結果を確認できるように、UseCodegenSolver オプションを true に設定します。

x0 = zeros(size(d));
options = optimoptions("lsqlin",...
    Algorithm="active-set",UseCodegenSolver=true);

問題を解くには、lsqlin を呼び出します。

[x,fv,~,ef,output,lam] = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0,options);
Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

lsqlin でこの問題が解けたら、各タイプの非ゼロラグランジュ乗数の数を確認します。非ゼロラグランジュ乗数の総数を減算することで、制約なしの解の成分がいくつあるか確認します。

nl = nnz(lam.lower);
nu = nnz(lam.upper);
ni = nnz(lam.ineqlin);
nunconstrained = N - nl - nu - ni;

fprintf('Number of solution components at lower bounds: %g\n',nl);
Number of solution components at lower bounds: 3
fprintf('Number of solution components at upper bounds: %g\n',nu);
Number of solution components at upper bounds: 2
fprintf('Number of solution components at inequality: %g\n',ni);
Number of solution components at inequality: 5
fprintf('Number of unconstrained solution components: %g\n',nunonstrained);
Number of unconstrained solution components: 5

コード生成手順

コード生成を使用して同じ問題を解くには、次の手順を実行します。

  1. これまでのステップをすべて組み込んだ関数を記述します。出力の数を減らすために、Display オプションを "off" に設定します。

    function [x,fv,lam] = solvelsqlin
    N = 15; % Number of variables
    rng default % For reproducibility
    A = randn([10*N,N]);
    b = 5*ones(size(A,1),1);
    Aeq = []; % No equality constraints
    beq = [];
    ub = ones(N,1);
    lb = -ub;
    C = 10*eye(N) + randn(N);
    C = (C + C.')/2; % Symmetrize the matrix
    d = 20*randn(N,1);
    x0 = zeros(size(d));
    options = optimoptions("lsqlin",Algorithm="active-set",...
        Display="off",UseCodegenSolver=true);
    [x,fv,~,ef,output,lam] = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0,options);
    
    nl = nnz(lam.lower);
    nu = nnz(lam.upper);
    ni = nnz(lam.ineqlin);
    nunconstrained = N - nl - nu - ni;
    fprintf("Number of solution components at lower bounds: %g\n",nl);
    fprintf("Number of solution components at upper bounds: %g\n",nu);
    fprintf("Number of solution components at inequality: %g\n",ni);
    fprintf("Number of unconstrained solution components: %g\n",nunconstrained);
    end
  2. コード生成用の構成を作成します。この場合は、"mex" を使用します。

    cfg = coder.config("mex");
  3. 関数 solvelsqlin 用のコードを生成します。

    codegen -config cfg solvelsqlin
    Code generation successful.
  4. solvelsqlin_mex.mexw64 のような名前の生成されたファイルを実行して、生成されたコードをテストします。

    [x2,fv2,lam2] = solvelsqlin_mex;
    Number of solution components at lower bounds: 3
    Number of solution components at upper bounds: 2
    Number of solution components at inequality: 5
    Number of unconstrained solution components: 5
  5. 各タイプのラグランジュ乗数の数は、2 つの解で同じになります。返された解の差異を確認します。

    disp([norm(x - x2), abs(fv - fv2)])
          1.0e-12 *
    
        0.0006    0.4547

    2 つの解に重要な差はありません。

参考

| | (MATLAB Coder) |

トピック