Optimization in Python

Geometric intro in Python

import matplotlib.pyplot as plt
import numpy as np


# The problem in array form
A = np.array([[2, 4], [3, 2]])
b = np.array([40, 50])
obj = np.array([3, 3])

# The problem in Mathtex
T = '''$Min. 3x + 3y$\n
$s.t.$ (Subject to)\n
$2x + 4y \geq 40$
Greater than 40 once projected to $\\vec {(2,4)}$\n
$3x + 2y \geq 50$
Greater than 50 once projected to $\\vec {(3,2)}$'''

# Intercept, feasiblity cutoff and axes angles
c = (b / A.T).T
f = (A.T * b / (A * A).sum(1)).T
r = np.arctan(np.divide(*A.T)) * 180 / np.pi

# Feasiblity region for both continous and integer cases
X = np.dstack(np.meshgrid(np.linspace(0, 20, 400), np.linspace(0, 25, 500)))
Xi = np.dstack(np.meshgrid(np.arange(20), np.arange(25)))
F = (np.dot(X, A.T) > b).all(2)
Fi = (np.dot(Xi, A.T) > b).all(2)


def plot_constraint0(ax):
    ax.annotate("",
                xy=f[0], xycoords='data',
                xytext=(0, 0), textcoords='data',
                arrowprops=dict(arrowstyle="->", color='r',
                                connectionstyle="arc3"),
                )

def plot_constraint1(ax):
    ax.annotate("",
                xy=f[1], xycoords='data',
                xytext=(0, 0), textcoords='data',
                arrowprops=dict(arrowstyle="->", color='b',
                                connectionstyle="arc3"),
                )

def plot_objective(ax):
    ax.annotate("",
                xy=(12, 12), xycoords='data',
                xytext=(0, 0), textcoords='data',
                arrowprops=dict(arrowstyle="->", color='g',
                                connectionstyle="arc3"),
                )

def plot_boundnlabel(ax):
    for i, color in enumerate('rb'):
        ax.plot([c[i, 0], 0], [0, c[i, 1]], c=color)
        ax.text(*(f[i] * 0.5), s='Constraint %s' % i,
                ha='center', va='center', rotation=90-r[i])

def plot_objectnlabel(ax):
    ax.text(6, 6, s='Object', ha='center', va='center', rotation=45)
    for i in np.linspace(15, 20, 10):
        ax.plot([0, i], [i, 0], 'g--', alpha=0.2)

def plot_lpfeasiblity(ax):
    ax.contourf(X[:, :, 0], X[:, :, 1], 1/F, alpha=0.1)

def plot_ipfeasiblity(ax):
    ax.scatter(*np.where(Fi.T), marker='+')

plt.rcParams['font.family'] = ['Arial']
fig, axes = plt.subplots(1, 5, sharex=True, sharey=True)
axes = axes.ravel()

axes[0].text(-2, 20, T, ha='left', va='top', fontsize=12)
axes[0].axis('off')

plot_constraint0(axes[1])
plot_constraint1(axes[1])
plot_boundnlabel(axes[1])
axes[1].set_title('''Legit solutions must cross the
cutoff lines once projected on to\ncontraint axes''')

plot_constraint0(axes[2])
plot_constraint1(axes[2])
plot_boundnlabel(axes[2])
plot_lpfeasiblity(axes[2])
axes[2].set_title('''Legit solutions must satisfy
both constraints\n''')

plot_constraint0(axes[3])
plot_constraint1(axes[3])
plot_boundnlabel(axes[3])
plot_objective(axes[3])
plot_objectnlabel(axes[3])
plot_lpfeasiblity(axes[3])
axes[3].set_title('''Need to find the minimus among
legit solutions once projected to\nobject axis''')

plot_constraint0(axes[4])
plot_constraint1(axes[4])
plot_boundnlabel(axes[4])
plot_objective(axes[4])
plot_objectnlabel(axes[4])
plot_lpfeasiblity(axes[4])
plot_ipfeasiblity(axes[4])
axes[4].set_title('''For integer problem\nlegit solutions are further
confined to the integer grid''')

for ax in axes:
    ax.set_aspect('equal')

plt.tight_layout()