Оперативное управление на основе решения задач линейного программирования



Введение

Характерной особенностью современных производственных процессов как объектов управления является наличие в них большого числа управляемых и управляющих величин, которые влияют друг на друга в том смысле, что изменение какой-либо одной из них вызывает изменение всех других. Это обстоятельство требует специального подхода к решению задач управления такими объектами.

Производственному процессу обязательно предшествует этап планирования. Производственные планы определяют, что должно производиться, в каком количестве, кем, когда и каким образом, т.е. устанавливают объект, номенклатуру, объем, сроки и технологию производства. Однако в ходе реализации производственных планов возникают случайные отклонения, обусловленные особенностями процесса производства.

Целью создания АСУ является оптимизация производственных процессов, т.е. получение максимального эффекта от затрачиваемых материальных и трудовых ресурсов. В функции подсистемы оперативного управления входит коррекция планов на основе оперативной информации о состоянии производственного объекта с целью уменьшения влияния возмущающих факторов на ход производственного процесса.

Задачи оптимизации управления, рассматриваемые в пособии, сводятся к задачам линейного и нелинейного программирования. Основное внимание уделяется методам решения, которые позволяют разработать несколько более сложный в логическом плане, но менее трудоемкий, т.е. оперативный метод решения. Пособие можно рассматривать как набор методов, органически друг с другом связанных и друг на друга опирающихся и позволяющих существенно повысить оперативность принятия управленческих решений.

Рассматриваемый подход оперативного решения задач математического программирования продемонстрирован рядом примеров оперативного решения задач линейного и нелинейного программирования. Уделено внимание в пособии вопросам оперативной информационной поддержки для принятия управленческих решений. Математические модели и методы, рассмотренные в пособии, могут найти применение и для автоматизации управления в условиях чрезвычайных ситуаций [1].

  1. Оперативное управление на основе решения задач

линейного программирования

Большой класс задач оптимизации [2,3] можно свести к задаче линейного программирования — В пособии используются многомерно-матричные обозначения и операции [4]: А[1,1] – матрица; А(1,0) – вектор-столбец; А(0,1) – вектор-строка; А(0,0) или просто А (с индексом или без индекса) – скаляр.

I =C(0,1)X(1,0)max                                                   (1)

при наличии ограничений

                                    A(1,1)X(1,0)B(1,0),                                                        (2)

X(1,0) 0                                .

ЗдесьI – критерий качества оптимального плана, С(0,1) – вектор-строка коэффициентов целевой функцииI, Х – вектор-столбец искомых переменных, А – матрица ограничений, определяемых конкретным производственным процессом, В – вектор-столбец правых частей ограничений. Можно считать, что вектор В(1,0) задает величину ресурсов, которые можно использовать в производственном процессе.

Составляющие векторов С(0,1) и В(1,0) в сильной степени зависят от текущих затрат на производство, транспортных расходов, от цен на комплектующие изделия и т.д.; поэтому эти составляющие являются менее стабильными, чем коэффициентыAij матрицы А, характеризующие расходi-ресурса на единицуj-продукции. В связи с этим в первую очередь будет уделено внимание оперативному управлению в условиях отклонения векторов С(0,1), В(1,0) от их «плановых» значений. Случай увеличения расходаi-ресурса наXj продукцию (например, смена производственного оборудования) также требует при ограниченных ресурсах Вj оперативного вмешательства в производственный процесс. Решение этих вопросов дано в этой главе.

1.1. Общая схема решения задачи линейного программирования

на  основе  определения  фундаментальной  системы  решений

для системы однородных линейных неравенств

На основе решений системы линейных неравенств, изложенных в пособии [5], можно построить алгоритм оперативного решения задачи линейного программирования. Например, для системы ограничений (2) можно получить решение в виде

X1=11P1+12P2+1mPm

X2=21P1+22P2+2mPm                                                (3)

  Xn+1=n+1,1P1+n+1,2P2+n+1,mPm          .

Здесьn – число элементов вектораX(1,0). Подставляя Х1, Х2,…,Хn в выражение (1) получим целевую функциюI, выраженную через переменныеPi, i=. Задача линейного программирования имеет вид

                                           I =                                                         (4)

при наличии ограниченияxn+1 = 1                                                                (5)

Время решения задачи линейного программирования грубо оценивается соотношением [6]:

Время =k(число ограничений)3.

Здесьk – константа пропорциональности.

Поэтому в задачах линейного программирования рекомендуют переходить к виду задачи с малым числом ограничений [7].

Важно отметить, что приведенная оценка является достаточно грубой и ее следует рассматривать лишь как рекомендацию по выбору метода решения задачи. Излагаемый в пособии подход позволит разработчику расширить арсенал применяемых методов, особенно для решения нелинейных задач оптимизации. Однако приведенные в пособии специально подобранные примеры, указывают на необходимость тщательного анализа оперативности решения в каждом конкретном случае.

1.2. Оперативное управление при изменении коэффициентов

целевой функции

Исходные данные для задач с целью методического согласования взяты из [7], т.е. в начале проводится этап планирования по методике [7], затем при изменении части переменных задачи осуществляется оперативное управление.

Пример 1. На этапе планирования решить задачу линейного программирования

I = 3x1 + 2x2max                                                                (6)

при наличии ограничений:

                    -5x1 + 4x2 20

   2x1+ 3x2 24                           (7)

x1   -  3x2 3

x1   + 0x2 5

Решение задачи осуществляется стандартными методами линейного программирования [7] и имеет вид х1 = 5,  х2= 4,7.

При этом величина критерияI = 24,4.

Пример2. На этапе управления оперативно скорректировать план, полученный в примере 1, при новых значениях коэффициентов целевой функции с1 = 2; с2 = 4.

Решение системы линейных неравенств (7) имеет вид

X1 = 0P1 + 0P2 + 1,565P3 + 12P4 + 37,644P5 + 52,174P6

X2 = 0P1 + 5P2 + 6,957P3 +  0P4  + 35,135P5 +   6,957P6                             (8)

X3 = 1P1 + 1P2 + 1P3        +  4P4  +   7,529P5 + 10,435P6 = 1         .

Подставляя х1, х2 в выражение для целевой функции при измененных значениях с1 и с2,приходим к задаче линейного программирования

F = 0P1 + 20P2 + 30,958P3 + 24P4 + 215,825P5 + 132,176P6 max            (9)

