2 Графичен метод

2.1 Стъпки на графичния метод

2.1.1 Построява се множеството 𝑋

Множеството X в контекста на линейното програмиране обикновено се дефинира чрез линейни неравенства от вида: ax_1 + bx_2 \leq c

За да представим графично това множество, първо трябва да начертаем неговата граница – правата линия, дефинирана от уравнението: ax_1 + bx_2 = c

Общ алгоритъм за начертаване Най-бързият начин да начертаем права в координатна система (x_1, x_2) е чрез намиране на пресечните ѝ точки с осите (т.нар. “отрези”):

Пресечна точка с оста x_2 (когато x_1 = 0): x_1 = 0 \implies bx_2 = c \implies x_2 = \frac{c}{b}
Пресечна точка с оста x_1 (когато x_2 = 0): x_2 = 0 \implies ax_1 = c \implies x_1 = \frac{c}{a}

Нека е дадена правата с уравнение 2x_1 + 3x_2 = 6.

Намиране на точка върху ординатата (x_2): Ако приемем x_1 = 0, уравнението става: 3x_2 = 6 \implies x_2 = \frac{6}{3} = 2 Получаваме точка с координати (0, 2).
Намиране на точка върху абсцисата (x_1): Ако приемем x_2 = 0, уравнението става: 2x_1 = 6 \implies x_1 = \frac{6}{2} = 3 Получаваме точка с координати (3, 0).
Графично представяне: Двете точки се отбелязват в координатната система и през тях се прекарва правата линия.

Точка	x_1 (Абсциса)	x_2 (Ордината)
A	0	2
B	3	0

Определяне на полуравнината (Множество X) Тъй като целта е да определим множеството X за неравенството 2x_1 + 3x_2 \leq 6:

Линията: Представлява границата, където стойността е точно 6.
Областта: За да разберем коя страна на правата да защриховаме, тестваме с точката (0,0): 2(0) + 3(0) = 0 \leq 6
Извод: Тъй като 0 \leq 6 е вярно, търсеното множество е полуравнината, която съдържа началото на координатната система.

2.1.2 Построява се нормалния вектор на целевата функция

Какво е нормален вектор?

В геометрията, нормала към даден обект (права или равнина) е обект, който е перпендикулярен на него.

В двумерното пространство (R^2), нормален вектор на дадена права наричаме всеки вектор, който е ортогонален (сключва ъгъл 90°) с нея. Обикновено за удобство го чертаем с начало в координатното начало (0,0).

Математическо определение:

Ако имаме права p, дефинирана с общото уравнение: p: ax_1 + bx_2 - c = 0

Тогава нормалният вектор \vec{n} се състои директно от коефициентите пред x_1 и x_2: \vec{n} = \begin{pmatrix} n_1 \\ n_2 \end{pmatrix} = \begin{pmatrix} a \\ b \end{pmatrix}

Пример: За правата p: 2x_1 - 3x_2 - 6 = 0, нормалният вектор е: \vec{n} = (2, -3)

Построяване на нормалния вектор

Когато разглеждаме целева функция или уравнение на права от вида z(x) = ax_1 + bx_2

Коефициентите (a, b) определят посоката на вектора.
Векторът \vec{n} = \begin{pmatrix} a \\ b \end{pmatrix} е винаги перпендикулярен на линията ax_1 + bx_2 = c.
Промяната на константата c не променя наклона на правата, а само я “плъзга” успоредно в равнината.

Взаимно положение на две прави

Защо правите ax_1 + bx_2 = c_1 и ax_1 + bx_2 = c_2 са успоредни? Защото имат един и същ нормален вектор, което означава, че имат еднакъв наклон.

Нека са дадени две прави в общ вид: g_1: A_1X + B_1Y + C_1 = 0 g_2: A_2X + B_2Y + C_2 = 0

Съществуват три основни случая за тяхното взаимно положение:

