39

Московский Авиационный Институт

(Технический Университет)

Кафедра 308

Курсовая работа

Выбор параметров контроля с использованием метода динамического программирования и метода ветвей и границ

Вариант II(2)

Выполнила

студентка

группы КТ-515

Принял

Москва

2008г.

Содержание

Задание

1. Метод динамического программирования

1.1 Теоретическая часть

2.2 Практическая часть

- ручной счёт

- листинг программы

2. Метод ветвей и границ

2.1 Теоретическая часть

2.2 Практическая часть

- ручной счёт

- листинг программы

Вывод

Литература

Задание

Вариант II(2)

Выбор параметров контроля с использованием метода динамического программирования и метода ветвей и границ при непересекающихся элементах объекта контроля и ограничениях по затратам на контроль С?16.

Исходные данные: вероятность отказов элементов и затраты на контроль параметров.

Выбрать такие параметры, чтобы С?16 при Q=Q_max.

1. Метод динамического программирования

1.1 Теоретическая часть

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

Пусть работоспособность объекта контроля характеризуется совокупностью n взаимосвязанных параметров, образующих множество S={x₁, x₂, …, x_n}. Проверка всех параметров из S влечет контроль всех N элементов системы и дает однозначный ответ: объект исправен, если все N элементов исправны, или неисправен, если по крайней мере один из элементов отказал. Для x_i определено подмножество R(x_i) элементов, проверяемых при контроле i-го параметра, причем предполагаем, что эти подмножества могут пересекаться, т.е. i, j: R(x_i)R(x_j). Пусть - некоторый набор параметров из множества S, т.е. S. Тогда и S. Значения x_i из S можно представить булевым вектором, причем

x_i = 1, если x_i,

0, если x_i.

Задача выбора параметров в этом случае формулируется двояко:

1) найти набор ?, для которого

P(?)=max

при ?x_i·c(x_i)?C; iЄ?

2) найти набор ?, для которого

?x_i·c(x_i)=min

при P(?)?P_з,

где P(?) - апостериорная вероятность работоспособного состояния объекта контроля при положительном исходе контроля выбранных параметров S; с(x_i) - затраты на контроль i-го параметра; Р_з - требуемая достоверность контроля; С - ограничение на общую стоимость контроля.

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

P(?)=Р₀/1-?Р_i,

iЄR(?)

где Р₀=?(1-р_i) - априорная вероятность безотказной работы объекта:

iЄR(S)

Р₀=1-?Р_i;

iЄR(S)

Р_i - нормированная вероятность отказа системы из-за отказа i-го элемента:

Р_i=(p_i/(1-p_i))/(1+? p_k/(1-p_k);

kЄR(S)

p_i - априорная вероятность отказа i-го элемента. Тогда вероятность того, что отказ будет обнаружен при проверке k-го параметра, можно вычислить по формуле:

Q_k=?P_k

kЄR(x_k)

При возможности наличия в ОК произвольного числа отказов

P(?)=?(1-p_i)/?(1-p_i)

iЄR(S) iЄR(?)

Можно использовать простой перебор вариантов, однако возникающие при этом вычислительные трудности не позволяют сделать этого даже для простых систем (при n>10). В связи с этим комплектование набора будем трактовать как многошаговый процесс, состоящий из последовательного выбора отдельных параметров.

В соответствии с общим принципом оптимальности разобьем весь имеющийся ресурс стоимости С на С отрезков единичной длины. (В практических случаях заданные положительные величины с(x_i) и С можно считать всегда целыми. Если это не так, то необходимо перейти к более мелким стоимостным единицам в зависимости от разрядности дробной части.). Рассмотрим наряду с интересующей нас исходной задачей множество аналогичных задач

f(Y)=max ?(x), Y Є [0,C],

xЄX_Y

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

Пусть ?₀=min c(x_i).

i=1,…,n

Тогда при всех ? Є [0,?₀] соответствующие множества ?_? состоят, из одного нулевого элемента и f(Y)=0 для всех таких ?. Для ресурса ? Є [?₀, С] согласно общей схеме динамического программирования справедливы следующие рекуррентные соотношения:

f(Y_k)=max [Q_i + f[Y_k - c(x_i)] - G_i (1)

iЄI_Y

где k=Y₀, Y₀+1, …, C; I_Y - множество тех i, для которых с(x_i)?Y_k, начиная с номера k=max c(x_i) уравнение (1) решается для всех i= 1,…,n;

G_i = ?P_i - сумма вероятностей элементов i-го параметра, которые пересекаются с

IЄR(x_i)??_l^*

элементами подмножества ?_l^*, образованного на шаге Y_k - c(x_i).

Если i, j; R(x_i)?R(x_j)= , то G_i=0 и

f(Y_k)=max {Q_i + f[Y_k - c(x_i)]} (2)

iЄI_Y

Для решения интересующей нас задачи опишем простой численный метод, не требующий предварительного определения всех допустимых наборов и основанный на рекуррентных соотношениях (1). Для всех целых ? = ?₀, С по формуле (1) вычисляются величины f(Y_k) и при этом фиксируются индексы i_Yk^*, на которых достигаются максимумы в (1). Искомый вектор ? формируется последовательно включением в набор параметра i_Yk и подмножества ?_l^*, зафиксированного на шаге Y_k - c(x_i). При этом, если Y_kЄ ?_l*, то на данном шаге этот параметр исключается из рассмотрения, так как каждый параметр может включаться в набор не более одного раза. Если на некотором ?-м шаге окажется, что f(Y_?)< f(Y_?-1), то в качестве ?_?* принимается подмножество ?_?-1* и фиксируется параметр i_{Y ?-1}, причем за f(Y_?)< принимается значение f(Y_?-1). Заметим, что если в задаче P(?)=max при

?x_i·c(x_i)?C

iЄ?

принять более жесткое ограничение, а именно ?c(x_i)=C, то последнее не допустимо, iЄ? так как в этом случае max f(Y_k) может быть меньше max f(Y_k-1) из-за того, что он достигается на другом подмножестве параметров.

Общая сложность метода, очевидно, ?(n) ? c(n+1), т.е. экспоненциальная функция при переборе заменена линейной функцией. При этом для запоминания промежуточных значений необходимо k?2c ячеек памяти. Если в качестве максимизируемого критерия использовать P(?)=?(1-p_i)/?(1-p_i), то необходимо решить задачу динамического iЄR(S) iЄR(?) программирования с мультипликативным критерием. Для этого достаточно прологарифмировать это выражение и обозначить

V=lgP(?)=lgР₀-?lg(1-p_i). (3)

iЄR(?)

Так как выражение, стоящее под знаком ? в (3), отрицательно, то, V= V_max тогда, когда максимальна величина суммы, т.е. в этом случае получим новую целевую функцию

V=??_i, где ?_i=lg (1-p_i),

iЄR(?)

обладающую свойством аддитивности и обращающуюся в максимум одновременно с P(?).

1.2 Практическая часть

Ручной счёт

Данные для расчета:

С?16

Таблица 1

Для удобства расчетов проранжируем таблицу1 следующим образом:

Таблица 2

Вычисления сведем в таблицу 3:

Таблица 3

Оптимальный набор включает параметры ?^*= {1,2,7,3,4,10,9} при этом

P(?) = 0,61+0,16 = 0,77 и С = 16.

Листинг программы

unit Unit1;

interface

uses

Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,

Dialogs, ToolWin, ComCtrls, mdCONTROLS, Grids, StdCtrls, ExtCtrls, Unit2,

Buttons;

type

TForm1 = class(TForm)

sgH: TStringGrid;

sgP: TStringGrid;

sgC: TStringGrid;

sgQ: TStringGrid;

lbC: TLabeledEdit;

BitBtn1: TBitBtn;

Label1: TLabel;

sgW: TStringGrid;

Label2: TLabel;

procedure FormCreate(Sender: TObject);

procedure BitBtn1Click(Sender: TObject);

procedure sgExit(Sender: TObject);

private

{ Private declarations }

public

H: TH;

P: TP;

C: TC;

W: TW;

end;

var

Form1: TForm1;

implementation

{$R *.dfm}

procedure TForm1.FormCreate(Sender: TObject);

var i,j: integer;

x: Byte;

f: TextFile;

begin

AssignFile(f, data.txt);

Reset(f);

sgW.Cells[0,0] := W;

// Ввод исходной матрицы

readln(f);

for j:=1 to 10 do

begin

sgH.Cells[0,j] := IntToStr(j);

sgW.Cells[0,j] := IntToStr(j);

for i:=1 to 10 do

begin

sgH.Cells[i,0] := IntToStr(i);

read(f, x);

sgH.Cells[i,j] := IntToStr(x);

if x = 1 then

H[i-1,j-1] := true

else

H[i-1,j-1] := false;

end;

readln(f);

end;

// Ввод вероятностей

readln(f);

readln(f);

sgP.Cells[0,0] := P;

for i:=1 to 10 do

begin

read(f, P[i-1]);

sgP.Cells[i,0] := FloatToStr(P[i-1]);

end;

readln(f);

// Ввод стоимостей

readln(f);

readln(f);

sgC.Cells[0,0] := C;

for j:=1 to 10 do

begin

read(f, C[j-1]);

sgC.Cells[0,j] := IntToStr(C[j-1]);

end;

CloseFile(f);

// Ввод вероятностей обнаружения отказа

sgQ.Cells[0,0] := Q;

for j:=1 to 10 do

sgQ.Cells[0,j] := FloatToStr(Q(j-1,H,P));

lbC.Text := 1;

end;

procedure TForm1.BitBtn1Click(Sender: TObject);

var i: integer;

begin

label1.Caption := FloatToStr(maxf(1, StrToInt(lbC.Text), H,P,C, W));

for i:=1 to 10 do

begin

sgW.Cells[2,i] := IntToStr(W[i-1].N);

if W[i-1].E then

sgW.Cells[1,i] := 1

else

sgW.Cells[1,i] := 0;

end;

end;

procedure TForm1.sgExit(Sender: TObject);

var i,j: integer;

begin

for j:=1 to 10 do

for i:=1 to 10 do

if sgH.Cells[i,j] = 1 then

H[i-1,j-1] := true

else

H[i-1,j-1] := false;

for i:=1 to 10 do

P[i-1] := StrToFloat(sgP.Cells[i,0]);

for j:=1 to 10 do

C[j-1] := StrToInt(sgC.Cells[0,j]);

// Ввод вероятностей обнаружения отказа

for j:=1 to 10 do

sgQ.Cells[0,j] := FloatToStr(Q(j-1,H,P));

end;

end.

unit Unit2;

interface

type

TH = array [0..9, 0..9] of boolean;

TP = array [0..9] of extended;

TC = array [0..9] of integer;

TDateW = record

E: boolean;

N: integer;

end;

TW = array [0..9] of TDateW;

function Q(j: integer; H: TH; P: TP): extended;

function maxf(n, Yk: integer; H: TH; P: TP; C: TC; var W: TW): extended;

implementation

function Q(j: integer; H: TH; P: TP): extended;

var i: integer;

begin

Result := 0;

for i:=0 to 9 do

if H[i,j] then

Result := Result + P[i];

end;

function G(j: integer; H: TH; P: TP; W: TW): extended;

var i,k: integer;

begin

Result := 0;

for i:=0 to 9 do

if H[i,j] then

for k:=0 to 9 do

if W[k].E and H[i,k] then

begin

Result := Result + P[i];

Break;

end;

end;

function f(n, Yk, j: integer; H: TH; P: TP; C: TC; var W: TW): extended;

begin

Result := Q(j,H,P) + maxf(n+1, Yk - C[j], H,P,C, W) - G(j,H,P,W);

end;

function maxf(n, Yk: integer; H: TH; P: TP; C: TC; var W: TW): extended;

var j,i: integer;

ft: extended;

Wt: TW;

begin

Result := 0;

for i:=0 to 9 do

begin

W[i].E := false;

W[i].N := 0;

end;

for j:=0 to 9 do

if C[j] <= Yk then

begin

for i:=0 to 9 do

begin

Wt[i].E := false;

Wt[i].N := 0;

end;

ft := f(n, Yk, j, H,P,C, Wt);

if Result < ft then

begin

Result := ft;

W := Wt;

W[j].E := true;

W[j].N := n;

end;

end;

end;

end.

2. Метод ветвей и границ

2.1 Теоретическая часть

Рассмотрим следующую задачу целочисленного программирования. Требуется максимизировать выражение:

n

L=?c_j•x_j (4)

j=1

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

n

?а_ij•x_j?b_i, i=1, …,m (5)

j=1

x_jЄ{0;1}, j=1, …,n

причем с_j?0, a_ij?0.

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

Эффективность вычислительных алгоритмов зависит от точности и простоты способа определения верхней границы возможных решений и точности определения приближенного решения. Чем точнее способ определения верхней границы целевой функции, тем больше бесперспективных ветвей отсекается в процессе оптимизации. Однако увеличение точности расчета верхних границ связано с возрастанием объема вычислений. Например, если для оценки верхней границы использовать симплекс-метод, то результат будет достаточно точным, но потребует большого объема вычислительной работы.

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

Обозначим: U - множество переменных x_j; S - множество фиксированных переменных, вошедших в допустимое решение; E_s - множество зависимых переменных, которые не могут быть включены в множество S, так как для них выполняется неравенство

а_ij> b_i - ?а_ij•x_j, i=1, …,m

x_jЄS

G_S - множество свободных переменных, из которых производится выбор для включения в S очередной переменной.

Рассмотрим первоначально одномерную задачу, когда m=1 и задача (4) имеет только одно ограничение вида (5).

Обозначим h_1j = c_j/a_1j и допустим, что x_jЄS (j=1, …,k<n) и выполняются условия

h_1,k+1? h_1,k+2? …? h_1l, l?n,

l

?a_1j>b₁ - ? a_1j•x_j,

j=k+1 x_jЄS

l-1

?a_1j? b₁ - ? a_1j•x_j,

j=k+1 x_jЄS

Условия означают, что во множество S без нарушения неравенства (5) можно дополнительно ввести элементы x_k+1, x_k+2, …, x_l-1. При введении во множество S элементов x_k+1, x_k+2, …, x_l неравенство (5) не выполняется.

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

H_s =?c_j•x_j + L_s,

x_jЄS

где

l-1

L_s = ?с_j + h_1l? b₁ ,

j=k+1

l-1

? b₁ = b₁ - ? a_1j•x_j - ?a_1j.

x_jЄS j=k+1

Из условий следует, что L_s не меньше максимального значения величины

?c_j•x_j

x_jЄG_S

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

? a_1j •x_j b₁ - ? a_1j•x_j = b₁,

x_jЄG_S x_jЄS

x_jЄ {0;1}, x_jЄG_S ,

Выбор очередной переменной для включения во множество S производится с помощью условия

h_1r(x_r) = max h_1j(x_j)

x_jЄG_S

Для выбранной переменной x_r определяются величины H_s(x_r) и H_s(x_r), т.е. в S включаются x_r = 1 или x_r = 0.

Если в процессе решения окажется, что во множестве G_S нет элементов, которые могут быть введены во множество S без нарушения ограничения (5), то полученное решение L_s =?c_j•x_j принимается в качестве первого приближенного x_jЄS решения L₀.

Все вершины дерева возможных вариантов, для которых выполняются условия

H_s? L₀, из дальнейшего рассмотрения исключаются.

Из оставшихся ветвей выбирается ветвь с максимальным значением H_s, и процесс поиска оптимального варианта продолжается. Если в процессе решения будет найдено L_s = ?c_j•x_j > L₀, то полученное решение принимается

x_jЄS

в качестве нового приближенного результата. Вычислительная процедура заканчивается, если для оставшихся ветвей выполняется условие H_s? L₀.

2.2 Практическая часть

Ручной счёт

Данные для расчета:

С?16

Таблица 4

Для удобства расчетов проранжируем таблицу 4 в порядке убывания h_j и присвоим новые номера элементам, следующим образом:

Таблица 5

Для формирования таблицы 6 произведем расчеты:

1)

8

?с_j=19>b₁ - ? c_j•x_j=16-0=16;

j=1 x_jЄS

7

?с_j=1616;

j=1

7

? с = с - ? с_j•x_j - ?с_j=16-0-16=0

x_jЄS j=1

H_s(x₁) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

8

?с_j=18>b₁ - ? c_j•x_j=16-0=16;

j=2 x_jЄS

7

?с_j=1516;

j=2

7

? с = с - ? с_j•x_j - ?с_j=16-0-15=1

x_jЄS j=2

H_s(x₁) = q2+q3+q4+q5+q6+q7+h₈? с = 0.5767

2)