при наличии ограничения

1P1 +1P2 +1P3 +4P4 +7,529P5 +10,435P6 =1 .(10)

Решение задачи (9), (10) имеет вид:

P1 = 0;P2 = 0;P3 = 1;P4 = 0;P5 = 0;P6 = 0 .

При этом х1 = 1,565, х2 = 6,957,F = 31,0, что означает необходимость введения оперативных изменений в начальный план производства.

С целью демонстрации возможности рассматриваемого подхода приведем решение с помощью симплексного метода без применения и с применением решения систем линейных неравенств.

Пример 1а. Найти максимум критерия (6) при наличии ограничений (7). Решение без применения решения системы линейных неравенств представлено в табл.1-4.

                                        Таблица 1

1

2

В

Y1

-5

4

20

Y2

2

3

24

Y3

1

-3

3

Y4

1

0

5

fmax

-3

-2

0

А(3,1) – разрешающий элемент.

                                        Таблица 2

-Y3

2

В

Y1

5

-11

35

Y2

-2

9

18

X1

1

-3

3

Y4

-1

3

2

fmax

3

-11

9

А(4,2) – разрешающий элемент.

Таблица 3

-Y3

-Y4

В

Y1

1,3

3,7

42,3

Y2

1,0

-3,0

12,0

X1

0

1

5,0

X2

-0,33

0,33

0,67

fmax

-0,67

3,7

16,3

Таблица4

-Y2

-Y4

В

Y1

-1,3

7,7

26,3

Y3

1,0

-3,0

12,0

X1

0

1

5,0

X2

0,33

-0,67

4,7

fmax

0,67

1,7

24,3

А(2,1) – разрешающий элемент.

Решение   найдено.

 Общее число вычисляемых элементов  матриц  равно  60.

Пример 1б. Найти максимум критерия (6) при наличии ограничений (7). Решение задачи с применением решения систем линейных неравенств представлено табл.5-7.

                                                                                                             Таблица 5

-X1

-X2

-X3

-X4

-X5

-X6

B

0

1

1

1

4

7,5

10,4

1,0

fmax

0

-10

-18,6

-36

-183,2

-170,4

0,0

A(1,1) – разрешающий элемент.

     Таблица 6

0

-X2

-X3

-X4

-X5

-X6

B

X1

1

1

1

4

7,5

10,4

1,0

fmax

0

-10

-18,6

-36

-183,2

-170,4

0

A(1,5) – разрешающий элемент.

                                                                                                             Таблица 7

0

-X2

-X3

-X4

-X1

-X6

B

X5

0,13

0,13

0,13

0,53

0,13

1,4

0,13

fmax

24,3

14,3

5,7

61,3

24,3

83,5

24,3

Решение найдено. Общее число вычисляемых элементов матриц уменьшилось и равно 42 (вместо 60 в предыдущем случае).

1.3. Оперативное управление при изменении величины ресурсов

Исходные условия задачи на этапе планирования заданы соотношениями (6), (7).

Пример 3. На этапе управления оперативно скорректировать план, полученный в примере 1, при новом значении ресурса В2 = 15 вместо В2 = 24.

Решение системы линейных неравенств исходной задачи имеет вид (8). Подставляя значения х1, х2 во второе неравенство системы (7), получаем новое ограничение

15Р2 + 24,171Р3 + 24Р4 + 180,693Р5 + 125,219Р6 15 .                              (11)

Таким образом, приходим к задаче линейного программирования с малым числом ограничений: найти максимум критерия (6) после подстановки в него х1, х2 из системы (8).

F = 0 P1 + 10P2 + 18,609P3 + 36P4 + 183,202P5 + 170,436P6 max         (12)

при наличии ограничений (10), (11).

Решение задачи имеет вид:

Р1 = 0; Р2 = 0; Р3 = 0; Р4 = 0; Р5 = 0,033207; Р6 = 0,072.

При этом х1 = 4,8; х2 = 1,54;

Fmax = 18,33,

что, как и в предыдущем примере, означает необходимость введения оперативных изменений в начальный план производства.

Пример 3а. Найти максимум критерия (6) при наличии ограничений (7) и дополнительного ограничения

  2x1+ 3x2 15.

Решение без применения решения системы линейных неравенств представлено табл.8-11.

                                       Таблица 8

1

2

В

Y1

-5

4

20

Y2

2

3

24

Y3

1

-3

3

Y4

1

0

5

Y5

2

3

15

fmax

-3

-2

0

А(3,1) – разрешающий элемент.

                                        Таблица 9

-Y3

2

В

Y1

5

-11

35

Y2

-2

9

18

X1

1

-3

3

Y4

-1

3

2,0

Y5

-2

9

9,0

fmax

3

-11

9

А(4,2) – разрешающий элемент.

                                                Таблица 10

-Y3

-Y4

В

Y1

1,3

3,7

42,3

Y2

1

-3

12

X1

0

1

5

X2

-0,33

0,33

0,67

Y5

1,0

-3,0

3,0

fmax

-0,67

3,7

16,3

А(5,1) – разрешающий элемент.

Таблица 11

-Y5

-Y4

В

Y1

-1,3

7,7

38,3

Y2

-1

0

9

X1

0

1

5,0

X2

0,33

-0,67

1,7

Y3

1

-3

3,0

fmax

0,67

1,7

18,3

Решение найдено.

Общее число вычисляемых элементов матриц равно 72.

Пример 3б. Найти максимум критерия (6) при наличии ограничений (7) и дополнительного ограничения

  2x1+ 3x2 15.

Решение системы линейных ограничений при этом будет иметь вид:

X1 = 12P3 + 52,174P4 + 69,565P5 ;

X2 =   5P2 +   6,957P4 + 23,188P5 ;

X2 = P1 + P2 + 4P3 + 10,435P4 + 13,913P5 = 1 .

Критерий (6) примет вид

F = 10P2 +36P3 + 170,436P4 + 264,071P5.

Решение задачи 3б с применением решения системы линейных неравенств представлено табл.12-13.

                                                                                           Таблица 12

-X1

-X2

-X3

-X4

-X5

B

0

1

1

1

10,4

13,9

1,0

fmax

0

-10

-36

-170,4

-264,1

0

A(1,5) – разрешающий элемент.

         Таблица 13

0

-X2

-X3

-X4