Случай	Условие (Коефициенти)	Описание
1. Съвпадат	\frac{A_1}{A_2} = \frac{B_1}{B_2} = \frac{C_1}{C_2}	Правите са идентични.
2. Успоредни	\frac{A_1}{A_2} = \frac{B_1}{B_2} \neq \frac{C_1}{C_2}	Имат еднакъв наклон, но различни пресечни точки.
3. Пресичат се	\frac{A_1}{A_2} \neq \frac{B_1}{B_2}	Имат различна посока и се пресичат в една точка.

Обобщение

Нормалният вектор е най-лесният начин да определим наклона на една права в аналитичната геометрия.
Коефициентите пред променливите в уравнението директно ни дават координатите на този вектор.
Ако две прави имат пропорционални коефициенти пред x_1 и x_2, те са или успоредни, или съвпадат.

2.1.3 Намиране на оптимален план

Свойство 1: Оптималните планове на множеството 𝑋 не могат да бъдат вътрешни точки за множеството. Ако съществуват такива планове, то те са контурни за множеството 𝑋. Контурът на множеството 𝑋 е образуван от части от правите (отсечки), които съответстват на ограниченията. Точката, в която се пресичат две такива прави се нарича опорен план (връх на множеството 𝑋).

Свойство 2: Двумерната задача може да няма решение (оптимален план), може да има един оптимален план или безброй много, които се изобразяват като отсечка или лъч в множеството 𝑋.

2.2 Задача за планиране на производство

Бутиково кафене предлага два вида еспресо: Супер еспресо и Делукс еспресо. За приготвянето на един килограм от първия вид еспресо са необходими по равни части арабика и робуста, а рецептата за Делукс предвижда смес от арабика и робуста в пропорция 1 към 3. Доставчиците са готови да осигурят 120 кг арабика и 160 кг. робуста. Заведението знае, че няма да може да продаде повече от 150 кг. Делукс еспресо. От всеки продаден килограм Супер еспресо заведението печели 40 EUR, докато печалбата от килограм Делукс възлиза на 50 EUR.

Колко от двата типа еспресо ще препоръчате на кафенето да произведе, за да максимизира печалбата си?

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

Целеви променливи:

\begin{align*} & x_1: \text{ Супер еспресо (кг.)}\\ & x_2: \text{ Делукс еспресо (кг.)} \end{align*}

\max z = 40 x_1 + 50 x_2 \text{ (целева функция)}

\begin{align*} 0.5 x_1 + 0.25 x_2 & \leq 120 \text{ (арабика)}\\ 0.5 x_1 + 0.75 x_2 & \leq 160 \text{ (робуста)} \\ 0 \cdot x_1 + x_2 & \leq 150 \text{ (търсене Делукс)}\\ x_1 & \geq 0 \\ x_2 & \geq 0 \end{align*}

Както и в предишната задача ще изобразим графично допустимото множество, като начертаем правите, към всяко от петте неравенства:

\begin{align} 0.5 x_1 & + 0.25 x_2 & = & 120 \\ 0.5 x_1 & + 0.75 x_2 & = & 160 \\ 0 \cdot x_1 & + x_2 & = & 150 \end{align}

Първо ще пресметнем пресечните точки на трите прави с осите x_1 и x_2?

Права арабика: (0, 120 / 0.25 = 480), (120 / 0.5 = 240, 0)
Права робуста: (0, 160 / 0.75), (160 / 0.5, 0)
Права търсене на Делукс: (0, 150), (100, 150) Тази права е успоредна на оста x_1.

Покажи

from IPython.display import display
import matplotlib.pyplot as plt
import numpy as np
from plotly import graph_objects as go

fig, ax = plt.subplots()

ax.plot([0, 240], [480, 0], label=r"$0.5 x_1 + 0.25 x_2 = 120$ (арабика)")
ax.plot([0, 320], [213.33, 0], label=r"$0.5 x_1 + 0.75 x_2 = 160$ (робуста)")
ax.plot([0, 400], [150, 150], label=r"$x_2 = 150$ (търсене на Делукс)")

