diff --git a/.ipynb_checkpoints/Untitled-checkpoint.ipynb b/.ipynb_checkpoints/Untitled-checkpoint.ipynb
new file mode 100644
index 0000000..4e4b89a
--- /dev/null
+++ b/.ipynb_checkpoints/Untitled-checkpoint.ipynb
@@ -0,0 +1,649 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "data (generic function with 2 methods)"
+ ]
+ },
+ "execution_count": 1,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function data(a=0.5)\n",
+ " a = a\n",
+ "\t\tna = 3\n",
+ "\t\tnc = 4\n",
+ "\t\tns = 4\n",
+ "\t\tnss = 16\n",
+ "\t\tz1 = ones(4)\n",
+ "\t\tz2 = [0.72, 0.92, 1.12, 1.32]\n",
+ "\t\tz3 = [0.86, 0.96, 1.06, 1.16]\n",
+ " z = [16, 3]\n",
+ " z = [[1, i, j] for i in z2 for j in z3] # All possible payoffs\n",
+ " z = vcat(z'...)\n",
+ "\t\tp = ones(3)\n",
+ "\t\tpi = 1/16\n",
+ "\t\te = [2, 0, 0]\n",
+ "\n",
+ "\t\treturn Dict(\"a\"=>a,\"na\"=>na,\"nc\"=>nc,\"ns\"=>ns,\"nss\"=>nss,\"e\"=>e,\"p\"=>p,\"z\"=>z,\"pi\"=>pi)\n",
+ "end\n",
+ "\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "Dict{String,Any} with 9 entries:\n",
+ " \"na\" => 3\n",
+ " \"e\" => [2, 0, 0]\n",
+ " \"ns\" => 4\n",
+ " \"nc\" => 4\n",
+ " \"nss\" => 16\n",
+ " \"z\" => [1.0 0.72 0.86; 1.0 0.72 0.96; … ; 1.0 1.32 1.06; 1.0 1.32 1.16]\n",
+ " \"a\" => 0.5\n",
+ " \"p\" => [1.0, 1.0, 1.0]\n",
+ " \"pi\" => 0.0625"
+ ]
+ },
+ "execution_count": 2,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "d = data()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\NLopt\\faRdv.ji for NLopt [76087f3c-5699-56af-9a33-bf431cd00edd]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using NLopt"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "u (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 4,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function u(x, data::Dict)\n",
+ " -exp(-data[\"a\"].*x)\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "obj (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 5,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function obj(x::Vector, grad::Vector, data::Dict)\n",
+ " a = data[\"a\"]\n",
+ " z = data[\"z\"]\n",
+ " if length(grad) > 0\n",
+ " grad[1] = -a * u(x[1], data)\n",
+ " for i in 1:3\n",
+ " grad[i+1] = -a .* sum(data[\"pi\"] .* z[:, i] .* u.(z * x[2:end], Ref(data)))\n",
+ " end\n",
+ "end \n",
+ " return u(x[1], data) + sum(data[\"pi\"] .* u.(z * x[2:end], Ref(data)))\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.0\n",
+ " 0.0\n",
+ " 0.0\n",
+ " 0.0"
+ ]
+ },
+ "execution_count": 6,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "grad = zeros(4)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "-0.8280607425183014"
+ ]
+ },
+ "execution_count": 7,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "obj(ones(4), grad, d)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.3032653298563167 \n",
+ " 0.11076504140283401\n",
+ " 0.11021903217359377\n",
+ " 0.11118090020589945"
+ ]
+ },
+ "execution_count": 8,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "grad"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "constr (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 9,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function constr(x::Vector,grad::Vector,data::Dict)\n",
+ " if length(grad) > 0\n",
+ " grad[1] = 1\n",
+ " grad[2:end] = d[\"p\"]\n",
+ " end\n",
+ " return x[1] + sum(d[\"p\"] .* (x[2:end] .- d[\"e\"]))\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "2.0"
+ ]
+ },
+ "execution_count": 10,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "constr(ones(4), zeros(4), d)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "max_NLopt (generic function with 2 methods)"
+ ]
+ },
+ "execution_count": 11,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function max_NLopt(a=0.5)\n",
+ " d = data(a)\n",
+ " opt = Opt(:LD_MMA, 4)\n",
+ " lower_bounds!(opt, [0, -Inf, -Inf, -Inf]) # As cannot sell more then initial endowment\n",
+ " max_objective!(opt, (x,g) -> obj(x,g, d))\n",
+ " inequality_constraint!(opt, (x, g) -> constr(x, g, d))\n",
+ " ftol_rel!(opt, 1e-9)\n",
+ " NLopt.optimize(opt, [0,2,0,0])\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 12,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "(-1.2082140916618629, [1.00797, -1.41217, 0.801299, 1.6029], :FTOL_REACHED)"
+ ]
+ },
+ "execution_count": 12,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "max_NLopt()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 13,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\DataFrames\\AR9oZ.ji for DataFrames [a93c6f00-e57d-5684-b7b6-d8193f3e46c0]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using DataFrames"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 14,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "table_NLopt (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 14,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function table_NLopt()\n",
+ "\t\td = DataFrame(a=[0.5;1.0;5.0],c = zeros(3),omega1=zeros(3),omega2=zeros(3),omega3=zeros(3),fval=zeros(3))\n",
+ "\t\tfor i in 1:nrow(d)\n",
+ "\t\t\txx = max_NLopt(d[i,:a])\n",
+ "\t\t\tfor j in 2:ncol(d)-1\n",
+ "\t\t\t\td[i,j] = xx[2][j-1]\n",
+ "\t\t\tend\n",
+ "\t\t\td[i,end] = xx[1]\n",
+ "\t\tend\n",
+ "\t\treturn d\n",
+ "\tend"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 15,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ " | a | c | omega1 | omega2 | omega3 | fval |
|---|
| Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
|---|
3 rows × 6 columns
| 1 | 0.5 | 1.00797 | -1.41217 | 0.801299 | 1.6029 | -1.20821 |
|---|
| 2 | 1.0 | 1.00401 | -0.205867 | 0.400689 | 0.801171 | -0.732819 |
|---|
| 3 | 5.0 | 1.0008 | 0.758793 | 0.0801445 | 0.160263 | -0.013422 |
|---|
"
+ ],
+ "text/latex": [
+ "\\begin{tabular}{r|cccccc}\n",
+ "\t& a & c & omega1 & omega2 & omega3 & fval\\\\\n",
+ "\t\\hline\n",
+ "\t& Float64 & Float64 & Float64 & Float64 & Float64 & Float64\\\\\n",
+ "\t\\hline\n",
+ "\t1 & 0.5 & 1.00797 & -1.41217 & 0.801299 & 1.6029 & -1.20821 \\\\\n",
+ "\t2 & 1.0 & 1.00401 & -0.205867 & 0.400689 & 0.801171 & -0.732819 \\\\\n",
+ "\t3 & 5.0 & 1.0008 & 0.758793 & 0.0801445 & 0.160263 & -0.013422 \\\\\n",
+ "\\end{tabular}\n"
+ ],
+ "text/plain": [
+ "3×6 DataFrame\n",
+ "│ Row │ a │ c │ omega1 │ omega2 │ omega3 │ fval │\n",
+ "│ │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │\n",
+ "├─────┼─────────┼─────────┼───────────┼───────────┼──────────┼───────────┤\n",
+ "│ 1 │ 0.5 │ 1.00797 │ -1.41217 │ 0.801299 │ 1.6029 │ -1.20821 │\n",
+ "│ 2 │ 1.0 │ 1.00401 │ -0.205867 │ 0.400689 │ 0.801171 │ -0.732819 │\n",
+ "│ 3 │ 5.0 │ 1.0008 │ 0.758793 │ 0.0801445 │ 0.160263 │ -0.013422 │"
+ ]
+ },
+ "execution_count": 15,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "table_NLopt()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m registry at `C:\\Users\\lorenzo\\.julia\\registries\\General`\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m git-repo `https://github.com/JuliaRegistries/General.git`\n",
+ "\u001b[2K\u001b[?25h[1mFetching:\u001b[22m\u001b[39m [========================================>] 100.0 %.0 %1 %.7 %92.7 %\u001b[32m\u001b[1m Resolving\u001b[22m\u001b[39m package versions...\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m ZipFile ────── v0.8.1\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m CSVFiles ───── v0.15.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m PyCall ─────── v1.91.1\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m FileIO ─────── v1.0.6\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m TextParse ──── v0.9.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m Clustering ─── v0.13.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m DoubleFloats ─ v0.7.18\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m Documenter ─── v0.22.0\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Project.toml`\n",
+ " \u001b[90m [49dc2e85]\u001b[39m\u001b[92m + Calculus v0.4.1\u001b[39m\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Manifest.toml`\n",
+ " \u001b[90m [5d742f6a]\u001b[39m\u001b[93m ↑ CSVFiles v0.14.0 ⇒ v0.15.0\u001b[39m\n",
+ " \u001b[90m [aaaa29a8]\u001b[39m\u001b[93m ↑ Clustering v0.12.4 ⇒ v0.13.0\u001b[39m\n",
+ " \u001b[90m [e30172f5]\u001b[39m\u001b[93m ↑ Documenter v0.21.5 ⇒ v0.22.0\u001b[39m\n",
+ " \u001b[90m [497a8b3b]\u001b[39m\u001b[93m ↑ DoubleFloats v0.7.17 ⇒ v0.7.18\u001b[39m\n",
+ " \u001b[90m [5789e2e9]\u001b[39m\u001b[93m ↑ FileIO v1.0.5 ⇒ v1.0.6\u001b[39m\n",
+ " \u001b[90m [14197337]\u001b[39m\u001b[91m - GenericLinearAlgebra v0.1.0\u001b[39m\n",
+ " \u001b[90m [01680d73]\u001b[39m\u001b[92m + GenericSVD v0.2.1\u001b[39m\n",
+ " \u001b[90m [438e738f]\u001b[39m\u001b[93m ↑ PyCall v1.91.0 ⇒ v1.91.1\u001b[39m\n",
+ " \u001b[90m [e0df1984]\u001b[39m\u001b[93m ↑ TextParse v0.8.0 ⇒ v0.9.0\u001b[39m\n",
+ " \u001b[90m [a5390f91]\u001b[39m\u001b[93m ↑ ZipFile v0.8.0 ⇒ v0.8.1\u001b[39m\n",
+ "\u001b[32m\u001b[1m Building\u001b[22m\u001b[39m PyCall ─→ `C:\\Users\\lorenzo\\.julia\\packages\\PyCall\\a5Jd3\\deps\\build.log`\n",
+ "\u001b[32m\u001b[1m Building\u001b[22m\u001b[39m ZipFile → `C:\\Users\\lorenzo\\.julia\\packages\\ZipFile\\YHTbb\\deps\\build.log`\n"
+ ]
+ }
+ ],
+ "source": [
+ "import Pkg; Pkg.add(\"Calculus\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 16,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\Calculus\\G3wEN.ji for Calculus [49dc2e85-a5d0-5ad3-a950-438e2897f1b9]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using Calculus"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 17,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.3032653298574901 \n",
+ " 0.11076504140290903\n",
+ " 0.11021903217009567\n",
+ " 0.11118090021276408"
+ ]
+ },
+ "execution_count": 17,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "Calculus.gradient(arg -> obj(arg, zeros(4), d), ones(4))"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 18,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\JuMP\\DmXqY.ji for JuMP [4076af6c-e467-56ae-b986-b466b2749572]\n",
+ "└ @ Base loading.jl:1184\n",
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\Ipopt\\yMQMo.ji for Ipopt [b6b21f68-93f8-5de0-b562-5493be1d77c9]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using JuMP\n",
+ "using Ipopt"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "\u001b[32m\u001b[1m Status\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Project.toml`\n",
+ " \u001b[90m [aae01518]\u001b[39m\u001b[37m BandedMatrices v0.8.1\u001b[39m\n",
+ " \u001b[90m [6e4b80f9]\u001b[39m\u001b[37m BenchmarkTools v0.4.2\u001b[39m\n",
+ " \u001b[90m [a134a8b2]\u001b[39m\u001b[37m BlackBoxOptim v0.4.0\u001b[39m\n",
+ " \u001b[90m [49dc2e85]\u001b[39m\u001b[37m Calculus v0.4.1\u001b[39m\n",
+ " \u001b[90m [324d7699]\u001b[39m\u001b[37m CategoricalArrays v0.5.2\u001b[39m\n",
+ " \u001b[90m [3da002f7]\u001b[39m\u001b[37m ColorTypes v0.7.5\u001b[39m\n",
+ " \u001b[90m [34da2185]\u001b[39m\u001b[37m Compat v1.4.0\u001b[39m\n",
+ " \u001b[90m [a93c6f00]\u001b[39m\u001b[37m DataFrames v0.17.0\u001b[39m\n",
+ " \u001b[90m [1313f7d8]\u001b[39m\u001b[37m DataFramesMeta v0.4.0\u001b[39m\n",
+ " \u001b[90m [5721bf48]\u001b[39m\u001b[37m DataVoyager v0.3.1\u001b[39m\n",
+ " \u001b[90m [01453d9d]\u001b[39m\u001b[37m DiffEqDiffTools v0.7.1\u001b[39m\n",
+ " \u001b[90m [31c24e10]\u001b[39m\u001b[37m Distributions v0.16.4\u001b[39m\n",
+ " \u001b[90m [63b64609]\u001b[39m\u001b[37m ExamplePackage v0.1.0 [`..\\..\\dev\\ExamplePackage`]\u001b[39m\n",
+ " \u001b[90m [2fe49d83]\u001b[39m\u001b[37m Expectations v1.0.2\u001b[39m\n",
+ " \u001b[90m [442a2c76]\u001b[39m\u001b[37m FastGaussQuadrature v0.3.2\u001b[39m\n",
+ " \u001b[90m [1a297f60]\u001b[39m\u001b[37m FillArrays v0.4.0\u001b[39m\n",
+ " \u001b[90m [9d5cd8c9]\u001b[39m\u001b[37m FixedEffectModels v0.7.1\u001b[39m\n",
+ " \u001b[90m [587475ba]\u001b[39m\u001b[37m Flux v0.7.1\u001b[39m\n",
+ " \u001b[90m [f6369f11]\u001b[39m\u001b[37m ForwardDiff v0.10.2\u001b[39m\n",
+ " \u001b[90m [38e38edf]\u001b[39m\u001b[37m GLM v1.1.0\u001b[39m\n",
+ " \u001b[90m [28b8d3ca]\u001b[39m\u001b[37m GR v0.37.0\u001b[39m\n",
+ " \u001b[90m [2aa259a0]\u001b[39m\u001b[37m HWconstrained v0.1.0 [`C:\\Users\\lorenzo\\.julia\\dev\\HWconstrained`]\u001b[39m\n",
+ " \u001b[90m [d8de2c54]\u001b[39m\u001b[37m HWunconstrained v0.1.0 [`C:\\Users\\lorenzo\\.julia\\dev\\HWunconstrained`]\u001b[39m\n",
+ " \u001b[90m [7073ff75]\u001b[39m\u001b[37m IJulia v1.16.0\u001b[39m\n",
+ " \u001b[90m [43edad99]\u001b[39m\u001b[37m InstantiateFromURL v0.2.1\u001b[39m\n",
+ " \u001b[90m [a98d9a8b]\u001b[39m\u001b[37m Interpolations v0.11.1\u001b[39m\n",
+ " \u001b[90m [b6b21f68]\u001b[39m\u001b[37m Ipopt v0.5.1\u001b[39m\n",
+ " \u001b[90m [4076af6c]\u001b[39m\u001b[37m JuMP v0.18.5\u001b[39m\n",
+ " \u001b[90m [e5e0dc1b]\u001b[39m\u001b[37m Juno v0.5.5\u001b[39m\n",
+ " \u001b[90m [5ab0869b]\u001b[39m\u001b[37m KernelDensity v0.5.1\u001b[39m\n",
+ " \u001b[90m [b964fa9f]\u001b[39m\u001b[37m LaTeXStrings v1.0.3\u001b[39m\n",
+ " \u001b[90m [23fbe1c1]\u001b[39m\u001b[37m Latexify v0.6.0\u001b[39m\n",
+ " \u001b[90m [5078a376]\u001b[39m\u001b[37m LazyArrays v0.5.1\u001b[39m\n",
+ " \u001b[90m [0fc2ff8b]\u001b[39m\u001b[37m LeastSquaresOptim v0.7.1\u001b[39m\n",
+ " \u001b[90m [093fc24a]\u001b[39m\u001b[37m LightGraphs v1.2.0\u001b[39m\n",
+ " \u001b[90m [76087f3c]\u001b[39m\u001b[37m NLopt v0.5.1\u001b[39m\n",
+ " \u001b[90m [2774e3e8]\u001b[39m\u001b[37m NLsolve v3.0.1\u001b[39m\n",
+ " \u001b[90m [429524aa]\u001b[39m\u001b[37m Optim v0.17.2\u001b[39m\n",
+ " \u001b[90m [d96e819e]\u001b[39m\u001b[37m Parameters v0.10.3\u001b[39m\n",
+ " \u001b[90m [14b8a8f1]\u001b[39m\u001b[37m PkgTemplates v0.4.1\u001b[39m\n",
+ " \u001b[90m [91a5bcdd]\u001b[39m\u001b[37m Plots v0.23.0\u001b[39m\n",
+ " \u001b[90m [f27b6e38]\u001b[39m\u001b[37m Polynomials v0.5.2\u001b[39m\n",
+ " \u001b[90m [92933f4c]\u001b[39m\u001b[37m ProgressMeter v0.9.0\u001b[39m\n",
+ " \u001b[90m [1fd47b50]\u001b[39m\u001b[37m QuadGK v2.0.3\u001b[39m\n",
+ " \u001b[90m [fcd29c91]\u001b[39m\u001b[37m QuantEcon v0.15.0\u001b[39m\n",
+ " \u001b[90m [1a8c2f83]\u001b[39m\u001b[37m Query v0.10.1\u001b[39m\n",
+ " \u001b[90m [612083be]\u001b[39m\u001b[37m Queryverse v0.2.0\u001b[39m\n",
+ " \u001b[90m [ce6b1742]\u001b[39m\u001b[37m RDatasets v0.6.1\u001b[39m\n",
+ " \u001b[90m [d519eb52]\u001b[39m\u001b[37m RegressionTables v0.2.0\u001b[39m\n",
+ " \u001b[90m [295af30f]\u001b[39m\u001b[37m Revise v1.0.2\u001b[39m\n",
+ " \u001b[90m [f2b01f46]\u001b[39m\u001b[37m Roots v0.7.4\u001b[39m\n",
+ " \u001b[90m [60ddc479]\u001b[39m\u001b[37m StatPlots v0.9.2\u001b[39m\n",
+ " \u001b[90m [90137ffa]\u001b[39m\u001b[37m StaticArrays v0.10.2\u001b[39m\n",
+ " \u001b[90m [2913bbd2]\u001b[39m\u001b[37m StatsBase v0.27.0\u001b[39m\n",
+ " \u001b[90m [3eaba693]\u001b[39m\u001b[37m StatsModels v0.5.0\u001b[39m\n",
+ " \u001b[90m [37b6cedf]\u001b[39m\u001b[37m Traceur v0.2.1\u001b[39m\n",
+ " \u001b[90m [0ae4a718]\u001b[39m\u001b[37m VegaDatasets v0.5.0\u001b[39m\n",
+ " \u001b[90m [112f6efa]\u001b[39m\u001b[37m VegaLite v0.5.2\u001b[39m\n",
+ " \u001b[90m [37e2e46d]\u001b[39m\u001b[37m LinearAlgebra \u001b[39m\n",
+ " \u001b[90m [3fa0cd96]\u001b[39m\u001b[37m REPL \u001b[39m\n",
+ " \u001b[90m [9a3f8284]\u001b[39m\u001b[37m Random \u001b[39m\n",
+ " \u001b[90m [2f01184e]\u001b[39m\u001b[37m SparseArrays \u001b[39m\n",
+ " \u001b[90m [10745b16]\u001b[39m\u001b[37m Statistics \u001b[39m\n",
+ " \u001b[90m [8dfed614]\u001b[39m\u001b[37m Test \u001b[39m\n"
+ ]
+ }
+ ],
+ "source": [
+ "] st"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 33,
+ "metadata": {},
+ "outputs": [
+ {
+ "ename": "UndefVarError",
+ "evalue": "UndefVarError: with_optimizer not defined",
+ "output_type": "error",
+ "traceback": [
+ "UndefVarError: with_optimizer not defined",
+ "",
+ "Stacktrace:",
+ " [1] max_JuMP(::Float64) at .\\In[33]:3",
+ " [2] max_JuMP() at .\\In[33]:2",
+ " [3] top-level scope at In[33]:11"
+ ]
+ }
+ ],
+ "source": [
+ "function max_JuMP(a=0.5)\n",
+ "\t\td = data(a)\n",
+ "\t\tm = Model(with_optimizer(Ipopt.Optimizer))\n",
+ "\t\t@variable(m, c >= 0)\n",
+ "\t\t@variable(m, x[1:3])\n",
+ "\t\t@NLconstraint(m, c + sum(d[\"p\"][i] * (x[i] - d[\"e\"][i]) for i in 1:3) == 0)\n",
+ "\t\t@NLobjective(m, Max, -exp(-a * c) + sum(-exp(-a * sum(d[\"za\"]'[j][i] * x[i] for i in 1:3)) * d[\"pi\"][j] for j in 1:16))\n",
+ "\t\tJuMP.optimize!(m)\n",
+ "\t\treturn Dict(\"obj\"=>objective_value(m),\"c\"=>value(c),\"omegas\"=>[value(x[i]) for i in 1:length(x)])\n",
+ "\tend\n",
+ "max_JuMP()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 38,
+ "metadata": {},
+ "outputs": [
+ {
+ "ename": "LoadError",
+ "evalue": "syntax: unexpected \"]\"",
+ "output_type": "error",
+ "traceback": [
+ "syntax: unexpected \"]\"",
+ ""
+ ]
+ }
+ ],
+ "source": []
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Julia 1.1.0",
+ "language": "julia",
+ "name": "julia-1.1"
+ },
+ "language_info": {
+ "file_extension": ".jl",
+ "mimetype": "application/julia",
+ "name": "julia",
+ "version": "1.1.0"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/Project.toml b/Project.toml
index a987f1f..b63824c 100644
--- a/Project.toml
+++ b/Project.toml
@@ -3,6 +3,13 @@ uuid = "2aa259a0-4d4c-11e9-0909-ad9b1ef308f5"
authors = ["florian oswald "]
version = "0.1.0"
+[deps]
+DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
+Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
+JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+NLopt = "76087f3c-5699-56af-9a33-bf431cd00edd"
+
[extras]
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
diff --git a/Untitled.ipynb b/Untitled.ipynb
new file mode 100644
index 0000000..4e4b89a
--- /dev/null
+++ b/Untitled.ipynb
@@ -0,0 +1,649 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "data (generic function with 2 methods)"
+ ]
+ },
+ "execution_count": 1,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function data(a=0.5)\n",
+ " a = a\n",
+ "\t\tna = 3\n",
+ "\t\tnc = 4\n",
+ "\t\tns = 4\n",
+ "\t\tnss = 16\n",
+ "\t\tz1 = ones(4)\n",
+ "\t\tz2 = [0.72, 0.92, 1.12, 1.32]\n",
+ "\t\tz3 = [0.86, 0.96, 1.06, 1.16]\n",
+ " z = [16, 3]\n",
+ " z = [[1, i, j] for i in z2 for j in z3] # All possible payoffs\n",
+ " z = vcat(z'...)\n",
+ "\t\tp = ones(3)\n",
+ "\t\tpi = 1/16\n",
+ "\t\te = [2, 0, 0]\n",
+ "\n",
+ "\t\treturn Dict(\"a\"=>a,\"na\"=>na,\"nc\"=>nc,\"ns\"=>ns,\"nss\"=>nss,\"e\"=>e,\"p\"=>p,\"z\"=>z,\"pi\"=>pi)\n",
+ "end\n",
+ "\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "Dict{String,Any} with 9 entries:\n",
+ " \"na\" => 3\n",
+ " \"e\" => [2, 0, 0]\n",
+ " \"ns\" => 4\n",
+ " \"nc\" => 4\n",
+ " \"nss\" => 16\n",
+ " \"z\" => [1.0 0.72 0.86; 1.0 0.72 0.96; … ; 1.0 1.32 1.06; 1.0 1.32 1.16]\n",
+ " \"a\" => 0.5\n",
+ " \"p\" => [1.0, 1.0, 1.0]\n",
+ " \"pi\" => 0.0625"
+ ]
+ },
+ "execution_count": 2,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "d = data()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\NLopt\\faRdv.ji for NLopt [76087f3c-5699-56af-9a33-bf431cd00edd]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using NLopt"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "u (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 4,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function u(x, data::Dict)\n",
+ " -exp(-data[\"a\"].*x)\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "obj (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 5,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function obj(x::Vector, grad::Vector, data::Dict)\n",
+ " a = data[\"a\"]\n",
+ " z = data[\"z\"]\n",
+ " if length(grad) > 0\n",
+ " grad[1] = -a * u(x[1], data)\n",
+ " for i in 1:3\n",
+ " grad[i+1] = -a .* sum(data[\"pi\"] .* z[:, i] .* u.(z * x[2:end], Ref(data)))\n",
+ " end\n",
+ "end \n",
+ " return u(x[1], data) + sum(data[\"pi\"] .* u.(z * x[2:end], Ref(data)))\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.0\n",
+ " 0.0\n",
+ " 0.0\n",
+ " 0.0"
+ ]
+ },
+ "execution_count": 6,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "grad = zeros(4)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "-0.8280607425183014"
+ ]
+ },
+ "execution_count": 7,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "obj(ones(4), grad, d)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.3032653298563167 \n",
+ " 0.11076504140283401\n",
+ " 0.11021903217359377\n",
+ " 0.11118090020589945"
+ ]
+ },
+ "execution_count": 8,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "grad"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "constr (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 9,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function constr(x::Vector,grad::Vector,data::Dict)\n",
+ " if length(grad) > 0\n",
+ " grad[1] = 1\n",
+ " grad[2:end] = d[\"p\"]\n",
+ " end\n",
+ " return x[1] + sum(d[\"p\"] .* (x[2:end] .- d[\"e\"]))\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "2.0"
+ ]
+ },
+ "execution_count": 10,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "constr(ones(4), zeros(4), d)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "max_NLopt (generic function with 2 methods)"
+ ]
+ },
+ "execution_count": 11,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function max_NLopt(a=0.5)\n",
+ " d = data(a)\n",
+ " opt = Opt(:LD_MMA, 4)\n",
+ " lower_bounds!(opt, [0, -Inf, -Inf, -Inf]) # As cannot sell more then initial endowment\n",
+ " max_objective!(opt, (x,g) -> obj(x,g, d))\n",
+ " inequality_constraint!(opt, (x, g) -> constr(x, g, d))\n",
+ " ftol_rel!(opt, 1e-9)\n",
+ " NLopt.optimize(opt, [0,2,0,0])\n",
+ "end"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 12,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "(-1.2082140916618629, [1.00797, -1.41217, 0.801299, 1.6029], :FTOL_REACHED)"
+ ]
+ },
+ "execution_count": 12,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "max_NLopt()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 13,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\DataFrames\\AR9oZ.ji for DataFrames [a93c6f00-e57d-5684-b7b6-d8193f3e46c0]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using DataFrames"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 14,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "table_NLopt (generic function with 1 method)"
+ ]
+ },
+ "execution_count": 14,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "function table_NLopt()\n",
+ "\t\td = DataFrame(a=[0.5;1.0;5.0],c = zeros(3),omega1=zeros(3),omega2=zeros(3),omega3=zeros(3),fval=zeros(3))\n",
+ "\t\tfor i in 1:nrow(d)\n",
+ "\t\t\txx = max_NLopt(d[i,:a])\n",
+ "\t\t\tfor j in 2:ncol(d)-1\n",
+ "\t\t\t\td[i,j] = xx[2][j-1]\n",
+ "\t\t\tend\n",
+ "\t\t\td[i,end] = xx[1]\n",
+ "\t\tend\n",
+ "\t\treturn d\n",
+ "\tend"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 15,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ " | a | c | omega1 | omega2 | omega3 | fval |
|---|
| Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
|---|
3 rows × 6 columns
| 1 | 0.5 | 1.00797 | -1.41217 | 0.801299 | 1.6029 | -1.20821 |
|---|
| 2 | 1.0 | 1.00401 | -0.205867 | 0.400689 | 0.801171 | -0.732819 |
|---|
| 3 | 5.0 | 1.0008 | 0.758793 | 0.0801445 | 0.160263 | -0.013422 |
|---|
"
+ ],
+ "text/latex": [
+ "\\begin{tabular}{r|cccccc}\n",
+ "\t& a & c & omega1 & omega2 & omega3 & fval\\\\\n",
+ "\t\\hline\n",
+ "\t& Float64 & Float64 & Float64 & Float64 & Float64 & Float64\\\\\n",
+ "\t\\hline\n",
+ "\t1 & 0.5 & 1.00797 & -1.41217 & 0.801299 & 1.6029 & -1.20821 \\\\\n",
+ "\t2 & 1.0 & 1.00401 & -0.205867 & 0.400689 & 0.801171 & -0.732819 \\\\\n",
+ "\t3 & 5.0 & 1.0008 & 0.758793 & 0.0801445 & 0.160263 & -0.013422 \\\\\n",
+ "\\end{tabular}\n"
+ ],
+ "text/plain": [
+ "3×6 DataFrame\n",
+ "│ Row │ a │ c │ omega1 │ omega2 │ omega3 │ fval │\n",
+ "│ │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │ \u001b[90mFloat64\u001b[39m │\n",
+ "├─────┼─────────┼─────────┼───────────┼───────────┼──────────┼───────────┤\n",
+ "│ 1 │ 0.5 │ 1.00797 │ -1.41217 │ 0.801299 │ 1.6029 │ -1.20821 │\n",
+ "│ 2 │ 1.0 │ 1.00401 │ -0.205867 │ 0.400689 │ 0.801171 │ -0.732819 │\n",
+ "│ 3 │ 5.0 │ 1.0008 │ 0.758793 │ 0.0801445 │ 0.160263 │ -0.013422 │"
+ ]
+ },
+ "execution_count": 15,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "table_NLopt()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m registry at `C:\\Users\\lorenzo\\.julia\\registries\\General`\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m git-repo `https://github.com/JuliaRegistries/General.git`\n",
+ "\u001b[2K\u001b[?25h[1mFetching:\u001b[22m\u001b[39m [========================================>] 100.0 %.0 %1 %.7 %92.7 %\u001b[32m\u001b[1m Resolving\u001b[22m\u001b[39m package versions...\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m ZipFile ────── v0.8.1\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m CSVFiles ───── v0.15.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m PyCall ─────── v1.91.1\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m FileIO ─────── v1.0.6\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m TextParse ──── v0.9.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m Clustering ─── v0.13.0\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m DoubleFloats ─ v0.7.18\n",
+ "\u001b[32m\u001b[1m Installed\u001b[22m\u001b[39m Documenter ─── v0.22.0\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Project.toml`\n",
+ " \u001b[90m [49dc2e85]\u001b[39m\u001b[92m + Calculus v0.4.1\u001b[39m\n",
+ "\u001b[32m\u001b[1m Updating\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Manifest.toml`\n",
+ " \u001b[90m [5d742f6a]\u001b[39m\u001b[93m ↑ CSVFiles v0.14.0 ⇒ v0.15.0\u001b[39m\n",
+ " \u001b[90m [aaaa29a8]\u001b[39m\u001b[93m ↑ Clustering v0.12.4 ⇒ v0.13.0\u001b[39m\n",
+ " \u001b[90m [e30172f5]\u001b[39m\u001b[93m ↑ Documenter v0.21.5 ⇒ v0.22.0\u001b[39m\n",
+ " \u001b[90m [497a8b3b]\u001b[39m\u001b[93m ↑ DoubleFloats v0.7.17 ⇒ v0.7.18\u001b[39m\n",
+ " \u001b[90m [5789e2e9]\u001b[39m\u001b[93m ↑ FileIO v1.0.5 ⇒ v1.0.6\u001b[39m\n",
+ " \u001b[90m [14197337]\u001b[39m\u001b[91m - GenericLinearAlgebra v0.1.0\u001b[39m\n",
+ " \u001b[90m [01680d73]\u001b[39m\u001b[92m + GenericSVD v0.2.1\u001b[39m\n",
+ " \u001b[90m [438e738f]\u001b[39m\u001b[93m ↑ PyCall v1.91.0 ⇒ v1.91.1\u001b[39m\n",
+ " \u001b[90m [e0df1984]\u001b[39m\u001b[93m ↑ TextParse v0.8.0 ⇒ v0.9.0\u001b[39m\n",
+ " \u001b[90m [a5390f91]\u001b[39m\u001b[93m ↑ ZipFile v0.8.0 ⇒ v0.8.1\u001b[39m\n",
+ "\u001b[32m\u001b[1m Building\u001b[22m\u001b[39m PyCall ─→ `C:\\Users\\lorenzo\\.julia\\packages\\PyCall\\a5Jd3\\deps\\build.log`\n",
+ "\u001b[32m\u001b[1m Building\u001b[22m\u001b[39m ZipFile → `C:\\Users\\lorenzo\\.julia\\packages\\ZipFile\\YHTbb\\deps\\build.log`\n"
+ ]
+ }
+ ],
+ "source": [
+ "import Pkg; Pkg.add(\"Calculus\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 16,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\Calculus\\G3wEN.ji for Calculus [49dc2e85-a5d0-5ad3-a950-438e2897f1b9]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using Calculus"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 17,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "4-element Array{Float64,1}:\n",
+ " 0.3032653298574901 \n",
+ " 0.11076504140290903\n",
+ " 0.11021903217009567\n",
+ " 0.11118090021276408"
+ ]
+ },
+ "execution_count": 17,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "Calculus.gradient(arg -> obj(arg, zeros(4), d), ones(4))"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 18,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\JuMP\\DmXqY.ji for JuMP [4076af6c-e467-56ae-b986-b466b2749572]\n",
+ "└ @ Base loading.jl:1184\n",
+ "┌ Info: Recompiling stale cache file C:\\Users\\lorenzo\\.julia\\compiled\\v1.1\\Ipopt\\yMQMo.ji for Ipopt [b6b21f68-93f8-5de0-b562-5493be1d77c9]\n",
+ "└ @ Base loading.jl:1184\n"
+ ]
+ }
+ ],
+ "source": [
+ "using JuMP\n",
+ "using Ipopt"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "\u001b[32m\u001b[1m Status\u001b[22m\u001b[39m `C:\\Users\\lorenzo\\.julia\\environments\\v1.1\\Project.toml`\n",
+ " \u001b[90m [aae01518]\u001b[39m\u001b[37m BandedMatrices v0.8.1\u001b[39m\n",
+ " \u001b[90m [6e4b80f9]\u001b[39m\u001b[37m BenchmarkTools v0.4.2\u001b[39m\n",
+ " \u001b[90m [a134a8b2]\u001b[39m\u001b[37m BlackBoxOptim v0.4.0\u001b[39m\n",
+ " \u001b[90m [49dc2e85]\u001b[39m\u001b[37m Calculus v0.4.1\u001b[39m\n",
+ " \u001b[90m [324d7699]\u001b[39m\u001b[37m CategoricalArrays v0.5.2\u001b[39m\n",
+ " \u001b[90m [3da002f7]\u001b[39m\u001b[37m ColorTypes v0.7.5\u001b[39m\n",
+ " \u001b[90m [34da2185]\u001b[39m\u001b[37m Compat v1.4.0\u001b[39m\n",
+ " \u001b[90m [a93c6f00]\u001b[39m\u001b[37m DataFrames v0.17.0\u001b[39m\n",
+ " \u001b[90m [1313f7d8]\u001b[39m\u001b[37m DataFramesMeta v0.4.0\u001b[39m\n",
+ " \u001b[90m [5721bf48]\u001b[39m\u001b[37m DataVoyager v0.3.1\u001b[39m\n",
+ " \u001b[90m [01453d9d]\u001b[39m\u001b[37m DiffEqDiffTools v0.7.1\u001b[39m\n",
+ " \u001b[90m [31c24e10]\u001b[39m\u001b[37m Distributions v0.16.4\u001b[39m\n",
+ " \u001b[90m [63b64609]\u001b[39m\u001b[37m ExamplePackage v0.1.0 [`..\\..\\dev\\ExamplePackage`]\u001b[39m\n",
+ " \u001b[90m [2fe49d83]\u001b[39m\u001b[37m Expectations v1.0.2\u001b[39m\n",
+ " \u001b[90m [442a2c76]\u001b[39m\u001b[37m FastGaussQuadrature v0.3.2\u001b[39m\n",
+ " \u001b[90m [1a297f60]\u001b[39m\u001b[37m FillArrays v0.4.0\u001b[39m\n",
+ " \u001b[90m [9d5cd8c9]\u001b[39m\u001b[37m FixedEffectModels v0.7.1\u001b[39m\n",
+ " \u001b[90m [587475ba]\u001b[39m\u001b[37m Flux v0.7.1\u001b[39m\n",
+ " \u001b[90m [f6369f11]\u001b[39m\u001b[37m ForwardDiff v0.10.2\u001b[39m\n",
+ " \u001b[90m [38e38edf]\u001b[39m\u001b[37m GLM v1.1.0\u001b[39m\n",
+ " \u001b[90m [28b8d3ca]\u001b[39m\u001b[37m GR v0.37.0\u001b[39m\n",
+ " \u001b[90m [2aa259a0]\u001b[39m\u001b[37m HWconstrained v0.1.0 [`C:\\Users\\lorenzo\\.julia\\dev\\HWconstrained`]\u001b[39m\n",
+ " \u001b[90m [d8de2c54]\u001b[39m\u001b[37m HWunconstrained v0.1.0 [`C:\\Users\\lorenzo\\.julia\\dev\\HWunconstrained`]\u001b[39m\n",
+ " \u001b[90m [7073ff75]\u001b[39m\u001b[37m IJulia v1.16.0\u001b[39m\n",
+ " \u001b[90m [43edad99]\u001b[39m\u001b[37m InstantiateFromURL v0.2.1\u001b[39m\n",
+ " \u001b[90m [a98d9a8b]\u001b[39m\u001b[37m Interpolations v0.11.1\u001b[39m\n",
+ " \u001b[90m [b6b21f68]\u001b[39m\u001b[37m Ipopt v0.5.1\u001b[39m\n",
+ " \u001b[90m [4076af6c]\u001b[39m\u001b[37m JuMP v0.18.5\u001b[39m\n",
+ " \u001b[90m [e5e0dc1b]\u001b[39m\u001b[37m Juno v0.5.5\u001b[39m\n",
+ " \u001b[90m [5ab0869b]\u001b[39m\u001b[37m KernelDensity v0.5.1\u001b[39m\n",
+ " \u001b[90m [b964fa9f]\u001b[39m\u001b[37m LaTeXStrings v1.0.3\u001b[39m\n",
+ " \u001b[90m [23fbe1c1]\u001b[39m\u001b[37m Latexify v0.6.0\u001b[39m\n",
+ " \u001b[90m [5078a376]\u001b[39m\u001b[37m LazyArrays v0.5.1\u001b[39m\n",
+ " \u001b[90m [0fc2ff8b]\u001b[39m\u001b[37m LeastSquaresOptim v0.7.1\u001b[39m\n",
+ " \u001b[90m [093fc24a]\u001b[39m\u001b[37m LightGraphs v1.2.0\u001b[39m\n",
+ " \u001b[90m [76087f3c]\u001b[39m\u001b[37m NLopt v0.5.1\u001b[39m\n",
+ " \u001b[90m [2774e3e8]\u001b[39m\u001b[37m NLsolve v3.0.1\u001b[39m\n",
+ " \u001b[90m [429524aa]\u001b[39m\u001b[37m Optim v0.17.2\u001b[39m\n",
+ " \u001b[90m [d96e819e]\u001b[39m\u001b[37m Parameters v0.10.3\u001b[39m\n",
+ " \u001b[90m [14b8a8f1]\u001b[39m\u001b[37m PkgTemplates v0.4.1\u001b[39m\n",
+ " \u001b[90m [91a5bcdd]\u001b[39m\u001b[37m Plots v0.23.0\u001b[39m\n",
+ " \u001b[90m [f27b6e38]\u001b[39m\u001b[37m Polynomials v0.5.2\u001b[39m\n",
+ " \u001b[90m [92933f4c]\u001b[39m\u001b[37m ProgressMeter v0.9.0\u001b[39m\n",
+ " \u001b[90m [1fd47b50]\u001b[39m\u001b[37m QuadGK v2.0.3\u001b[39m\n",
+ " \u001b[90m [fcd29c91]\u001b[39m\u001b[37m QuantEcon v0.15.0\u001b[39m\n",
+ " \u001b[90m [1a8c2f83]\u001b[39m\u001b[37m Query v0.10.1\u001b[39m\n",
+ " \u001b[90m [612083be]\u001b[39m\u001b[37m Queryverse v0.2.0\u001b[39m\n",
+ " \u001b[90m [ce6b1742]\u001b[39m\u001b[37m RDatasets v0.6.1\u001b[39m\n",
+ " \u001b[90m [d519eb52]\u001b[39m\u001b[37m RegressionTables v0.2.0\u001b[39m\n",
+ " \u001b[90m [295af30f]\u001b[39m\u001b[37m Revise v1.0.2\u001b[39m\n",
+ " \u001b[90m [f2b01f46]\u001b[39m\u001b[37m Roots v0.7.4\u001b[39m\n",
+ " \u001b[90m [60ddc479]\u001b[39m\u001b[37m StatPlots v0.9.2\u001b[39m\n",
+ " \u001b[90m [90137ffa]\u001b[39m\u001b[37m StaticArrays v0.10.2\u001b[39m\n",
+ " \u001b[90m [2913bbd2]\u001b[39m\u001b[37m StatsBase v0.27.0\u001b[39m\n",
+ " \u001b[90m [3eaba693]\u001b[39m\u001b[37m StatsModels v0.5.0\u001b[39m\n",
+ " \u001b[90m [37b6cedf]\u001b[39m\u001b[37m Traceur v0.2.1\u001b[39m\n",
+ " \u001b[90m [0ae4a718]\u001b[39m\u001b[37m VegaDatasets v0.5.0\u001b[39m\n",
+ " \u001b[90m [112f6efa]\u001b[39m\u001b[37m VegaLite v0.5.2\u001b[39m\n",
+ " \u001b[90m [37e2e46d]\u001b[39m\u001b[37m LinearAlgebra \u001b[39m\n",
+ " \u001b[90m [3fa0cd96]\u001b[39m\u001b[37m REPL \u001b[39m\n",
+ " \u001b[90m [9a3f8284]\u001b[39m\u001b[37m Random \u001b[39m\n",
+ " \u001b[90m [2f01184e]\u001b[39m\u001b[37m SparseArrays \u001b[39m\n",
+ " \u001b[90m [10745b16]\u001b[39m\u001b[37m Statistics \u001b[39m\n",
+ " \u001b[90m [8dfed614]\u001b[39m\u001b[37m Test \u001b[39m\n"
+ ]
+ }
+ ],
+ "source": [
+ "] st"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 33,
+ "metadata": {},
+ "outputs": [
+ {
+ "ename": "UndefVarError",
+ "evalue": "UndefVarError: with_optimizer not defined",
+ "output_type": "error",
+ "traceback": [
+ "UndefVarError: with_optimizer not defined",
+ "",
+ "Stacktrace:",
+ " [1] max_JuMP(::Float64) at .\\In[33]:3",
+ " [2] max_JuMP() at .\\In[33]:2",
+ " [3] top-level scope at In[33]:11"
+ ]
+ }
+ ],
+ "source": [
+ "function max_JuMP(a=0.5)\n",
+ "\t\td = data(a)\n",
+ "\t\tm = Model(with_optimizer(Ipopt.Optimizer))\n",
+ "\t\t@variable(m, c >= 0)\n",
+ "\t\t@variable(m, x[1:3])\n",
+ "\t\t@NLconstraint(m, c + sum(d[\"p\"][i] * (x[i] - d[\"e\"][i]) for i in 1:3) == 0)\n",
+ "\t\t@NLobjective(m, Max, -exp(-a * c) + sum(-exp(-a * sum(d[\"za\"]'[j][i] * x[i] for i in 1:3)) * d[\"pi\"][j] for j in 1:16))\n",
+ "\t\tJuMP.optimize!(m)\n",
+ "\t\treturn Dict(\"obj\"=>objective_value(m),\"c\"=>value(c),\"omegas\"=>[value(x[i]) for i in 1:length(x)])\n",
+ "\tend\n",
+ "max_JuMP()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 38,
+ "metadata": {},
+ "outputs": [
+ {
+ "ename": "LoadError",
+ "evalue": "syntax: unexpected \"]\"",
+ "output_type": "error",
+ "traceback": [
+ "syntax: unexpected \"]\"",
+ ""
+ ]
+ }
+ ],
+ "source": []
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Julia 1.1.0",
+ "language": "julia",
+ "name": "julia-1.1"
+ },
+ "language_info": {
+ "file_extension": ".jl",
+ "mimetype": "application/julia",
+ "name": "julia",
+ "version": "1.1.0"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/src/HWconstrained.jl b/src/HWconstrained.jl
index 3b6faa4..79afdc3 100644
--- a/src/HWconstrained.jl
+++ b/src/HWconstrained.jl
@@ -7,25 +7,37 @@ greet() = print("Hello World!")
export data, table_NLopt, table_JuMP
- function data(a=0.5)
-
-
-
-
-
-
-
-
-
+function data(a=0.5)
+ a = a
+ na = 3
+ nc = 4 # Budget constraint + 3 nonnegativity constraints
+ ns = 4
+ nss = 16
+ z1 = ones(4)
+ z2 = [0.72, 0.92, 1.12, 1.32]
+ z3 = [0.86, 0.96, 1.06, 1.16]
+ z = [(1, i, j) for i in z2 for j in z3] # All possible payoffs
+ p = ones(3)
+ pi = repeat([1/16], nss)
+ e = [2, 0, 0]
+
return Dict("a"=>a,"na"=>na,"nc"=>nc,"ns"=>ns,"nss"=>nss,"e"=>e,"p"=>p,"z"=>z,"pi"=>pi)
- end
+end
function max_JuMP(a=0.5)
-
- return Dict("obj"=>objective_value(m),"c"=>value(c),"omegas"=>[value(omega[i]) for i in 1:length(omega)])
+ d = data(a)
+ m = Model(with_optimizer(Ipopt.Optimizer))
+
+ @variable(m, c>=0)
+ @variable(m, x[1:3])
+ @NLconstraint(m, c + sum(d["p"][i] * (x[i] - d["e"][i]) for i in 1:3) == 0 )
+ @NLobjective(m, Max, -exp(-a * c) + sum(-exp(-a * sum(d["z"][j][i] * x[i] for i in 1:3)) * d["pi"][j] for j in 1:16))
+ JuMP.optimize!(m)
+ return Dict("obj"=>objective_value(m),"c"=>value(c),"omegas"=>[value(x[i]) for i in 1:length(x)])
end
+
function table_JuMP()
d = DataFrame(a=[0.5;1.0;5.0],c = zeros(3),omega1=zeros(3),omega2=zeros(3),omega3=zeros(3),fval=zeros(3))
for i in 1:nrow(d)
@@ -39,67 +51,42 @@ greet() = print("Hello World!")
return d
end
-
-
-
-
-
-
-
- function obj(x::Vector,grad::Vector,data::Dict)
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- end
- function constr(x::Vector,grad::Vector,data::Dict)
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- end
- function max_NLopt(a=0.5)
-
-
-
-
-
-
-
-
- end
+
+
+function obj(x::Vector, grad::Vector, data::Dict)
+ a = data["a"]
+ z = data["z"]
+ u(x) = -exp(-a * x)
+ if length(grad) > 0
+ grad[1] = -a * u(x[1], data)
+ for i in 1:3
+ grad[i+1] = -a .* sum(data["pi"] .* z[:, i] .* u.(z * x[2:end], Ref(data)))
+ end
+end
+ return u(x[1], data) + sum(data["pi"] .* u.(z * x[2:end], Ref(data)))
+end
+
+
+function constr(x::Vector,grad::Vector,data::Dict)
+ if length(grad) > 0
+ grad[1] = 1
+ grad[2:end] = d["p"]
+ end
+ return x[1] + sum(d["p"] .* (x[2:end] .- d["e"]))
+end
+
+function max_NLopt(a=0.5)
+ d = data(a)
+ opt = Opt(:LD_MMA, 4)
+ lower_bounds!(opt, [0, -Inf, -Inf, -Inf]) # As cannot sell more then initial endowment
+ max_objective!(opt, (x,g) -> obj(x,g, d))
+ inequality_constraint!(opt, (x, g) -> constr(x, g, d))
+ ftol_rel!(opt, 1e-9)
+ NLopt.optimize(opt, [0,2,0,0])
+end
+
function table_NLopt()
d = DataFrame(a=[0.5;1.0;5.0],c = zeros(3),omega1=zeros(3),omega2=zeros(3),omega3=zeros(3),fval=zeros(3))
diff --git a/test/runtests.jl b/test/runtests.jl
index 1793ffc..5476e99 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -1,26 +1,33 @@
using HWconstrained
using Test
+using Calculus
@testset "HWconstrained.jl" begin
@testset "testing components" begin
@testset "tests gradient of objective function" begin
-
-
-
-
-
-
+ d = data()
+ grad = zeros(4)
+ obj(ones(4), zeros(4), d)
+ @test Calculus.gradient(arg -> obj(arg, zeros(4), d), ones(4)) == grad
+
+
+
+
+
end
@testset "tests gradient of constraint function" begin
-
-
-
-
-
-
+ d = data()
+ grad = zeros(4)
+ constr(ones(4), zeros(4), d)
+ @test Calculus.gradient(arg -> constr(arg, zeros(4), d), ones(4)) == grad
+
+
+
+
+
end
end