Lecture 4: Non-linear Programing, Constrained Maximization¶

Items for reiteration:¶

  1. The lecture notes are posted on my website. I've added link to the location at the top of the canvas page.
  2. The lectures are not recorded.
  3. Please post your questions about course content on the canvas discussions for each lecture (they go live after the lecture).
  4. Thanks to those of me that have reached out with questions! I enjoy answering them, and it helps me a lot to see what is working and what is not.
  5. The review lectures will address everyone's favorite question: "Dr. Morris, please tell me exactly what will be on the exam."

A note on incremental cost and marginal cost:¶

  • Marginal Cost: as it's used in economics (i.e. MR=MC under competition) is a mathmatical term: the derivative of the cost function w.r.t. quantity.
  • Incremental Cost is a practical term: the change in cost for a one unit change in output.

The two are often interchangable, but we should know that they are different. They are the same when the cost curve is linear, they are intechangable whenever we decide that a linear approximation is suitable for our purposes.

Let's go through the homework:¶

You manage a two-product firm. The production technology requires a mixture of capital and labor to produce each product. Capital is shared, while labor is specific to each of the two products. For the first product, capital ($K \geq 0$) and labor ($L_1 \geq 0$) must satisfy the following, in order to produce $q_1$ units: $q_{1} \leq \sqrt{KL_{1}}$. Likewise, producing $q_2$ units of the second product requires capital ($K$) and labor (now $L_2 \geq 0$) such that: $q_2 \leq\sqrt{KL_2}$. In addition, total capital is limited to a maximum of 200 units. (So $K \leq 200$.) Naturally, $K$, $L_1$, $L_2$ are all required to be non-negative. Capital costs 100 per unit, labor for the first product costs 140 per unit and labor for the second product costs 175 per unit. The first product sells for 275 per unit, and the second sells for 300 per unit.

  1. Initially suppose only the first product is present. Determine and interpret your optimal production plan.
  2. Next suppose only the second product is present. Determine and interpret your optimal production plan.
  3. Now assume both products are present. Determine and interpret your optimal production plan.
  4. Repeat the three parts above assuming the first product sells for 200 per unit.
  1. Fill in the following table:
$q_1$ $q_2$ MC of $q_1$ with $q_2$ fixed. MC of $q_2$ with $q_1$ fixed.
50 50
50 100
50 150
100 50
100 100
100 150
150 100
150 150

(MC is the marginal cost)

1. Prod 1: Determine and interpret your optimal production plan.¶

What are $q_{1} \leq \sqrt{KL_{1}}$ and $q_2 \leq\sqrt{KL_2}$?¶

  • A name for this sort of function is a "production function"
  • The expression $q_{1} \leq \sqrt{KL_{1}}$ is that rate at which our production process changes inputs into outputs.
  • Notice that this even though there is a $\leq$ sign, that this is not a constraint per se because we get to choose $K$ and $L_1$, and it is only that choice that constrains the amount of $q_1$ that we produce.
  • Also notice that if we are trying to optimize we will never produce less than $q_1=\sqrt{KL_{1}}$

Let's look at some plots to see the implications of this. Note that I'll use python for the examples today, this is because it makes plotting easy, and I want you to think of these as list of steps.

In [2]:
# quick set up for our plot
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

# create the framework for our plot
l1 = np.linspace(0, 300, 300)
k = np.linspace(0, 200, 300)
L1, K = np.meshgrid(l1, k)

# CALCULATE THE PRODUCTION FUNCTION 
Q1 = np.sqrt(L1*K)
In [3]:
# Lets look at the objects we created:
l1[0:10] # first 10 elements of l1
# k[0:10] 
# K
# L1
Q1
Out[3]:
array([[  0.        ,   0.        ,   0.        , ...,   0.        ,
          0.        ,   0.        ],
       [  0.        ,   0.81922734,   1.15856241, ...,  14.11830831,
         14.14205653,  14.16576494],
       [  0.        ,   1.15856241,   1.63845468, ...,  19.96630309,
         19.99988814,  20.0334169 ],
       ...,
       [  0.        ,  14.11830831,  19.96630309, ..., 243.3105196 ,
        243.71978906, 244.12837239],
       [  0.        ,  14.14205653,  19.99988814, ..., 243.71978906,
        244.12974694, 244.53901755],
       [  0.        ,  14.16576494,  20.0334169 , ..., 244.12837239,
        244.53901755, 244.94897428]])
