Book Example: Economic equilibria#

Description: economic model using complementarity conditions from Chapter 19 AMPL book

Tags: ampl-book, finance, complementarity-problem, cbc, gurobi

Notebook author: Marcos Dominguez Velad <marcos@ampl.com>

Model author: N/A

# Install dependencies
%pip install -q amplpy pandas

# Google Colab & Kaggle integration
from amplpy import AMPL, ampl_notebook

ampl = ampl_notebook(
    modules=["coin", "gurobi"],  # modules to install
    license_uuid="default",  # license to use
)  # instantiate AMPL object and register magics

Economic equilibria#

This model is based on the models showed in the Chapter 19 of the AMPL book, related to complementarity problems.

Classic model#

The following model does not use complementarity conditions.

Sets and parameters:
- Products PROD and demand[i] for each product i.
- Activities ACT and cost[j] for each activity j.
- Amount of product i produced by activity j, io[i,j].
Decision variables Level[j]: levels of production activities (for each activity j).
Objective: minimize total production cost,

\[ \sum \limits_{j \in ACT} cost_j \cdot Level_j; \]

So the AMPL formulation is implemented as:

%%writefile econmin.mod
set PROD;   # products
set ACT;    # activities

param cost {ACT} > 0;      # cost per unit of each activity
param demand {PROD} >= 0;  # units of demand for each product
param io {PROD,ACT} >= 0;  # units of each product from
                           # 1 unit of each activity

var Level {j in ACT} >= 0;

minimize Total_Cost:  sum {j in ACT} cost[j] * Level[j];

subject to Demand {i in PROD}:
      sum {j in ACT} io[i,j] * Level[j] >= demand[i];

Overwriting econmin.mod

import pandas as pd

# ACT cost data
act_cost_data = {
    "ACT": ["P1", "P1a", "P2", "P2a", "P2b", "P3", "P3c", "P4"],
    "cost": [2450, 1290, 1850, 3700, 2150, 2200, 2370, 2170],
}
act_cost_df = pd.DataFrame(act_cost_data).set_index("ACT")

# PROD demand data
prod_demand_data = {
    "PROD": ["AA1", "AC1", "BC1", "BC2", "NA2", "NA3"],
    "demand": [70000, 80000, 90000, 70000, 400000, 800000],
}
prod_demand_df = pd.DataFrame(prod_demand_data).set_index("PROD")

# IO matrix data
io_data = {
    "AA1": [60, 8, 8, 40, 15, 70, 25, 60],
    "AC1": [20, 0, 10, 40, 35, 30, 40, 20],
    "BC1": [10, 20, 15, 35, 15, 15, 30, 15],
    "BC2": [15, 20, 10, 10, 15, 15, 30, 10],
    "NA2": [938, 1180, 945, 278, 1182, 896, 1029, 1397],
    "NA3": [295, 770, 440, 430, 315, 400, 370, 450],
}
io_index = ["P1", "P1a", "P2", "P2a", "P2b", "P3", "P3c", "P4"]
io_df = pd.DataFrame(io_data, index=io_index)

ampl.reset()
ampl.read("econmin.mod")

ampl.set_data(prod_demand_df, "PROD")
ampl.set_data(act_cost_df, "ACT")
ampl.param["io"] = io_df.T

ampl.solve(solver="cbc")
ampl.display("Level, Demand.dual")

cbc 2.10.12: cbc 2.10.12: optimal solution; objective 6808640.553
0 simplex iterations
:      Level   Demand.dual    :=
AA1       .      16.7051
AC1       .       5.44585
BC1       .      57.818
BC2       .       0
NA2       .       0
NA3       .       0
P1       0           .
P1a   1555.3         .
P2       0           .
P2a      0           .
P2b      0           .
P3     147.465       .
P3c   1889.4         .
P4       0           .
;

assert ampl.solve_result == "solved", ampl.solve_result

Complementarity model#

Consider the new variables Price[i] for each product, and solve the problem to find an equilibrium instead of an optimum solution. The equilibrium is subject to two conditions:

First, for each product, total output must meet demand and price must be nonnegative, and in addition there must be a complementarity between these relationships, where production exceeds demand the price must be zero, or equivalently, if the price is positive, the production must equal the demand.

subject to Pri_Compl {i in PROD}:
   Price[i] >= 0 complements
      sum {j in ACT} io[i,j] * Level[j] >= demand[i];

Second, for each activity j, the value of resulting product i is Price[i]*io[i,j], so activity j would produce a total value of

\[\sum \limits_{i \in PROD} Price[i] \cdot io[i,j]\]

When equilibrium happens, the previous value must not exceed the activity’s cost per unit cost[j]. Moreover, there is a complementarity between this relationship and the level of activity j, where cost exceeds total value the activity must be zero, or, where the activity cost is positive then the total value must be equal to the cost. This can be expressed as:

subject to Lev_Compl {j in ACT}:
   Level[j] >= 0 complements
      sum {i in PROD} Price[i] * io[i,j] <= cost[j];

%%writefile econ.mod
set PROD;   # products
set ACT;    # activities

param cost {ACT} > 0;      # cost per unit of each activity
param demand {PROD} >= 0;  # units of demand for each product
param io {PROD,ACT} >= 0;  # units of each product from
                           # 1 unit of each activity
var Price {i in PROD};
var Level {j in ACT};

subject to Pri_Compl {i in PROD}:
   Price[i] >= 0 complements
      sum {j in ACT} io[i,j] * Level[j] >= demand[i];

subject to Lev_Compl {j in ACT}:
   Level[j] >= 0 complements
      sum {i in PROD} Price[i] * io[i,j] <= cost[j];

Overwriting econ.mod

Remark: the model is square, as there are \(n+m\) variables and \(n+m\) constraints (\(n\) number of products and \(m\) number of activities). Some solvers can take advantage of this to use different solving techniques.

Finally, let’s see that the equilibrium model gives also the optimal solution for the first classic formulation.

ampl.reset()
ampl.read("econ.mod")

ampl.set_data(prod_demand_df, "PROD")
ampl.set_data(act_cost_df, "ACT")
ampl.param["io"] = io_df.T

ampl.solve(solver="gurobi")
ampl.display("Level, sum {j in ACT} cost[j] * Level[j]")

Gurobi 12.0.1: Gurobi 12.0.1: optimal solution
0 simplex iterations
Objective = find a feasible point.
Level [*] :=
 P1     0
P1a  1555.3
 P2     0
P2a     0
P2b     0
 P3   147.465
P3c  1889.4
 P4     0
;

sum{j in ACT} cost[j]*Level[j] = 6808640

assert ampl.solve_result == "solved", ampl.solve_result