-X1

B

X5

0,07

0,07

0,29

0,75

0,07

0,07

fmax

19

9

39,9

27,6

19

19

Решение найдено. Общее число вычисляемых элементов матриц при решении задачи линейного программирования уменьшилось и равно 24 (вместо 72 в предыдущем случае).

1.4. Оперативное управление при изменении расхода ресурса

Исходные условия задачи на этапе планирования заданы соотношениями (6), (7).

Пример 4. На этапе управления оперативно скорректировать план, полученный в примере 1, при новом значении расхода ресурса А21 = 4 вместо запланированного ранее в (7) значения А21 = 2.

Решение системы линейных неравенств исходной задачи (7) имеет вид (8). Подставляя значения х1, х2 во второе измененное неравенство системы (7), получаем новое ограничение

0 Р1 + 15Р2 + 27,161Р3 + 48Р4 + 255,981Р5 + 229,567Р6 24 .                (13)

Таким образом, как и ранее, приходим к задаче линейного программирования с малым числом ограничений: найти максимум критерия (12) при наличии ограничений (10), (13).

Решение задачи имеет вид:

Р1 = 0; Р2 = 0; Р3 = 0,39; Р4 = 0; Р5 = 0; Р6 = 0,06.

При этом    х1=3,74 ,   х2=3,13 ,F=17,5,

т.е. изменение расхода ресурса В2 должно сопровождаться изменениями производственных заданий с учетом сложившейся производственной ситуации.

1.5. Подсистема формирования оперативных производственных           заданий

Автоматизированная система оперативного управления будет соответствовать прямому назначению только в том случае, если она сможет обеспечить регулирование хода производства путем оперативного решения комплекса задач.

Для примера пусть изменение плана х1, х2 потребует решения 5 взаимосвязанных задач. Решение задач необходимо для внесения изменений в технологический процесс и осуществления соответствующих организационно-технических изменений.

Граф взаимосвязи задач представлен на рис.1.

1                           3                          5

  1. 4

                    Рис.1. Граф взаимосвязи задач.

Примем для конкретности время решенияi-й задачи равным 10i.

Используя метод, основанный на вычислении длительности критического пути Ткр в графе, определяем, что наибольшее количество процессоров в вычислительной системе для реализации заданного набора задач за наименьшее время равно

  ,    где      ,

гдеi – время реализации задачи хi на одном процессоре. Для нашего примера Мкр = 1,67, т.е. число процессоров М = 2. Решение задачи календарного планирования для заданного комплекса задач можно выполнить с помощью решения задачи линейного программирования.

Введем обозначения:ti – время начала решенияi-й задачи,i=;t6 – время окончания решения всего комплекса задач.

Тогда задачу линейного программирования можно представить в следующем виде:  найти минимум критерияI =t6  при наличии ограничений

t3 – t1 10 ;   t5 – t3 30;

t5 – t2 20 ;   t4 – t2 20;

t6 – t5 50 ;   t6 – t4 40.

Решение задачи имеет вид

t1 = 0;t2= 20;t3 = 10;t4 = 40;t5 = 40;t6 = 90.

Распределение задач по процессорам осуществляется последовательно, начиная с наименьшего значенияti на свободные в данный момент процессоры, при этом приоритет имеет тот процессор, на который была распределена задача-аргумент.

Таким образом, задачи 1, 3, 5 решаются на первом процессоре, 2 и 4 – на втором процессоре.

  1. Оперативное управление на основе решения задач

нелинейного программирования

Рассматриваемый в пособии класс задач представим в следующем виде: найти минимум некоторой дифференцируемой нелинейной функцииn  переменных

J =f(x1,x2, …xn)min                                                                    (14)

при наличии ограничений вида (2).

Широкое распространение при решении подобных задач получил метод Давидона-Флетчера-Пауэлла (ДФП) [8] и выпуклый симплексный метод Зангвилла [8]. Рассмотрим кратко вначале алгоритм ДФП при отсутствии ограничений (2).

2.1. Метод ДФП

Начальный этап. Пусть - константа для остановки. Выбрать начальную точку х0(1,0) и начальную симметрическую положительно определенную матрицу Н(1,1), например единичную матрицу.

Положить х1(1,0) = х0(1,0),j = 1 и перейти к основному этапу.

Основной этап.Шаг 1. Если норма градиента функцииf(xj), то остановиться, в противном случае положитьdj = -Hjf(xj) и взять в качествеxj оптимальное решение задачи минимизацииf(xj +dj) при 0. Положитьxj+1 =xj +xjdj. Еслиjn, то перейти к шагу 2. Еслиj =n, то положитьx =xj+1 и далее повторить шаг 1 с новыми найденными начальными условиями. Таким образом, через каждыеn шагов производится восстановление начальной матрицыH1(1,1), что повышает устойчивость алгоритма.

Иногда в целях упрощения алгоритма матрицу Н(1,1) не пересчитывают.

Шаг 2. ПостроитьHj+1 следующим образом:

 ,                                                       (15)

гдеpj =jdj                                                                                                   (16)

qj =f(xj+1) -f(xj).                                                                               (17)

Заменитьj наj+1 и перейти к шагу 1.

При наличии ограничений типа равенств

                  А(1,1) Х(1,0) =b(1,0)                                                                       (18)

в методе ДФП вектор х0(1,0) и матрица Н1(1,1) определяются по формулам:

,                                (19)

                 .                  (20)

Отметим, что для построения Н1(1,1) и х0(1,0) требуется достаточно сложная операция – обращение матрицы [A(1,1)AT(1,1)].

В остальном алгоритм решения полностью соответствует алгоритму ДФП без ограничений. Иногда в целях упрощения алгоритма матрицу Н(1,1), как и в случае минимизации функции без ограничений, не пересчитывают, т.е. Н(1,1) = Н1(1,1).

Для случая ограничений (18) при дополнительном условии Х(1,0) 0 решение задачи поиска минимума (2) можно проводить с помощью метода ДФП. При этом нужно учитывать, что спускаясь по очередному направлениюdk,можно выйти из допустимой области. Чтобы этого не произошло, надо ограничить максимально возможное перемещение в данном направленииxj+1 =xj +jdj. Величинаj должна быть такой, чтобы обеспечивалось условиеxj+1 0. Проверка выполнения этого ограничения проще, чем проверка выполнения всех условий (2).