In [4]:
# Create the figure and add a 3D axis
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
# Plot the data
ax.plot_surface(L1, K, Q1)
# Set axis labels and show the plot
ax.set_xlabel('L1')
ax.set_ylabel('K')
ax.set_zlabel('Q1')
plt.show()

What if we hold K fixed? What will the graph look like?¶

In Q1 = np.sqrt(L1*K) we will plug in values for K and let L1 vary.

In [5]:
l1 = np.linspace(0, 300, 300)
Q1_10 = np.sqrt(L1*10)    # k = 10
Q1_100 = np.sqrt(L1*100)  # k = 100
Q1_200 = np.sqrt(L1*200)  # k = 200
In [6]:
# Lets look inside Q1 and Q1_10
Q1[0:5]
Q1_10[0:5]
Out[6]:
array([[ 0.        ,  3.16756134,  4.4796082 , ..., 54.58876359,
        54.68058664, 54.77225575],
       [ 0.        ,  3.16756134,  4.4796082 , ..., 54.58876359,
        54.68058664, 54.77225575],
       [ 0.        ,  3.16756134,  4.4796082 , ..., 54.58876359,
        54.68058664, 54.77225575],
       [ 0.        ,  3.16756134,  4.4796082 , ..., 54.58876359,
        54.68058664, 54.77225575],
       [ 0.        ,  3.16756134,  4.4796082 , ..., 54.58876359,
        54.68058664, 54.77225575]])

What if we hold K fixed? What will the graph look like?¶

In [7]:
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(L1, K, Q1_10)
ax.plot_surface(L1, K, Q1_100)
ax.plot_surface(L1, K, Q1_200)
ax.set_xlabel('L1')
ax.set_ylabel('K')
ax.set_zlabel('Q1')
plt.show()

Optimal production plan for $q_1$ alone:¶

Let's refer to the prices as "P", and the labor and capital costs as "C" $$ \Pi(P,C) \equiv 275 \times q_1 - 140 \times L_1 - 100 \times K $$ $$ q_1 \leq \sqrt{L_1 \times K} $$ $$ K \leq \bar{K} = 200 $$

This makes it seem like we have more choices than we actually have. The constraint (production function) makes it so that we can only choose two of K, q1, and l1.

$$ \Pi(P,C) \equiv 275 \times q_1 - 140 \times \frac{q_1^2}{K} - 100 \times K $$$$ K \leq \bar{K} = 200 $$

So what do we do with this?¶

$$ \Pi(P,C) \equiv 275 \times q_1 - 140 \times \frac{q_1^2}{K} - 100 \times K $$$$ K \leq \bar{K} = 200 $$

How do we choose the optimal production plan?

There are a number of approaches, I'm going to show you a very general python solver. This will be useful for those of you that are using python, I want the rest of you to focus on the steps we are taking so you can apply them to whatever solver you like.

iPRS: How did you solve this?

Since this is not a linear problem let's use GEKKO¶

In [8]:
# set up the solver
from gekko import GEKKO
m = GEKKO(remote=False) # here we create a gekko object called 'm'

Step 1: Choice Variables and their constraints¶

In [9]:
# Initialize the decision variables
q1 = m.Var(
    name="q1", 
    lb=0 # the lower bound
)
q1.value=1 # most solvers run faster when you give a starting point

k =  m.Var(
    name="k", 
    lb=0, # lower bound
    ub=200 # upper bound
) 
k.value=200 # this is our first guess to speed up the solution

Step 2: Write down the objective function (what we want to maximize), any remaining constraints, and solve.¶

$$ \Pi(P,C) \equiv 275 \times q_1 - 140 \times \frac{q_1^2}{K} - 100 \times K $$

What happened to $ K \leq \bar{K} = 200 $?

In [10]:
m.Maximize(
    275*q1 - 140*((q1**2)/k) - 100*k
)
m.solve(disp=False) # silencing the output because it is diagnostic

The output from the model is just the choice variables: we need to calculate profit and labor:

In [11]:
q1, k
Out[11]:
([196.42857143], [200.0])