ax.set_xlim((0, 400))
ax.set_ylim((0, 500))
ax.set_xlabel(r'$x_1$ (Супер, кг.)')
ax.set_ylabel(r'$x_2$ (Делукс, кг.)')
ax.legend(loc=0)

Покажи

vertices = [
    [0, 0],
    [0, 150],
    [95, 150],
    [200, 80],    
    [240, 0]
]

vertices_x = [v[0] for v in vertices]
vertices_y = [v[1] for v in vertices]

for v in vertices:
    ax.annotate(
        f"({v[0]}, {v[1]})",
        (v[0], v[1]),
        textcoords="offset points",
        xytext=(0, 5)
    )

ax.fill(vertices_x, vertices_y, color='grey', alpha=0.3)

display(fig)

Целева функция и оптимален план

За да определим оптималния план графично ще начертаем прави, съответстващи на различни нива на печалба.

Целевата функция зависи от две променливи, което означава, че за да я изобразим графично ще трябва да използваме трето измерение.

Покажи

# Create a 3d scatter plot for the vertices of the feasible region with plotly

fig = go.Figure(data=[go.Scatter3d(
    x=vertices_x,
    y=vertices_y,
    z=[0] * len(vertices_x),
    mode='markers+text',
    marker=dict(size=2, color='red'),
    text=[f"({v[0]}, {v[1]}, 0)" for v in vertices],
    textposition="top center"
)])

# Draw lines connecting the vertices

for i in range(len(vertices)):
    x_line = [vertices[i][0], vertices[(i + 1) % len(vertices)][0]]
    y_line = [vertices[i][1], vertices[(i + 1) % len(vertices)][1]]
    z_line = [0, 0]
    fig.add_trace(go.Scatter3d(
        x=x_line,
        y=y_line,
        z=z_line,
        mode='lines',
        line=dict(color='blue', width=2)
    ))

fig.update_layout(
    title="Целевата функция на проблема с кафето",
    scene=dict(
        xaxis_title=r"Super (kg)",
        yaxis_title=r"Deluxe (kg)",
        zaxis_title=r"Печалба (Euro)",
        zaxis=dict(range=[0, 12000]),
        aspectmode='cube'
    ),
    showlegend=False
)

# Draw the plane of the objective function

X1 = np.linspace(0, 250, 20)
X2 = np.linspace(0, 250, 20)

X1, X2 = np.meshgrid(X1, X2)
Z = 40 * X1 + 50 * X2

fig.add_trace(go.Surface(
    x=X1,
    y=X2,
    z=Z,
    opacity=0.5,
    colorscale='Viridis'
))

# Draw a plane at level z = 6000

z_plane_1 = 6000

fig.add_trace(go.Surface(
    x=X1,
    y=X2,
    z=np.full(X1.shape, z_plane_1),
    opacity=0.2,
    colorscale='Viridis',
    showscale=False
))

# Draw a line at the intersection of the objective plane and the z = 10000 plane
intersection_x_1 = np.linspace(0, 250, 20)
intersection_y_1 = (z_plane_1 - 40 * intersection_x_1) / 50

fig.add_trace(go.Scatter3d(
    x=intersection_x_1,
    y=intersection_y_1,
    z=[z_plane_1] * len(intersection_x_1),
    mode='lines',
    opacity=0.4,
    line=dict(color='red', width=3),
    name='Intersection'
))

# Draw the intersection line on the xy plane

fig.add_trace(go.Scatter3d(
    x=intersection_x_1,
    y=intersection_y_1,
    z=[0] * len(intersection_x_1),
    mode='lines',
    opacity=0.4,
    line=dict(color='red', width=2),
    name='Intersection'
))
fig.update_layout(
    title="Целева функция",
    scene=dict(
        xaxis_title=r"Super (kg)",
        yaxis_title=r"Deluxe (kg)",
        zaxis_title=r"Gewinn (Euro)",
        zaxis=dict(range=[0, 13000]),
        yaxis=dict(range=[0, 250]),
        xaxis=dict(range=[0, 250]),
        aspectmode='cube'
    ),
    showlegend=False
)