2.2. Оперативное управление при изменении целевой функции

На основе решений (3) метод ДФП позволяет повысить оперативность управления.

Пример 5. Минимизировать

I = 2x12 + 2x22 – 2x1x2 – 4x1 – 6x2(21)

при условиях

x1 +x2 +x3 = 2

x1 + 5x2 +x4 = 5                                                                        (22)

x1,x2,x3,x4 0                   .

Решим данную задачу методом ДФП.

По формуле (19), (20) находим

H1= .

X0T = [0,6; 0,867; 0,533; 0,067] .

Тогда  градиент функции и вектор направления будут равны

fT = [-3,33; -3,73; 0; 0];

d1T = [1,5; -0,196; -1,306; -0,52].

Вектор положения запишется в виде

X1 = 0,6 + 1,502222 ;

X2 = 0,866667 - 0,1955556 ;

X3 = 0,5333 - 1,306667 ;

X4 = 0,06667 - 0,524444 .

Из условия положительностиxj 0,j = 1,4.

Xmax = 0,1271203.

Минимальное значение функции (21) при данном направлении движения получается в точке = 0,1271203;x1 = 0,7905;x2 = 0,8418. В найденной точке минимума находим новые значения градиента функции и вектора направления

= [-2,5016; -4,1938; 0; 0];

2=-H1 ;d2T = [0,942; -0,054; -0,888; -0,672].

Новый вектор положения имеет вид

X1 = 0,7905 + 0,942;

X2 = 0,8418 - 0,054;

X3 = 0,367 - 0,888;

X4 = 0 - 0,672 .

Из последнего соотношения следует, что движение в направленииd2 невозможно, т.к. переменная х4 должна быть положительной величиной.

Примем х4 = 0 и произведем пересчет направления. При этом система ограничений (22) примет вид

х1 + х2 + х3 = 2                                                                   (22 а)

х1 + 5х2 + 0 х4 = 5              .

Вектор положения при этом запишется в виде

X1 = 0,7905 + 0,999 ;

X2 = 0,8418 - 0,2 ;

X3 = 0,367 - 0,799 ;

X4 = 0 + 0 .

Из условия положительностиxj 0,  находим .

Минимальное значение  функции (21) при данном направлении движения получится в точке = 0,459;X1 = 1,2495;X2 = 0,75;X3 = 0;X4 = 0.

ПримемX4 = 0 и произведем пересчет направления для найденной текущей точки минимума. Система ограничений остается прежней – (22 а).

Вектор положения при этом имеет вид

X1 = 1,2495 - 0,356;

X2 = 0,75 + 0,071;

X3 = 0 + 0,285;

X4 = 0.

Минимальное значение функции (21) получается в точке = 0,345855

При этом

X1 = 1,126376;

X2 = 0,774556;

X3 = 0,086;

X4 = 0;Imin = -7,1604,

что близко к точному оптимальному значениюI*=-7,16.

Для повышения оперативности управления матрицы Н1(1,1) для различных направлений, в том числе и при нулевых значениях,xi,xj могут быть рассчитаны заранее по формуле (19).

При измененииI подстройка начинается с последнего набораxi,i= . Матрица Н1, полученная по формуле (19), для данного направления позволит оперативно произвести подстройку оптимальной точкиxi,i=.

  1. Выпуклый симплексный метод Зангвилла

Рассмотрим алгоритм выпуклого симплексного метода Зангвилла для решения задачи минимизацииf(x) при условии А(1,1)Х(1,0) =b(1,0); Х(1,0) 0.

Начальный этап. Выбрать точку Х1, для которой АХ1 =b и Х1 0. Положитьk = 1 и перейти к основному этапу.

Основной этап.Шаг 1. При заданномXk определитьIk – множество индексовm наибольших компонент вектораXk;

B(1,1) = {ajjIk};N(1,1) = {ajjIk) ;                                        (23)

r(0,1) =f(Xk)T -f(Xk)TB-1 A ;                                                  (24)

= max{-rjrj 0};                                                                            (25)

= max{xjrjrj 0};                                                                           (26)

=        индекс, для которого = -r;                                              (27)

                       индекс, для которого =xr.(28)

Если = = 0, то остановиться. Если, то вычислитьdN в соответствии с (27) и (29):

             0, еслиjIk;j

dj =                                                                                                        (29)

            1, еслиjIk;j =           .

Если , то вычислитьdN по формулам (28), (30) :

               0, еслиjIk;j

dj =                                                                                                        (30)

              -1, еслиjIk;j =         .

Если 0, то вычислитьdN либо по формулам (27), (29), либо по формулам (28), (30). В любом случае определитьdB в соответствии с (31)

dB =B-1NdN                                                                                (31)

и перейти к шагу 2.

Шаг 2. Минимизироватьf(xk +dk) при условии 0max

где                                                                   (32)

Пустьk -  оптимальное решение этой задачи. Положитьxk+1 =xk +kdk, заменитьk наk+1 и перейти к шагу 1.

Пример 6. Решим задачу из примера 5 выпуклым симплексным методом Зангвилла. Возьмем в качестве начальной точкиx1 = (0;0;2;5)T. Заметим, что

f(x) = (4x1 – 2x2 – 4; 4x2 – 2x1 – 6; 0; 0)T.

Итерация 1

Поиск направления минимизации. В точке х1 имеемf(x1)={-4; -6; 0; 0}T.

Соответственно алгоритмуI1 = {3,4}, поэтому

,      а            .

Приведенный градиент вычисляется в соответствии с (24)

rT = (-4; -6; 0; 0)

В соответствии с (25) получаем

=max{-r1, -r2, -r3, -rj} = -r2= 6.   Из  (26) имеем =max{x3r3,x4r4} = 0 и, следовательно, из (27) получаем, что = 2.Построение направления проводится по формулам (29) и (31).

Из (29) имеемdNT = (d1,d2) = (0;1), а из (31) следует, что

dBT = {d3,d4} =  =(-1; -5).

Линейный поиск из точки х1(0,0,2,5)Т осуществляется по направлениюd1 = {0;1;-1;-5}T. Максимальное значение, для которого точкаx1 +d1 допустима, определяется по формуле (32):

Кроме того, имеемf(x1 +d1) = 22 - 6.. Следовательно, нужно решить задачу минимизации 22 - 6  при условии 0 1. Оптимальным решением является х =1, так что х2 = х1 +d1 = (0,1,1,0)T.