8

?с_j=18>b₁ - ? c_j•x_j=16-1=15;

j=2 x_jЄS

7

?с_j=1515;

j=2

7

? с = с - ? с_j•x_j - ?с_j=16-1-15=0

xjЄS j=2

H_s(x₂) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

8

?с_j=16>b₁ - ? c_j•x_j=16-1=15;

j=3 xjЄS

7

?с_j=1315;

j=3

7

? с = с - ? с_j•x_j - ?с_j=16-1-13=2

xjЄS j=3

H_s(x₂) = q1+q3+q4+q5+q6+q7+h₈? с = 0.5734

3)

8

?с_j=16>b₁ - ? c_j•x_j=16-1-2=13;

j=3 xjЄS

7

?с_j=1313;

j=3

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-13=0

xjЄS j=3

H_s(x₃) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

8

?с_j=15>b₁ - ? c_j•x_j=16-1-2=13;

j=4 xjЄS

7

?с_j=1213;

j=4

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-12=1

xjЄS j=4

H_s(x₃) = q1+q2+q4+q5+q6+q7+h₈? с = 0.5967

4)

8

?с_j=15>b₁ - ? c_j•x_j=16-1-2-1=12;

j=4 xjЄS

7

?с_j=1212;

j=4

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-12=0

xjЄS j=4