# Draw a plane at the optimal solution z = 12000

z_plane_2 = 12000

fig.add_trace(go.Surface(
    x=X1,
    y=X2,
    z=np.full(X1.shape, z_plane_2),
    opacity=0.2,
    colorscale='Viridis',
    showscale=False
))

intersection_x = np.linspace(0, 250, 20)
intersection_y = (z_plane_2 - 40 * intersection_x) / 50

# Draw the intersection line on the xy plane

fig.add_trace(go.Scatter3d(
    x=intersection_x,
    y=intersection_y,
    z=[z_plane_2] * len(intersection_x),
    mode='lines',
    opacity=0.4,
    line=dict(color='red', width=3),
    name='Intersection'
))
fig.add_trace(go.Scatter3d(
    x=intersection_x,
    y=intersection_y,
    z=[0] * len(intersection_x),
    mode='lines',
    opacity=0.4,
    line=dict(color='red', width=3),
    name='Intersection'
))

fig.show()

Фигура 2.3: Die Gewinnfunktion, die zulässige Menge und Isolinien des Gewinns bei 6000 und 12000 Euro.

Нека да начертаем тези прави на графиката на допустимото множество. Фигура 2.4 показва графиката на допустимото множество и правите на печалбата от 3000 лв., 50000 лв. и 12000 лв.

Покажи

fig, ax = plt.subplots()

ax.fill(vertices_x, vertices_y, color='grey', alpha=0.3)
for v in vertices:
    ax.annotate(
        f"({v[0]}, {v[1]})",
        (v[0], v[1]),
        textcoords="offset points",
        xytext=(0, 5)
    )

ax.plot([0, 3000 / 40], [3000 / 50, 0], label=r'$40 x_1 + 50 x_2 = 3000$')
ax.plot([0, 5000 / 40], [5000 / 50, 0], label=r'$40 x_1 + 50 x_2 = 5000$')
ax.plot([0, 12000 / 40], [12000 / 50, 0], label=r'$40 x_1 + 50 x_2 = 12000$')

ax.quiver(0, 0, 40, 50, angles='xy', scale_units='xy', scale=1, color='grey', label='Нормален вектор')

ax.set_xlim((0, 260))
ax.set_ylim((0, 260))

ax.legend(loc=0)

2.3 Задачa за максимизиране на печалба

Управителят на компания иска да определи какво количество от двата продукта на компанията да произведе за определен период, така че да осигури максимална печалба.

В таблицата са посочени общия брой работни часове, с които разполага всеки отдел, и часовете, нужни за производството на всеки продукт:

Отдел	Продукт 1 (раб. часове за 1 бр.)	Продукт 2 (раб. часове за 1 бр.)	Общ брой налични работни часове
A	1.00	0.35	100
Б	0.30	0.20	36
В	0.20	0.50	50
Печалба за единица	30.00 euro	15.00 euro	-

Съставяне на математическия модел

Първата стъпка е да дефинираме променливите, целта и ограниченията на базата на таблицата.

Неизвестни:
- x_1 - количество от Продукт 1
- x_2 - количество от Продукт 2
Целева функция (Максимизиране на печалбата): z(x) = 30x_1 + 15x_2 \rightarrow \max
Ограничения (спрямо наличните часове):
- 1. Отдел А: x_1 + 0.35x_2 \leq 100
- 1. Отдел Б: 0.3x_1 + 0.2x_2 \leq 36
- 1. Отдел В: 0.2x_1 + 0.5x_2 \leq 50
- 1. Неотрицателност на неизвестните: x_1 \geq 0, x_2 \geq 0

Опростяване на модела и намиране на координати

За по-лесно пресмятане, умножаваме ограниченията за отделите по 10, за да премахнем десетичните дроби. Получаваме следната система:

z(x) = 30x_1 + 15x_2 \rightarrow \max

Ограничения: 10x_1 + 3.5x_2 \leq 1000 3x_1 + 2x_2 \leq 360 2x_1 + 5x_2 \leq 500 x_1 \geq 0, x_2 \geq 0