Итерация 2

Поиск направления минимизации.   В    точке   х2 = (0, 1, 1, 0)   множествоI2 = {2;3},    т.е.  ,      ,f(x2) = (-6; -2; 0; 0)T. При этом приведенный градиент (24)   .

В соответствии с   (25)   и   (26)   имеем =max{-r1, -r2, -r3} = -r1 = ; =max{x2r2;x3r3;x4r4} = 0,   т.е. = 1.

Согласно   (29),   (31)   имеемdNT = (d1,d4) = (1; 0) .

Тогдаd2 = (1;  -0,2; -0,8; 0).

Линейный поиск. Из точкиx2 = (0,1,1,0)T по направлениюd2 = (1, -0,2; -0,8; 0)T. Максимальное значение, для которого точка (x2 +d2) допустима, определяется по формуле (32):

.

Кроме того,f(x2 +xd2) = 2,48x2 – 5,6 - 4.

Следовательно, нужно решить следующую задачу: минимизировать 2,482 – 5,6 - 4 при условии 0

Оптимальное решение равноx2= , т.е.x3 =x2 +2d2 =.

Итерация 3

Поиск направления. В точкеx3 множествоI3 = (1,2) , т.е.  ,          ,         .

Из (24) получаемrT = (0;0;0;1). В этом случае =max{-r1, -r2, -r3} = 0.

= max(x1r1; x2r2; x3r3; x4r4) = 0.

Таким образом, точкаx3 оптимальна. Оптимальное значение целевой функции

I = -7,16.

Для практических расчетов целесообразно сочетать оба рассмотренных метода. Метод Зангвилла более прост для определения направления минимизации, а для поиска минимума в выбранном направлении удобен метод ДФП.

3. Оперативный расчет компромиссных планов

по нескольким критериям.

Многокритериальная оптимизация.

Производственные задачи, решаемые математическими методами, часто являются многокритериальными задачами. Так, производственный план предприятия должен обеспечивать максимально возможный выход продукции, малые издержки, высокую производительность труда, низкую себестоимость продукции, рентабельность производства и т.д.

Теоретические вопросы многокритериальной оптимизации рассмотрены в методических указаниях [9]. Покажем, как на основе решений системы линейных неравенств, изложенных в главе 1, можно повысить оперативность решения многокритериальных задач оптимизации.

3.1. Метод равных и наименьших относительных отклонений

Заданы два критерия  оптимизации

J1 = 3x1 + 2x2max,                                                                      (33)

J2 = 2x1 + 4x2max,                                                                      (34)

при наличии ограничений, общих для двух критериев, в виде системы (7). Используя систему (8), решаем последовательно для первого и второго критерия задачи на основе соотношений (4), (5).

При этом получимJ1max = 24,333;J2max = 30,958. Уравнение равнозначности [9, формула (4)] при этом имеет вид

0,0585х1 – 0,047х2 = 0 .                                                                    (35)

С учетом (8) это уравнение можно записать в виде

-2,335Р2 – 0,235Р3 + 0,702Р4 + 0,551Р5 + 2,725Р6 = 0 .                   (36)

Тогда компромиссное решение получим, решив задачу определения максимального значения, напримерJ2(34) при наличии ограничений (10), (36). Компромиссные значенияJ1 = 23;J2 = 29,2. Сокращение числа ограничений в последовательно решаемых задачах линейного программирования позволяет  повысить оперативность решения.

3.2. Метод обобщенной целевой функции

Решение для равнозначных критериев можно получить и на основе нахождения экстремума суммы взвешенных критериев

J =1J1 +2J2                                                                                  (37)

при наличии ограничения (10).

Важным элементом  при такой оптимизации является назначение коэффициентов веса каждого оптимизируемого критерия.

3.3. Критерии, различные по значимости.   Метод последовательных уступок

Исходные критерии ранжируются. В начале решается задача оптимизации по самому важному критерию без учета других критериев оптимизации. Пусть таким критерием будетJ1(33). Решаем задачу поиска максимумаJ1 при наличии ограничения (10). ПолучимJ1max = 24,3. Затем решаем задачу оптимизации по второму критериюJ2(34) при наличии ограничений (10) иJ1 J1уст, гдеJ1уст – величинаJ1, до которой может быть уменьшена величина первого критерия. Например,J1уст = 20. При этомJ2max = 30,4.

Учитывая повышение оперативности решения, можно одновременно повысить уровень информационной поддержки для принятия управленческого решения, решив задачу при различных величинах уступок. В табл.14 представлена зависимостьJ2 = f(J1уст).

                                                                                   Таблица 14

J2

20

22

23

J1уст

30,4

29,6

29,1

Прогнозирование подобных зависимостей можно осуществлять с помощью регрессионных моделей [10]. Данный вопрос для повышения его практической значимости рассмотрен в разделе 4 для случая, когда в результатах измерений присутствуют аномалии. В разделах 5,6 рассмотрены вопросы обработки измерительной информации для оперативного управления.

4. Оперативное оценивание параметров модели регрессии при наличии аномальных результатов измерений

В реальных системах обработки информации оценки  вектора регрессионных коэффициентовb для модели

       ,                                 (38)

гдеX[m] – вектор наблюдаемых линейно независимых факторов;

b – вектор неизвестных и подлежащих оценке параметров;

[m] – помеха типа белого шума, приходится проводить в условиях аномальных измерений (АИ)Y(l),l ().

Наибольшее распространение при решении поставленной задачи получили метод максимального правдоподобия при известном законе распределения ошибки [10] и метод наименьших модулей (МНМ), обеспечивающий устойчивое решение в условиях отклонения реального закона распределения ошибки от постулируемого априори закона распределения [11,12]. Однако все эти методы требуют в случае обнаружения АИ значительных вычислительных затрат для исключения влияния самих измерений на искомую регрессионную зависимость.

В пособии предлагается организовать процесс вычислений таким образом, чтобы не получать каждый раз результат заново, а корректировать процесс вычислений с учетом аномальности результата измерений. Подобный подход можно применять для обработки как накопленных данных, так и последовательно поступающей порции информации.

4.1. Метод выделения результата АИ

В основу его алгоритма положен итеративный метод решения на основе МНМ [11,13]. При этом оценка , полученная на основеN результатов измерения, имеет вид

                                                                   (39)

