Этот пример показывает значение использования разреженной арифметики, когда у вас есть разреженная проблема. Матрица имеет n
строки, где вы выбираете n
быть большим значением и несколькими ненулевыми диагональными полосами. Полная матрица размера n
- n
может израсходовать всю доступную память, но разреженная матрица не представляет проблемы.
Проблема состоит в том, чтобы минимизировать x'*H*x/2 + f'*x
при ограничениях
x(1) + x(2) + ... + x(n) <= 0
,
где f = [-1;-2;-3;...;-n]
H
разреженная симметричная ленточная матрица.
Создайте симметричную циркулянтную матрицу H
на основе сдвигов векторного [3,6,2,14,2,6,3]
, с 14 находиться на основной диагонали. Имейте матрицу быть n
- n
, где n = 30,000
.
n = 3e4;
H2 = speye(n);
H = 3*circshift(H2,-3,2) + 6*circshift(H2,-2,2) + 2*circshift(H2,-1,2)...
+ 14*H2 + 2*circshift(H2,1,2) + 6*circshift(H2,2,2) + 3*circshift(H2,3,2);
Просмотрите структуру разреженной матрицы.
spy(H)
Создайте переменную x
оптимизации и проблема
qprob
.
x = optimvar('x',n);
qprob = optimproblem;
Создайте целевую функцию и ограничения. Поместите цель и ограничения в qprob
.
f = 1:n; obj = 1/2*x'*H*x - f*x; qprob.Objective = obj; cons = sum(x) <= 0; qprob.Constraints = cons;
Решите задачу квадратичного программирования с помощью 'interior-point-convex'
по умолчанию алгоритм и линейная алгебра для разреженных матриц. Чтобы помешать решателю останавливаться преждевременно, установите
StepTolerance
опция к 0
.
options = optimoptions('quadprog','Algorithm','interior-point-convex',... 'LinearSolver','sparse','StepTolerance',0); [sol,fval,exitflag,output,lambda] = solve(qprob,'Options',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. <stopping criteria details>
Просмотрите значение целевой функции, количество итераций и множитель Лагранжа, сопоставленный с линейным ограничением неравенства.
fprintf('The objective function value is %d.\nThe number of iterations is %d.\nThe Lagrange multiplier is %d.\n',... fval,output.iterations,lambda.Constraints)
The objective function value is -3.133073e+10. The number of iterations is 7. The Lagrange multiplier is 1.500050e+04.
Оцените ограничение, чтобы видеть, что решение находится на контуре.
fprintf('The linear inequality constraint sum(x) has value %d.\n',sum(sol.x))
The linear inequality constraint sum(x) has value 7.599738e-09.
Сумма компонентов решения является нулем к в допусках.
Решение x
имеет три области: начальный фрагмент, итоговый фрагмент и приблизительно линейный фрагмент по большей части решения. Постройте эти три области.
subplot(3,1,1) plot(sol.x(1:60)) title('x(1) Through x(60)') subplot(3,1,2) plot(sol.x(61:n-60)) title('x(61) Through x(n-60)') subplot(3,1,3) plot(sol.x(n-59:n)) title('x(n-59) Through x(n)')