Этот пример показывает основанные на проблеме подходы для решения проблемы
C*x = d
в смысле наименьших квадратов. Здесь C и d дают массивы, и x является неизвестным массивом, и x ограничивается, чтобы быть неотрицательным. Другими словами,
минимизируйте sum((C*x - d).^2)
подвергните x >= 0.
Чтобы начаться, загрузите массивы C и d в вашу рабочую область.
load particleПросмотрите размеры массивов.
sizec = size(C)
sizec = 1×2
2000 400
sized = size(d)
sized = 1×2
2000 1
Создайте переменную x оптимизации соответствующего размера для умножения C. Наложите нижнюю границу 0 на элементах x.
x = optimvar('x',sizec(2),'LowerBound',0);
Создайте выражение целевой функции.
residual = C*x - d; obj = sum(residual.^2);
Создайте задачу оптимизации под названием nonneglsq и включайте целевую функцию в проблему.
nonneglsq = optimproblem('Objective',obj);Найдите решатель по умолчанию для проблемы.
opts = optimoptions(nonneglsq)
opts =
lsqlin options:
Options used by current Algorithm ('interior-point'):
(Other available algorithms: 'trust-region-reflective')
Set properties:
No options set.
Default properties:
Algorithm: 'interior-point'
ConstraintTolerance: 1.0000e-08
Display: 'final'
LinearSolver: 'auto'
MaxIterations: 200
OptimalityTolerance: 1.0000e-08
StepTolerance: 1.0000e-12
Show options not used by current Algorithm ('interior-point')
Решите проблему с помощью решателя по умолчанию.
[sol,fval,exitflag] = solve(nonneglsq)
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.
sol = struct with fields:
x: [400x1 double]
fval = 22.5795
exitflag =
OptimalSolution
Решатель по умолчанию, lsqlin, работал быстро. Но существуют другие решатели, которые решают проблемы этой природы. lsqnonneg решает проблемы этого типа, как делает quadprog. Попробуйте каждого.
[sol2,fval2,exitflag2] = solve(nonneglsq,'Solver','lsqnonneg'); [sol3,fval3,exitflag3] = solve(nonneglsq,'Solver','quadprog');
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.
Исследуйте первые десять записей каждого решения.
solutions = table(sol.x(1:10),sol2.x(1:10),sol3.x(1:10),'VariableNames',{'lsqlin','lsqnonneg','quadprog'})
solutions=10×3 table
lsqlin lsqnonneg quadprog
__________ _________ __________
0.049416 0.049416 0.049416
0.082985 0.082985 0.082985
0.078469 0.078469 0.078469
0.11682 0.11682 0.11682
0.07165 0.07165 0.07165
0.060525 0.060525 0.060525
1.2196e-09 0 1.2196e-09
2.703e-10 0 2.703e-10
1.5017e-10 0 1.5017e-10
1.2331e-10 0 1.2331e-10
Несколько записей lsqnonneg, кажется, точно нуль, в то время как другие решения не являются точно нулевыми. Например, восьмая запись в каждом решении:
fprintf('The entries are %d, %d, and %d.\n',sol.x(8),sol2.x(8),sol3.x(8))The entries are 2.702960e-10, 0, and 2.702960e-10.
Возможно, другой quadprog и алгоритм lsqlin были бы немного более точными.
options = optimoptions('lsqlin','Algorithm','trust-region-reflective'); [sol,fval,exitflag] = solve(nonneglsq,'Solver','lsqlin','Options',options);
Local minimum possible. lsqlin stopped because the relative change in function value is less than the function tolerance.
disp(sol.x(8))
2.2800e-15
options = optimoptions('quadprog','Algorithm','trust-region-reflective'); [sol3,fval3,exitflag3] = solve(nonneglsq,'Solver','quadprog','Options',options);
Local minimum possible. quadprog stopped because the relative change in function value is less than the sqrt of the function tolerance, the rate of change in the function value is slow, and no negative curvature was detected.
disp(sol3.x(8))
6.7615e-13
Решения немного ближе, чтобы обнулить использование алгоритма 'trust-region-reflective', чем алгоритмы внутренней точки по умолчанию.