где

                       (40)

                                                       (41)

где   - оценкаm-го измерения выходного сигнала.

Начальные значенияR[m]=1;  соответствуют определению параметров  по методу наименьших квадратов (МНК). Далее вычисления оценок по формулам (39) – (41) проводятся итерационно до тех пор, пока изменения оценок за одну итерацию не достигнут заданной малой величины. При этом наименьший весовой коэффициентR[l] указывает на наиболее грубоеl-измерение.

4.2. Рекуррентная процедура исключения АИ

Матрицу  и векторZn-1[n] приR[m]=1;  можно определить из ,ZN[n], полученных на первой итерации, путем исключения аномальных составляющих

                                         (42)

                                                                   (43)

На основании леммы об обращении матриц можно записать

                                       (44)

Используя выражения (43) и (44), получаем

                            (45)

4.3. Рекуррентная процедура учета текущего неаномального измерения

Эта процедура соответствует обычному рекуррентному МНК с нарастающей памятью.

Рекуррентные расчетные соотношения имеют вид:

При большом числе измерений итерационную процедуру проверки на аномальность очередного измерения следует начинать с весовой функции

В ряде случаев достаточно проверки

Величина определяется статистическими характеристиками помехи [10].

4.4. Построение модели регрессии при наличии сильной корреляции независимых факторов

Наличие корреляции независимых факторов в модели (38) и определение оценок вектора  на основе соотношения (39) приводят к неустойчивости решения.

В этом случае при выполнении операции обращения матрицы в выражении (39) следует применять псевдообращение [14]:

                                                                      (46)

Дальнейшая итерационная процедура выделения результата АИ, описанная выше, остается без изменений. Процедура требует настройки алгоритма псевдообращения на величину грубости псевдообращения. Практическая проверка алгоритма показала, что некоторое увеличение степени загрубления способствует сокращению числа циклов итерационного процесса выявления АИ.

5. Оперативное оценивание параметров входного сигнала измерительной системы при наличии помех

   В настоящее время в связи с повышением требований к оперативности управления предельно актуальной является проблема оперативного наиболее точного восстановления и преобразования искаженного сигнала измерительной системы [15]. Даже если исключить внешние помехи, сами измерительные системы искажают входной сигнал.

5.1. Оптимальное восстановление входного сигнала

В данном разделе рассматриваются помехоустойчивое восстановление  и преобразование сигнала, искаженного линейной инерционной системой. Примером такой системы может служить любой инерционный датчик. А поскольку именно датчики являются основным источником информации для работы автоматических и автоматизированных систем, проблема оперативного восстановления и преобразования сигнала встает особенно остро.

Одним из возможных способов повышения устойчивости и точности оценок входного сигнала является сужение класса функций восстанавливаемого входного сигнала.

В данном пособии принята на интервале наблюдения модель входного сигнала полиномиального типа

.                                                                                  (47)

В качестве моделей линейной динамической измерительной системы приняты передаточные функции вида [16]:

                                       (48)

     ,

с помощью которых могут быть описаны различные классы измерительных систем. Расчетная схема формирования наблюдаемого сигнала системы оптимальной коррекции динамических погрешностей измерений системой W(p) представлена на рис. 2.

      Рис. 2. Схема оптимальной фильтрации при наличии помехи типа                    «белый шум»

Необходимо синтезировать линейные фильтры Wi, i=1,2, выделяющие с наименьшей дисперсией ошибки из выходного сигнала измерительной системы YH(t) составляющие, пропорциональные контролируемым параметрам b1,b2,…,bn. Сигналы Z1(t), Z2(t) (рис. 1) являются случайными из-за случайности параметров b1,b2. Желаемым выходным сигналомj(t) синтезируемого линейного фильтра Wj является сигнал

.                                                                                          (49)

Сформируем ошибку при определении параметра bj с помощью фильтра Wj в виде

                                                        (50)

и определим оператор Wj из условия минимума дисперсии этой ошибки

.                                  (51)

Рассмотрим вначале решение данной задачи для случая помехи n(t) типа белого шума с корреляционной функцией Rn() =*(). Затем на основе данного решения произведем обобщение на случай произвольной помехи.

Оптимальные несмещенные оценки  в случае помехи типа белого шума определяются соотношениями:

.                                                                                        (52)

Элементы матрицы А и вектора  равны

(53)

.                                                                           (54)

Точность оценок характеризуется дисперсионной матрицей ошибок

,                                              (55)

где – интенсивность белого шума.

В случае дискретных измерений выходного сигнала измерительной системы  в моменты t=nT уравнения (53), (54) перепишутся в виде

,                                                                       (56)

.                                                                       (57)

Синтезируем на основе рассмотренной схемы рис. 1 оптимальный фильтр в случае помехи, отличной от белого шума.

Преобразуем расчетную схему, представленную на рис.2. Помеха n(t) характеризуется нулевым математическим ожиданием и корреляционной функцией Rn(). Будем искать оптимальный фильтр в виде последовательного соединения двух звеньев1 и2j, j=1,2. Звено1такое, что преобразует стационарную произвольную помеху n(t) в стационарный белый шум. Эквивалентная расчетная схема представлена на рис. 3.

Рис. 3. Схема оптимальной фильтрации при помехе, отличной от белого шума

Здесь1-1 представляет собой формирующий фильтр для помехи n(t) c корреляционной функцией Rn(). Окончательно схему (рис.3) можно представить в виде рис. 4.

Рис. 4. Схема оптимальной фильтрации при помехе, отличной от белого шума

Таким образом, предлагаемая выше методика синтеза оптимальных фильтров на основе соотношений (52)-(57) может быть применена и для случая синтеза фильтров при произвольной помехе.

5.2. Спектрально–временной метод восстановления сигнала.

Данный метод разработан научной школой В.В. Солодовникова [17]. Сущность этого метода сводится к нахождению импульсной переходной функции восстанавливающего фильтра. Так как в чистом виде дифференциатор – нереализуемое звено, то в данном методе используется приближенная реализация дифференцирующего звена, например в виде передаточной функции:

.                                                                        (58)

Подбирая, можноW(p) приблизить к дифференцирующему звену в заданной полосе частот. Отсюда получим

.                                                                      (59)

РаскладываемW(p) в ряд, например, по экспонентам

.                                                     (60)

