|
Документ взят из кэша поисковой машины. Адрес
оригинального документа
: http://num-anal.srcc.msu.ru/lib_na/cat/mn_htm_p/mnp1r_p.htm
Дата изменения: Thu Nov 5 15:28:26 2015 Дата индексирования: Sun Apr 10 04:10:44 2016 Кодировка: Windows-1251 |
|
Текст подпрограммы и версий mnp1r_p.zip |
Тексты тестовых примеров tmnp1r_p.zip |
Решение задачи квадратичного программирования при наличии линейных ограничений на переменные методом покоординатного спуска.
Решается задача квадратичного программирования
min { (p, x) + (Cx, x) | Ax ≤ b } ,
где p и x - векторы длины n, C - симметричная положительно определенная матрица размерности n * n, A - матрица размерности m * n, b - вектоp длины m.
Матрица C задается в компактной форме, т.е. представляется в виде вектоpа длины n (n + 1)/2, который состоит из элементов нижнего треугольника матрицы, выписанных последовательно по строкам.
Для получения решения исходной задачи строится двойственная ей задача: найти
min { φ(u) = (h, u) + (Gu, u) | u ≥ 0 } ,
где h = (1/2)AC - 1 p + b, G = (1/4)AC - 1 A*.
Для решения двойственной задачи используется метод покоординатного спуска. Для одномерной минимизации функции φ (u) вдоль направления спуска используется метод квадратичной аппроксимации.
Некоторая вычисленная точка uk ≥ 0 считается точкой минимума функции φ (u), если выполнено хотя бы одно из следующих условий:
| 1. | | uik - uik - 1 | ≤ XE для всех i =1, ..., m, где uk = (u1k, ..., umk) - точка, полученная на k - ой итерации метода, а XE - заданная точность решения задачи по аргументу; |
| 2. | | φ (uk) - φ (uk - 1) | ≤ FE, где uk - точка, вычисленная на k - ой итерации метода, а FE - заданная точность решения задачи по функционалу. |
Решение исходной задачи вычисляется по формуле
x(u) = - (1/2) C -1 (A*u + p) .
Кюнци Г.П., Крелле B., Нелинейное программирование, Изд - во "Советское радио", M., 1965.
procedure MNP1R(var M :Integer; var N :Integer; var A :Array of Real;
var B :Array of Real; var U :Array of Real;
var P :Array of Real; FSTEP :Real; IPAR :Integer;
var MAXK :Integer; XE :Real; var FE :Real;
var KOUNT :Integer; var CO :Array of Real;
var G :Array of Real; var I0 :Array of Integer;
var IERR :Integer);
Параметры
| M - | число стpок матрицы A (тип: целый); |
| N - | число столбцов матрицы A (тип: целый); |
| A - | вещественный двумерный массив размерности M * N; |
| B - | вещественный вектоp длины M; |
| U - | вещественный вектоp длины max {6M + N + 11; N (N + 1)/2 }; на входе первые N (N + 1)/2 компонент вектоpа содержат компактную запись матрицы C; |
| P - | вещественный вектоp длины N; на входе содержит компоненты вектоpа p, а на выходе содержит компоненты решения исходной задачи; |
| FSTEZ - | начальная длина шага (тип: вещественный); |
| IPAR - | параметр, управляющий вариантом покоординатного спуска: если IPAR = 1, то используется вариант с ускоренным дроблением шага (тип: целый); |
| MAXK - | целая переменная, на входе равная максимально допустимому числу вычислений функции, а на выходе - числу фактически выполненных вычислений функций; |
| XE - | заданная точность решения задачи по аргументу (тип: вещественный); |
| FE - | заданная точность решения задачи по функционалу (тип: вещественный); |
| KOUNT - | целая переменная, содержащая число фактически выполненных итераций метода; |
| CO - | вещественный вектоp длины N (N + 1)/2, используемый как рабочий; |
| G - | вещественный вектоp длины M (M + 1)/2, используемый как рабочий; |
| I0 - | целочисленный вектоp длины M, используемый как рабочий; |
| IERR - | целая переменная, указывающая причину окончания процесса: |
| IERR= 0 - | если найден минимум функции с заданной точностью по аргументу или по функционалу; |
| IERR= 4 - | если выполнено заданное максимальное число вычислений функции, но точность не была достигнута; |
| IERR=65 - | если матрица C не является положительно определенной. |
Версии: нет
Вызываемые подпрограммы
| AIH1R - | подпрограмма обращения положительно определенной симметричной матрицы с компактной формой представления методом квадратного корня. |
Замечания по использованию
|
Матрица C должна быть положительно определенной. Используются служебные подпрограммы: MNP11, MNP12, MNP13, MNP14, MNP15, MNP16. |
Найти
min [ Q(x) = 0.5x12 + 0.5x22 - x1 - 2x2 ]
при ограничениях
2x1 + 3x2 ≤ 6 ,
x1 + 4x2 ≤ 5 ,
x1 ≥ 0 ,
x2 ≥ 0 .
Точка минимума x* = (13/17, 18/17) .
Unit TMNP1R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, MNP1R_p;
function TMNP1R: String;
implementation
function TMNP1R: String;
var
KOUNT,IERR,_i :Integer;
IO :Array [0..3] of Integer;
СО :Array [0..2] of Real;
G :Array [0..9] of Real;
const
ХЕ :Real = 1.0E-8;
FE :Real = 1.0E-9;
M :Integer = 4;
N :Integer = 2;
FSTEZ :Real = 1.0;
IPAR :Integer = 1;
МАХК :Integer = 250;
A :Array [0..7] of Real = ( 2.0,1.0,-1.0,0.0,3.0,4.0,0.0,-1.0 );
B :Array [0..3] of Real = ( 6.0,5.0,0.0,0.0 );
P :Array [0..1] of Real = ( -1.0,-2.0 );
U :Array [0..36] of Real = ( 0.5,0.0,0.5,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0 );
begin
Result := ''; { результат функции }
MNP1R(M,N,A,B,U,P,FSTEZ,IPAR,MAXK,
XE,FE,KOUNT,CO,G,IO,IERR);
Result := Result + Format('%s',[' IERR=']);
Result := Result + Format('%2d ',[IERR]);
Result := Result + Format('%s',[' KOUNT=']);
Result := Result + Format('%3d ',[KOUNT]);
Result := Result + Format('%s',[' MAXK=']);
Result := Result + Format('%4d ',[MAXK]);
Result := Result + Format('%s',[' P(1)=']);
Result := Result + #$0D#$0A;
for _i:=0 to 1 do
begin
Result := Result + Format('%20.16f ',[P[_i]]);
if ( ((_i+1) mod 2)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
UtRes('TMNP1R',Result); { вывод результатов в файл TMNP1R.res }
exit;
end;
end.
Результаты:
IERR = 0
KOUNT = 95
MAXK = 117
P(1) = 7.646729 - 01
P(2) = 1.058691 + 00