Намираме пресечните точки (отрезите) на правите с координатните оси, за да дефинираме границите на множеството:

Права p_1: (100, 0) ; (0, 285.71)
Права p_2: (120, 0) ; (0, 180)
Права p_3: (250, 0) ; (0, 100)

Намиране на оптималния план

Точката, в която се намира оптималният план, е пресечната точка на правите p_1 и p_2. Нека я наречем точка A.

т. A = p_1 \cap p_2

За да намерим точните координати на т. A, решаваме системата уравнения, съставена от двете прави: 10x_1 + 3.5x_2 = 1000 3x_1 + 2x_2 = 360

Решението на системата ни дава координатите на оптималния план: x_1 = 77.89, \quad x_2 = 63.157 \implies \text{т. } A = (77.89, 63.157)

Намираме стойността на целевата функция в т. A: Заместваме получените координати в целевата функция z(x): z(A) = 30x_1 + 15x_2 = 30 \cdot 77.89 + 15 \cdot 63.157 = 3284.1

https://www.desmos.com/calculator/pvm4jqfzso

Краен резултат: \max_x z(x) = z(77.89, 63.157) = 3284.1 \text{ euro}

2.4 Задача, решена с GUROBI модел

Предприятие произвежда два вида боя: екстериорна (за външни повърхности) и интериорна (за вътрешни помещения). Производството двете бои изисква три суровини: разтворител, багрило и смола. За производство на един тон екстериорна боя са нужни 6 тона разтворител, един тон багрило и 2 тона смола. За производство на един тон интериорна боя са нужни 4 тона разтворител, 2 тона багрило и 7 тона смола. На разположение са следните количества суровини: 24 тона разтворител, 6 тона багрило и 4 тона смола. Колко тона от всяка боя да произведе предприятието, за да постигне максимална печалба? От продажбата на един тон екстериорна боя предприятието печели 5 000 лв., а от продажбата на един тон интериорна боя - 4 000 лв.

Таблица 2.1: Консумация на материали за производството на интериорна и екстериорна боя и печалба за продаден тон

	Екстериорна боя	Интериорна боя	Дневна наличност на ресурс
Разтворител	6	4	24
Багрило	1	2	6
Смола	2	7	4
Печалба (1 000 лв./тон)	5	4

Съставете математически модел, който да намери оптималния план за производство на боите.

Математически модел

Целеви променливи: \begin{align*} & x_1: \text{ Екстериорна боя (тонове)}\\ & x_2: \text{ Интериорна боя (тонове)} \end{align*}

\max z = 5 x_1 + 4 x_2 \text{ (целева функция)}

Ограничения:

\begin{align*} 6x_1 + 4x_2 & \leq 24 \text{ (разтворител)}\\ x_1 + 2x_2 & \leq 6 \text{ (багрило)}\\ 2x_1 + 7x_2 & \leq 4 \text{ (смола)}\\ x_1 & \geq 0 \\ x_2 & \geq 0 \end{align*}

Покажи

import gurobipy as gp
from gurobipy import GRB

# Create a new model
m = gp.Model("Paint")
m.Params.LogToConsole = 0

# Create variables
p_exterior = m.addVar(lb=0, name="exterior")
p_interior = m.addVar(lb=0, name="interior")

m.setObjective(5 * p_exterior + 4 * p_interior, GRB.MAXIMIZE)

m.addConstr(6 * p_exterior + 4 * p_interior <= 24, "Solvent")
m.addConstr(1 * p_exterior + 2 * p_interior <= 6, "Dye")
m.addConstr(2 * p_exterior + 7 * p_interior <= 4, "Resin")

m.optimize()

print(f"Optimal solution: p_exterior = {p_exterior.x}, p_interior = {p_interior.x}")

Restricted license - for non-production use only - expires 2027-11-29
Set parameter LogToConsole to value 0
Optimal solution: p_exterior = 2.0, p_interior = 0.0