Обозначим     .                                                               (61)

Из формулы (60) найдем коэффициенты С. Пусть , . Тогда,перейдя к матричной форме записи, С найдется как

.                                                                     (62)

Зная коэффициенты С, можно найти импульсную переходную функцию восстанавливающего фильтра по формуле (63):

.                                                                        (63)

На  рис.5 изображен график импульсной переходной функции восстанавливающего фильтра.

Рис. 5. Импульсная переходная функция восстанавливающего фильтра

Из рисунка видно, что импульсная переходная функция получилась конечной. Поэтому интеграл свертки (64) можно брать в реальном времени, т.е. на конечном интервале времени, меньшем времени переходного процесса.

.                                                         (64)

Это большое достоинство метода.

Однако, спектрально–временной метод основан на приближенной реализации дифференцирующего звена. Операция дифференцирования при наличии большого уровня помех приводит к сильным зашумлениям выходного сигнала по сравнению с методом, изложенным в разделе 5.1.

6. Повышение оперативности решения задачи восстановления сигнала на основе применения БПФ-свертки

Решение задачи восстановления входных сигналов [18] сводится к решению интегрального уравнения вида

,                                                                (65)

где К(t) – импульсная переходная функция стационарного оператора,x(t) – искомый входной сигнал. В большинстве случаев линейные операторы измерительных систем имеют интегрирующий характер [19].

Динамические характеристики восстанавливающих операторов, применяемых для  преобразования сигналаy(t), приближаются к характеристикам неидеальных дифференцирующих операторов, имеющих укороченные по длине импульсные переходные функции. Входной сигнал может быть восстановлен по формуле:

         .                                                                     (66)

Еслиn изменяется в больших пределах, то вычисление по формуле (66) может оказаться очень трудоемкой операцией по числу операций сложения и умножения. Некоторое упрощение можно получить, если применить БПФ-свертку без двоичных инверсий. Однако, импульсная переходная функцияK[n] дифференцирующих звеньев принципиально не может быть реализована точно. Попытку избежать этой принципиальной нереализуемостиK[n], а с ней и точного восстановления входного сигнала можно осуществить на основе итерационных процедур Ван-Циттерта [20]. Для уравнения (65) итерационная процедура может быть представлена в виде

                             (67)

Выбор величин обеспечивается сходимостью в среднем к решению уравнения (65). В выражении (67) импульсная переходная функция принципиально может быть реализована с высокой точностью. Некоторого упрощения итерационной процедуры можно добиться, если, как и раньше в случае восстановления по формуле (66), применять БПФ-свертку без двоичных инверсий. Ниже приведены примеры применения БПФ-свертки для восстановления входного сигнала с помощью соотношения (66) путем перехода в частотную область и примеры восстановления входного сигнала на основе дискретизации соотношения (67).

Пример 7. Восстановление входного сигнала интегратора.

Входной сигнал в виде единичного скачка:Xi=1, гдеi=N – число отсчетов. Для нашего примераN=16. Тогда на выходе интегратора будем иметь сигнал вида

YT=[ 0 1 2 3 4 5 6 7 8 9 10 0 0 0 0 0 ].

Дифференцирующий оператор:

VT=[ 1 –1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ].

Находим коэффициенты разложения в ряд Фурье для сигналовV иY. Матрица дискретного преобразования Фурье имеет вид:

где         .

Для восстановления входного сигнала на интервале наблюдения, необходимо наблюдаемый входной сигнал и дифференцирующий оператор дополнить соответствующим количеством нулей [18].

Коэффициенты разложения сигнала имеют вид:

Тогда спектр восстановленного сигнала имеет вид:XFj=VFjYFj,

Применяя обратное преобразование Фурье, получаем:  гдеB(l,m)=W-lm.

=[ 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 ].

Для сокращения числа операций целесообразно применять БПФ, представленный в матричной форме. Матрица преобразований имеет вид:

,                        (68)

гдеIn - единичная матрица размером (nn),Aj=A1A2An- кронекеровское произведение,Ai=A1A2An - прямая сумма.

     ,

где ,N = 2n - число отсчетов.

Формула (68) дляn = 4 представляется в виде :

                               (69)

где        .

i = 2,4,6,8j = 1,3,5,7

Результаты представленыFKP с двоично-инверсными номерами. С помощью формул (68) и (69) находим спектр выходного сигнала и спектр восстанавливающего фильтра. Перемножаем  их с двоично-инверсными номерами (не переставляя индексов). Входной сигнал находим по формуле :

                        (70)

ФункцииFKR в формулах (69) и (70) имеют двоично-инверсный порядок. Благодаря отсутствию операции с преобразованием индексов получим экономию в вычислительных операциях.

Пример 8. Восстановление входного сигнала интегрирующего звена  с помощью итерационных процедур.

Дифференцирующие звенья принципиально не могут быть реализованы с большой точностью, интегрирующие инерционные звенья принципиально могут быть реализованы с любой заданной точностью. Для восстановления входного сигнала целесообразно использовать итерационную процедуру:

                    (71)

где Т (шаг интегрирования) и определяются точностью интегрирования, К – матрица преобразования интегрирующего оператора,S = 0…7,M = 0…7,I = 1…2000 (число итераций). За счет выбора малого значения= 0,1     обеспечиваем сходимость итерационной процедуры.

Начальные значения входного сигнала принимаются равными выходному сигналу.

Входной сигнал:XT = [1  2  3  4  3  2  1  1].

Выходной сигнал:YT = [0,1  0,3  0,6  1  1,3  1,5  1,6  1,7].T=0,1,=0,1.

,

Пример 9. Восстановление входного сигнала апериодического звена

Для восстановления входного сигнала используем итерационную процедуру (71). Начальные значения входного сигнала принимаются равными выходному сигналу. Импульсная переходная функция рассчитывается по формулеK(t) =e-t.     Матрица преобразования имеет вид

Входной сигнал ХТ = [1   2   3   4   3   2   1   1].

Выходной сигнал:YT = [0,1  0,29  0,563  0,909  1,123  1,216  1,2  1,186].T=0,1, = 0,1.

Восстановленный входной сигнал после 2000 итерации

= [1   2   3   4   3   2   1   1].

Библиографический список

1. Косяченко С .А., Кузнецов Н.А., Кульба В.В., Шелков А.Б. Модели, методы и автоматизация управления в условиях чрезвычайных ситуаций//АиТ. 1998. № 6. С.3-15 .