We need to calculate profit and labor:

In [16]:
profit = (275*q1.value[0] - 140*((q1.value[0]**2)/k.value[0]) - 100*k.value[0])
l1 = q1.value[0]**2/k.value[0]

This gives us the following solutions¶

(note that I'm rounding here)

In [13]:
print('q1     ', int(q1.value[0])) # q1 is a list q1[0] is the first element of the list.
print('l1     ',int(l1))
print('K      ', int(k.value[0]))
print('profit ', int(profit))

q1      196
l1      192
K       200
profit  7008

Optimal production plan for $q_2$ alone:¶

Let's refer to the prices as "P", and the labor and capital costs as "C" $$ \Pi(P,C) \equiv 300 \times q_2 - 175 \times L_2 - 100 \times K $$ $$ q_2 \leq \sqrt{L_2 \times K} $$ $$ K \leq \bar{K} = 200 $$

We make a similar substitution here:

$$ \Pi(P,C) \equiv 300 \times q_2 - 175 \times \frac{q_2^2}{K} - 100 \times K $$$$ K \leq \bar{K} = 200 $$
In [14]:
# same set up repeated to clean out 'm'
m = GEKKO(remote=False)
q2 = m.Var(name="q2", lb=0)
q2.value=1 
k =  m.Var(name="k", lb=0, ub=200)
k.value=200

Write down the objective function again and solve:¶

$$ \Pi(P,C) \equiv + 300 \times q_2 - 175 \times \frac{q_2^2}{K} - 100 \times K $$
In [15]:
m.Maximize(
    300*q2 - 175*((q2**2)/k) - 100*k
)
m.solve(disp=False)
profit = (300*q2.value[0] - 175*((q2.value[0]**2)/k.value[0]) - 100*k.value[0])
l2 = q2.value[0]**2/k.value[0]

This gives us the following solutions¶

(note that I'm rounding here)

In [17]:
print('q2     ', int(q2.value[0]))
print('l2     ',int(l2))
print('K      ', int(k.value[0]))
print('profit ', int(profit))

q2      171
l2      146
K       200
profit  7008

Optimal production plan for both:¶

Let's refer to the prices as "P", and the labor and capital costs as "C" $$ \Pi(P,C) \equiv 275 \times q_1 + 300 \times q_2 - 140 \times L_1 - 175 \times L_2 - 100 \times K $$ $$ q_1 \leq \sqrt{L_1 \times K} $$ $$ q_2 \leq \sqrt{L_2 \times K} $$ $$ K \leq \bar{K} = 200 $$

Again both constraints allow us to eliminate $L_1$: $$ \Pi(P,C) \equiv 275 \times q_1 + 300 \times q_2 - 140 \times \frac{q_1^2}{K} - 175 \times \frac{q_2^2}{K} - 100 \times K $$ $$ K \leq \bar{K} = 200 $$

In [18]:
# reset gekko
m = GEKKO(remote=False)
# Initialize the decision variables
q1 = m.Var(name="q1", lb=0)
q1.value=1
q2 = m.Var(name="q2", lb=0)
q2.value=1
k =  m.Var(name="k", lb=0, ub=200) 
k.value=200
In [19]:
# write the objective function  and solve
m.Maximize(
    275*q1 - 140*((q1**2)/k)
    + 300*q2 - 175*((q2**2)/k)
    - 100*k
)
m.solve(disp=False)
In [20]:
# quick calcs to clean things up
profit = (
    275*q1.value[0] - 140*((q1.value[0]**2)/k.value[0]) 
    + 300*q2.value[0] - 175*((q2.value[0]**2)/k.value[0])  
    - 100*k.value[0]
)
l1 = q1.value[0]**2/k.value[0]
l2 = q2.value[0]**2/k.value[0]

This gives us the following solutions¶

(note that I'm rounding here)

In [21]:
print('q1     ', int(q1.value[0]))
print('l1     ',int(l1))
print('q2     ', int(q2.value[0]))
print('l2     ',int(l2))
print('K      ', int(k.value[0]))
print('profit ', int(profit))

q1      196
l1      192
q2      171
l2      146
K       200
profit  32723

We always pick k=200. What does this mean?¶

  • This means that the capital constraint is preventing us from making the next unit of output.
  • When a constraint does this we often call it a 'binding constraint'.

The q1,L1 and q2,L2 pairs that we pick are the same when we consider the products separately and together. What does this mean?¶

  • It means that the profit function is separable when the capital constraint is binding.

Repeat all three steps with $P_1=200$¶

In [22]:
m = GEKKO(remote=False)
q1 = m.Var(name="q1", lb=0)
q1.value=1
q2 = m.Var(name="q2", lb=0)
q2.value=1
k =  m.Var(name="k", lb=0, ub=200)
k.value=200 
m.Maximize(
    200*q1 - 140*((q1**2)/k) - 100*k
)
m.solve(disp=False) # silencing the out put because it is diagnostic
profit = (
    200*q1.value[0] 
    - 140*((q1.value[0]**2)/k.value[0]) 
    - 100*k.value[0]
)
l1 = q1.value[0]**2/k.value[0]

$Q_1$ when $P_1=200$¶

In [23]:
print('q1     ', int(q1.value[0]))
print('l1     ', int(l1))
print('K      ', int(k.value[0]))
print('profit ', int(profit))

q1      0
l1      0
K       0
profit  0

The optimal amount of Q1 to produce is 0.¶

What does this mean?

What does this mean? This means that a firm that only produces Q1 should not produce it at this price, given this cost structure. What does this mean for a firm that produces Q2? What should they do? enter the market for Q1 or not?

Lets look at the profit functions for the two single product firms with k=200¶

In [24]:
# set up
def Pi1(q1):
    return 200*q1 - 140*((q1**2)/200) 
def Pi2(q2):
    return 300*q2 - 175*((q2**2)/200) 
q1=np.linspace(0,300,300)
q2=np.linspace(0,300,300)

Pq1=Pi1(q1)- 100*200
Pq2=Pi2(q2)- 100*200
In [25]:
# plot so we can see what is goingon
plt.plot(q1,Pq1,label='Profit (q1)')
plt.plot(q1,Pq2,label='Profit (q2)')
plt.xlabel('Quantity')
plt.legend()
plt.grid()
plt.axhline(0, color='grey')
plt.title('Profit with k=200 and P1=200')
plt.show()

$Q_1$ and $Q_2$ when $P_1=200$¶

Now we are back to the multiproduct firm

In [26]:
# Let's look at the two together
# reset gekko
m = GEKKO(remote=False)
# Initialize the decision variables
q1 = m.Var(name="q1", lb=0)
q1.value=1
q2 = m.Var(name="q2", lb=0)
q2.value=1
k =  m.Var(name="k", lb=0, ub=200) # this is the constraint
k.value=200 # this is our first guess to speed up the solution
m.Maximize(
    200*q1 - 140*((q1**2)/k)
  + 300*q2 - 175*((q2**2)/k)
  - 100*k
)
m.solve(disp=False)
profit = (
    200*q1.value[0] - 140*((q1.value[0]**2)/k.value[0]) 
  + 300*q2.value[0] - 175*((q2.value[0]**2)/k.value[0])  
  - 100*k.value[0]
)
l1 = q1.value[0]**2/k.value[0]
l2 = q2.value[0]**2/k.value[0]
In [27]:
print('q1     ', int(q1.value[0]))
print('l1     ', int(l1))
print('q2     ', int(q2.value[0]))
print('l2     ', int(l2))
print('K      ', int(k.value[0]))
print('profit ', int(profit))

q1      142
l1      102
q2      171
l2      146
K       200
profit  20000

Lets plot this and see what it's telling us.¶

In [28]:
# Create data for the plot
def bigPi(q1,q2,k=200):
    return 200*q1 - 140*((q1**2)/k)+ 300*q2 - 175*((q2**2)/k) - 100*k
q1=np.linspace(0,300,300)
q2=np.linspace(0,300,300)
Q1, Q2 = np.meshgrid(q1, q2)
Pi12 = bigPi(Q1,Q2)
In [29]:
# Create the figure and add a 3D axis
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
# Plot the data
ax.plot_surface(Q1, Q2, Pi12)
# Set axis labels and show the plot
ax.set_xlabel('Q1')
ax.set_ylabel('Q2')
ax.set_zlabel('Profit')
plt.show()