Simple nonlinear complementarity problem
Contents
Simple nonlinear complementarity problem¶
Randall Romero Aguilar, PhD
This demo is based on the original Matlab demo accompanying the Computational Economics and Finance 2001 textbook by Mario Miranda and Paul Fackler.
Original (Matlab) CompEcon file: demslv13.m
Running this file requires the Python version of CompEcon. This can be installed with pip by running
!pip install compecon --upgrade
Last updated: 2022-Sept-05
About¶
The problem is to solve
(3)¶\[\begin{equation}
f(x, y) = \begin{bmatrix}1+xy -2x^3-x\\ 2x^2-y\end{bmatrix}
\end{equation}\]
subject to \(0 \leq x, y \leq 1\)
from compecon import MCP, jacobian
import numpy as np
x0 = [0.5, 0.5]
Solving the problem without the Jacobian¶
To solve this problem we create a MCP object using a lambda function.
def func(z):
x, y = z
return [1 + x*y - 2*x**3 - x, 2*x**2 - y]
F = MCP(func, [0, 0], [1,1])
Solve for initial guess \(x_0 = [0.5, 0.5]\)
#x = F.broyden(x0, transform='minmax', show=True) # FIXME: generates error
#x = F.broyden(x0, transform='ssmooth', show=True) # FIXME: generates error
#print(f'Solution is {x=}.')
Solving the problem with the Jacobian¶
def func2(z):
x, y = z
f = [1 + x*y - 2*x**3 - x, 2*x**2 - y]
J = [[y-6*x**2-1, x],[4*x, -1]]
return f, J
F2 = MCP(func2, [0, 0], [1,1])
x = F2.zero(x0, transform='minmax', show=True)
print(f'Solution is {x=}.')
Using the MINMAX transformation
Solving nonlinear equations by Newton's method
it bstep change
--------------------
0 1 1.56e-01
1 0 9.84e-03
2 0 3.19e-05
3 0 3.40e-10
Solution is x=array([0.7937, 1. ]).