2. Первозванский А.А. Математические модели в управлении производством. М.: Физматгиз, 1975. 616 с.

3 Курицкий Б.Я. Поиск оптимальных решений средствамиExcel 7.0. Санкт-Петербург.BHV, 1977. 384 с.

4. Милов Л.Т. Многомерно-матричные производные и анализ чувствительности систем автоматического управления//АиТ.1979.№ 9.С.15-25.

5. Кабанов А.Н. Математические методы повышения устойчивости алгоритмов АСУ. Учебное пособие. Рязань: РГРТА, 1985. 80 с.

6. Муртаф Б. Современное линейное программирование. М.: Мир, 1984. 224 с.

7. Линейное и целочисленное программирование: Методические указания /Рязан.гос.радиотех.акад.; Е.М. Дондик и др. Рязань, 1998. 44 с.

8. Базара М., Шетти К. Нелинейное программирование. Теория и алгоритмы. М.: Мир, 1988. 583 с.

9. Многокритериальная оптимизация: Методические указания/ Рязан.гос.радиотех.акад.; А.Н. Кабанов , Н.И. Тарасова  Рязань, 1997. 12 с.

10. Тихонов А.Н., Уфимцев М.В. Статистическая обработка результатов экспериментов. М.: Изд-во Моск. ун-та, 1988. 174 с.

11. Мудров В.И., Кушко В.Л. Методы обработки измерений. М.: Радио и связь, 1983. 304 с.

12. Гильбо Е.П., Челпанов И.Б. Обработка сигналов на основе упорядоченного выбора. М.: Сов. радио, 1976. 344 с.

13. Перельман И.И. Оперативная идентификация объектов управления. М.: Энергоиздат, 1982. 272 с.

14. Альберт А. Регрессия, псевдорегрессия и рекуррентное оценивание. М.: Наука,1977. 244 с.

   15. Клейман Е.Г. Идентификация входных сигналов в динамических системах// АиТ. 1999. №12.  С. 3-15.

16. Аш Ж. И др.  Датчики измерительных систем. В 2-х книгах. Кн.2: Пер с франц. М.: Мир, 1992. С. 292.

  1. Солодовников В.В., Семенов В.П., Пешель Н., Недо Д. Спектральный и интерполяционный методы. М.:Машиностроение,1978. 664 с.
  2. Василенко Г.И. Теория восстановления сигналов: о редукции к идеальному прибору в физике и технике. М, 1979. 272 с.
  3. Грабовский В.А. Динамические измерения:Основы метрологического обеспечения. Л.: Энегроатомиздат,1984. 224 с.
  4. Преображенский Н.Г., Пикалов В.В. Неустойчивые задачи диагностики плазмы. Новосибирск: Наука, 1982. 236 с.

Оглавление

Введение ………………………………………………………..………..…  1

  1. Оперативное управление на основе решения задач

     линейного программирования ………………………………………….…  1

1.1. Общая схема решения задачи линейного программирования

      на основе определения фундаментальной системы решений

      для системы однородных линейных неравенств ………………………..2

1.2. Оперативное управление при изменении коэффициентов

      целевой функции ……………………………………………………..……3

1.3. Оперативное управление при изменении величины ресурсов ………....6

1.4. Оперативное управление при изменении расхода ресурса …………….8

  1. Подсистема формирования оперативных производственных заданий...9
  1. Оперативное управление на основе решения задач

      нелинейного программирования …………………………………………10

2.1. Метод ДФП ………………………………………………………………. 10

2.2. Оперативное управление при изменении целевой функции …………..12

2.3. Выпуклый симплексный метод Зангвилла …………………………….. 14

  1. Оперативный расчет компромиссных планов по нескольким

критериям. Многокритериальная оптимизация .…………………… …17

3.1. Метод равных и наименьших относительных отклонений ………… .17

3.2. Метод обобщенной целевой функции …………………………………18

3.3. Критерии различны по значимости. Метод последовательных уступок……………………………………………………………….…..18

4.   Оперативное оценивание параметров модели регрессии при наличии   аномальных результатов измерений ……………………………………19

4.1. Метод выделения результата АИ ………………………………………19

4.2. Рекуррентная процедура исключения АИ ……………………………..  20

4.3. Рекуррентная процедура учета текущего неаномального  измерения ..20

4.4. Построение модели регрессии при наличии сильной   корреляции независимых факторов ………………………………………………….21

   5.   Оперативное оценивание параметров входного сигнала        измерительной  системы при наличии помех …………………………22

5.1. Оптимальное восстановление входного сигнала ………………………22

5.2. Спектрально–временной метод восстановления сигнала …………..…26

  1. Повышение оперативности решения задачи восстановления  сигнала

на основе применения БПФ-свертки …………………………………….27

Заключение …………………………………………………………………….31

Бибилиографический список …………………………………………………  31




Похожие работы, которые могут быть Вам интерестны.

1. Разработка модели и решение задачи линейного программирования (на примере задачи о составлении графика работы персонала)

2. Создание макросов для решения задач VBA

3. ГИБРИДНАЯ ИММУННАЯ СЕТЬ ДЛЯ РЕШЕНИЯ ЗАДАЧ СТРУКТУРНОЙ ИДЕНТИФИКАЦИИ

4. Разработка приближенно-аналитических методов решения задач теплообмена

5. РАЗВИТИЕ ЛОГИЧЕСКОГО МЫШЛЕНИЯ МЛАДШИХ ШКОЛЬНИКОВ НА УРОКАХ МАТЕМАТИКИ ПОСРЕДСТВОМ РЕШЕНИЯ ЛОГИЧЕСКИХ ЗАДАЧ

6. ОЦЕНКА КОЛЛЕКТОРНЫХ СВОЙСТВ ТОНКОСЛОИСТЫХ НЕФТЯНЫХ ПЛАСТОВ НА ОСНОВЕ РЕШЕНИЯ ОБРАТНОЙ ЗАДАЧИ ДИНАМИЧЕСКОЙ ИНВЕРСИИ

7. Оперативное лечение опухолей молочных желез

8. Цикличность планирования. Планирование: оперативное, тактическое, стратегическое

9. Решение линейного уравнения численным методом

10. Анализ линейного, пассивного проходного четырехполюсника