H_s(x₄) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

9

?с_j=17>b₁ - ? c_j•x_j=16-1-2-1=12;

j=5 xjЄS

8

?с_j=1112;

j=5

8

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-11=1

xjЄS j=5

H_s(x₄) = q1+q2+q3+q5+q6+q7+q8+h₉? с = 0.56167

5)

8

?с_j=11>b₁ - ? c_j•x_j=16-1-2-1-4=8;

j=5 xjЄS

7

?с_j=88;

j=5

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-4-8=0

xjЄS j=5

H_s(x₅) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

8

?с_j=9>b₁ - ? c_j•x_j=16-1-2-1-4=8;

j=6 xjЄS

7

?с_j=68;

j=6

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-4-6=2

xjЄS j=6

H_s(x₅) = q1+q2+q3+q4+q6+q7+h₈? с = 0.5934

6)

8

?с_j=9>b₁ - ? c_j•x_j=16-1-2-1-4-2=6;

j=6 xjЄS

7

?с_j=66;

j=6

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-4-2-6=0

xjЄS j=6

H_s(x₆) = q1+q2+q3+q4+q5+q6+q7+h₈? с = 0.61

9

?с_j=10>b₁ - ? c_j•x_j=16-1-2-1-4-2=6;

j=7 xjЄS

8

?с_j=46;

j=7

7

? с = с - ? с_j•x_j - ?с_j=16-1-2-1-4-2-4=2

xjЄS j=6

H_s(x₆) = q1+q2+q3+q4+q5+q7+q8+h₉? с = 0.56334

7) C_s = ? c_j•x_j, проверяем условие С-С_s. Если с (x_j)> С-С_s, то эти

xjЄS

элементы вводятся в множество E_s.

E_s = {x8,x9,x10}

с₇=1>b₁ - ? c_j•x_j=16-1-2-1-4-2-5=1;

xjЄS

с₇=11;

Условие не выполняется.

8) L_s = q1+q2+q3+q4+q5+q6+q7 = 0.61

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

Таблица 6

Для контроля системы необходимо использовать следующие параметры контроля: x1,x2,x3,x4,x5,x6,x7.

Перейдем на старую нумерацию элементов и построим дерево возможных вариантов решения. Оно представлено на рис.1. Около каждой вершины указана верхняя граница решения. Так как все эти оценки не превышают величины 0,61, то, следовательно, полученное решение L₀ = 0,61 является ...........