Skip to content

Commit

Permalink
Merge pull request #2 from Foggalong/reframing
Browse files Browse the repository at this point in the history 
Reframing as return maximization
  • Loading branch information
@Foggalong
Foggalong authored Apr 25, 2024
2 parents 83ee395 + c82f19e commit f21e77f
Show file tree
Hide file tree
Showing 4 changed files with 1,139 additions and 277 deletions.
139 changes: 67 additions & 72 deletions Example2/n2-test.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@
"np.set_printoptions(threshold=np.inf)\n",
"\n",
"# only show numpy output to five decimal places\n",
"np.set_printoptions(formatter={'float_kind':\"{:.5f}\".format})"
"np.set_printoptions(formatter={'float_kind':\"{:.8f}\".format})"
]
},
{
Expand Down Expand Up @@ -112,12 +112,12 @@
"Thread count: 4 physical cores, 8 logical processors, using up to 8 threads\n",
"\n",
"Optimize a model with 1 rows, 2 columns and 2 nonzeros\n",
"Model fingerprint: 0x37f20dc3\n",
"Model fingerprint: 0xd1ce5e64\n",
"Model has 2 quadratic objective terms\n",
"Coefficient statistics:\n",
" Matrix range [1e+00, 1e+00]\n",
" Objective range [5e-01, 1e+00]\n",
" QObjective range [1e+00, 1e+00]\n",
" Objective range [1e+00, 2e+00]\n",
" QObjective range [5e-01, 5e-01]\n",
" Bounds range [0e+00, 0e+00]\n",
" RHS range [1e+00, 1e+00]\n",
"Presolve time: 0.01s\n",
Expand All @@ -133,19 +133,19 @@
"\n",
" Objective Residual\n",
"Iter Primal Dual Primal Dual Compl Time\n",
" 0 9.99499500e+05 -1.00100025e+06 2.00e+03 0.00e+00 1.00e+06 0s\n",
" 1 -4.99749639e-01 -9.99502094e+02 1.00e+00 3.41e-13 9.99e+02 0s\n",
" 2 -5.00116624e-01 -2.65997154e+02 1.00e-06 6.99e-15 1.33e+02 0s\n",
" 3 -5.00816469e-01 -1.01361177e+00 9.33e-10 0.00e+00 2.56e-01 0s\n",
" 4 -5.58837435e-01 -5.96696571e-01 3.43e-12 0.00e+00 1.89e-02 0s\n",
" 5 -5.62490657e-01 -5.64212774e-01 0.00e+00 0.00e+00 8.61e-04 0s\n",
" 6 -5.62500000e-01 -5.62501707e-01 0.00e+00 0.00e+00 8.53e-07 0s\n",
" 7 -5.62500000e-01 -5.62500002e-01 0.00e+00 0.00e+00 8.53e-10 0s\n",
" 0 -4.97498625e+05 5.00500125e+05 2.00e+03 0.00e+00 9.99e+05 0s\n",
" 1 2.66553322e+00 9.98756178e+02 1.11e+00 0.00e+00 1.05e+03 0s\n",
" 2 1.40062213e+00 2.93466464e+02 1.11e-06 0.00e+00 1.46e+02 0s\n",
" 3 1.40285882e+00 2.15287433e+00 1.35e-09 0.00e+00 3.75e-01 0s\n",
" 4 1.67958671e+00 1.77570144e+00 1.33e-15 0.00e+00 4.81e-02 0s\n",
" 5 1.74615363e+00 1.75032001e+00 0.00e+00 2.01e-16 2.08e-03 0s\n",
" 6 1.74999618e+00 1.75000030e+00 0.00e+00 0.00e+00 2.06e-06 0s\n",
" 7 1.75000000e+00 1.75000000e+00 0.00e+00 0.00e+00 2.06e-09 0s\n",
"\n",
"Barrier solved model in 7 iterations and 0.02 seconds (0.00 work units)\n",
"Optimal objective -5.62500000e-01\n",
"Optimal objective 1.75000000e+00\n",
"\n",
"Standard: [0.25000 0.75000]\n"
"Standard: [0.00000001 0.99999999]\n"
]
}
],
Expand All @@ -164,8 +164,8 @@
"\n",
"# define the objective functions for standard problem\n",
"model_std.setObjective(\n",
" 0.5*w_std@(sigma@w_std) - lam*w_std.transpose()@mubar,\n",
"GRB.MINIMIZE)\n",
" w_std.transpose()@mubar - (lam/2)*w_std.transpose()@(sigma@w_std),\n",
"GRB.MAXIMIZE)\n",
"\n",
"# add sum-to-one constraint\n",
"model_std.addConstr(np.ones([1,2]) @ w_std == 1, name=\"sum-to-one\")\n",
Expand Down Expand Up @@ -196,44 +196,43 @@
"CPU model: Intel(R) Core(TM) i5-8350U CPU @ 1.70GHz, instruction set [SSE2|AVX|AVX2]\n",
"Thread count: 4 physical cores, 8 logical processors, using up to 8 threads\n",
"\n",
"Optimize a model with 2 rows, 3 columns and 3 nonzeros\n",
"Model fingerprint: 0xb83c03fb\n",
"Optimize a model with 1 rows, 3 columns and 2 nonzeros\n",
"Model fingerprint: 0xea4b60b6\n",
"Model has 2 quadratic objective terms\n",
"Model has 2 quadratic constraints\n",
"Model has 1 quadratic constraint\n",
"Coefficient statistics:\n",
" Matrix range [1e+00, 1e+00]\n",
" QMatrix range [1e-01, 4e+00]\n",
" Objective range [5e-01, 1e+00]\n",
" QObjective range [1e+00, 1e+00]\n",
" Objective range [1e+00, 2e+00]\n",
" QObjective range [5e-01, 5e-01]\n",
" Bounds range [0e+00, 0e+00]\n",
" RHS range [1e+00, 1e+00]\n",
"Presolve removed 1 rows and 0 columns\n",
"Presolve time: 0.01s\n",
"Presolved: 8 rows, 8 columns, 13 nonzeros\n",
"Presolved model has 3 second-order cone constraints\n",
"Presolved: 6 rows, 7 columns, 10 nonzeros\n",
"Presolved model has 2 second-order cone constraints\n",
"Ordering time: 0.00s\n",
"\n",
"Barrier statistics:\n",
" AA' NZ : 1.900e+01\n",
" Factor NZ : 3.600e+01\n",
" Factor Ops : 2.040e+02 (less than 1 second per iteration)\n",
" AA' NZ : 1.500e+01\n",
" Factor NZ : 2.100e+01\n",
" Factor Ops : 9.100e+01 (less than 1 second per iteration)\n",
" Threads : 1\n",
"\n",
" Objective Residual\n",
"Iter Primal Dual Primal Dual Compl Time\n",
" 0 -6.30136986e-01 -6.30136986e-01 1.61e+00 8.66e-01 7.30e-01 0s\n",
" 1 -1.14616869e+00 -1.60124001e+00 4.19e-01 9.52e-07 1.64e-01 0s\n",
" 2 -6.51626451e-01 -8.34836457e-01 4.61e-07 1.05e-12 1.67e-02 0s\n",
" 3 -6.72758817e-01 -6.83613315e-01 2.69e-08 3.87e-14 9.87e-04 0s\n",
" 4 -6.73595613e-01 -6.73729333e-01 2.55e-14 2.00e-15 1.22e-05 0s\n",
" 5 -6.73606055e-01 -6.73616513e-01 1.05e-12 9.30e-13 9.51e-07 0s\n",
" 6 -6.73611064e-01 -6.73611270e-01 6.32e-11 1.09e-11 1.88e-08 0s\n",
" 0 -5.73852130e-02 1.22656250e+00 8.97e-01 1.41e-01 3.51e-01 0s\n",
" 1 5.91715790e-01 9.52821847e-01 2.13e-01 1.51e-07 6.38e-02 0s\n",
" 2 5.15372316e-01 5.96285470e-01 2.34e-07 6.30e-10 8.99e-03 0s\n",
" 3 5.49225634e-01 5.53202935e-01 4.81e-10 1.05e-15 4.42e-04 0s\n",
" 4 5.52127256e-01 5.52382277e-01 1.12e-13 5.55e-16 2.83e-05 0s\n",
" 5 5.52206847e-01 5.52232334e-01 3.52e-12 3.66e-15 2.83e-06 0s\n",
" 6 5.52216328e-01 5.52218236e-01 1.38e-11 2.96e-14 2.12e-07 0s\n",
"\n",
"Barrier solved model in 6 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective -6.73611064e-01\n",
"Optimal objective 5.52216328e-01\n",
"\n",
"Standard: [0.25000 0.75000]\n",
"Robust: [0.41660 0.58340]\n"
"Standard: [0.00000001 0.99999999]\n",
"Robust: [0.83155134 0.16844866]\n"
]
}
],
Expand All @@ -245,18 +244,17 @@
"model_rbs = gp.Model(\"n02robust\")\n",
"\n",
"# initialise w for both models, z for robust model \n",
"w_rbs = model_rbs.addMVar(shape=n, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z_rbs = model_rbs.addVar(name=\"z\")\n",
"w_rbs = model_rbs.addMVar(shape=n, lb=0.0, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z_rbs = model_rbs.addVar(lb=0.0, name=\"z\")\n",
"\n",
"model_rbs.setObjective(\n",
" 0.5*w_rbs@(sigma@w_rbs) - lam*w_rbs.transpose()@mubar - kappa*z_rbs,\n",
"GRB.MINIMIZE)\n",
" w_rbs.transpose()@mubar - (lam/2)*w_rbs.transpose()@(sigma@w_rbs) - kappa*z_rbs,\n",
"GRB.MAXIMIZE)\n",
"\n",
"# add sum-to-half constraints to both models\n",
"model_rbs.addConstr(np.ones([1,2]) @ w_std == 1, name=\"sum-to-one\")\n",
"# add quadratic uncertainty constraint to the robust model\n",
"model_rbs.addConstr(z_rbs**2 <= np.inner(w_rbs, omega@w_rbs), name=\"uncertainty\")\n",
"model_rbs.addConstr(z_rbs >= 0, name=\"z positive\")\n",
"model_rbs.addConstr(z_rbs**2 >= w_rbs.transpose()@omega@w_rbs, name=\"uncertainty\")\n",
"\n",
"# solve problem with Gurobi, print alongside standard solution for comparison\n",
"model_rbs.optimize()\n",
Expand All @@ -268,7 +266,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that if $\\kappa = 0$, the standard solution not _exactly_ recovered, only accurate to four decimal places."
"Note that if $\\kappa = 0$, the standard solution not _exactly_ recovered, only accurate to seven decimal places."
]
},
{
Expand All @@ -285,44 +283,42 @@
"CPU model: Intel(R) Core(TM) i5-8350U CPU @ 1.70GHz, instruction set [SSE2|AVX|AVX2]\n",
"Thread count: 4 physical cores, 8 logical processors, using up to 8 threads\n",
"\n",
"Optimize a model with 2 rows, 3 columns and 3 nonzeros\n",
"Model fingerprint: 0x48b37c4b\n",
"Optimize a model with 1 rows, 3 columns and 2 nonzeros\n",
"Model fingerprint: 0xcd1e4542\n",
"Model has 2 quadratic objective terms\n",
"Model has 2 quadratic constraints\n",
"Model has 1 quadratic constraint\n",
"Coefficient statistics:\n",
" Matrix range [1e+00, 1e+00]\n",
" QMatrix range [1e-01, 4e+00]\n",
" Objective range [5e-01, 1e+00]\n",
" QObjective range [1e+00, 1e+00]\n",
" Objective range [1e+00, 2e+00]\n",
" QObjective range [5e-01, 5e-01]\n",
" Bounds range [0e+00, 0e+00]\n",
" RHS range [1e+00, 1e+00]\n",
"Presolve removed 1 rows and 0 columns\n",
"Presolve time: 0.01s\n",
"Presolved: 8 rows, 8 columns, 13 nonzeros\n",
"Presolved model has 3 second-order cone constraints\n",
"Presolved: 6 rows, 7 columns, 10 nonzeros\n",
"Presolved model has 2 second-order cone constraints\n",
"Ordering time: 0.00s\n",
"\n",
"Barrier statistics:\n",
" AA' NZ : 1.900e+01\n",
" Factor NZ : 3.600e+01\n",
" Factor Ops : 2.040e+02 (less than 1 second per iteration)\n",
" AA' NZ : 1.500e+01\n",
" Factor NZ : 2.100e+01\n",
" Factor Ops : 9.100e+01 (less than 1 second per iteration)\n",
" Threads : 1\n",
"\n",
" Objective Residual\n",
"Iter Primal Dual Primal Dual Compl Time\n",
" 0 -6.30136986e-01 -6.30136986e-01 1.67e+00 4.15e-01 3.76e-01 0s\n",
" 1 -7.63463502e-01 -9.90451466e-01 3.34e-01 4.57e-07 6.75e-02 0s\n",
" 2 -5.29769678e-01 -6.77289143e-01 3.67e-07 5.02e-13 1.34e-02 0s\n",
" 3 -5.47218963e-01 -5.75552518e-01 4.04e-13 3.33e-16 2.58e-03 0s\n",
" 4 -5.61928380e-01 -5.63204951e-01 4.66e-15 1.11e-16 1.16e-04 0s\n",
" 5 -5.62460343e-01 -5.62501958e-01 9.44e-14 2.22e-16 3.78e-06 0s\n",
" 6 -5.62499955e-01 -5.62500028e-01 2.26e-11 6.10e-14 6.60e-09 0s\n",
" 0 1.22656250e+00 1.22656250e+00 1.08e+00 6.94e-01 4.95e-01 0s\n",
" 1 1.61470687e+00 1.61496484e+00 5.67e-02 1.55e-01 6.97e-02 0s\n",
" 2 1.67968985e+00 1.70168242e+00 6.24e-08 2.02e-02 1.10e-02 0s\n",
" 3 1.74763426e+00 1.75020504e+00 7.06e-14 2.52e-05 2.97e-04 0s\n",
" 4 1.74999548e+00 1.75001176e+00 1.87e-14 2.84e-08 1.82e-06 0s\n",
" 5 1.74999978e+00 1.75000067e+00 4.92e-13 9.33e-10 9.90e-08 0s\n",
"\n",
"Barrier solved model in 6 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective -5.62499955e-01\n",
"Barrier solved model in 5 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective 1.74999978e+00\n",
"\n",
"Standard: [0.2500000031 0.7499999969]\n",
"Robust: [0.2499847190 0.7500152810]\n"
"Standard: [0.0000000076 0.9999999924]\n",
"Robust: [0.0000001531 0.9999998469]\n"
]
}
],
Expand All @@ -334,18 +330,17 @@
"model_rbs = gp.Model(\"n02robust_k0\")\n",
"\n",
"# initialise w for both models, z for robust model \n",
"w_rbs = model_rbs.addMVar(shape=n, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z_rbs = model_rbs.addVar(name=\"z\")\n",
"w_rbs = model_rbs.addMVar(shape=n, lb=0.0, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z_rbs = model_rbs.addVar(lb=0.0, name=\"z\")\n",
"\n",
"model_rbs.setObjective(\n",
" 0.5*w_rbs@(sigma@w_rbs) - lam*w_rbs.transpose()@mubar - kappa*z_rbs,\n",
"GRB.MINIMIZE)\n",
" w_rbs.transpose()@mubar - (lam/2)*w_rbs.transpose()@(sigma@w_rbs) - kappa*z_rbs,\n",
"GRB.MAXIMIZE)\n",
"\n",
"# add sum-to-half constraints to both models\n",
"model_rbs.addConstr(np.ones([1,2]) @ w_std == 1, name=\"sum-to-one\")\n",
"# add quadratic uncertainty constraint to the robust model\n",
"model_rbs.addConstr(z_rbs**2 <= np.inner(w_rbs, omega@w_rbs), name=\"uncertainty\")\n",
"model_rbs.addConstr(z_rbs >= 0, name=\"z positive\")\n",
"model_rbs.addConstr(z_rbs**2 >= w_rbs.transpose()@omega@w_rbs, name=\"uncertainty\")\n",
"\n",
"# solve problem with Gurobi, print alongside standard solution for comparison\n",
"model_rbs.optimize()\n",
Expand Down
37 changes: 23 additions & 14 deletions Scratch/numerical-error.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,16 @@
"source": [
"# Numerical Error\n",
"\n",
"The purpose of this notebook is to illustrate the numerical error which comes from setting $\\kappa = 0$."
"The purpose of this notebook is to illustrate the numerical error which comes from setting $\\kappa = 0$ and using \n",
"\n",
"```python3\n",
"model.addConstr(z**2 <= np.inner(w, omega@w), name=\"uncertainty\")\n",
"```\n",
"as the cone constraint rather than\n",
"```python3\n",
"model.addConstr(z**2 <= w.transpose()@omega@w, name=\"uncertainty\")\n",
"```\n",
"This has been fixed in the latest version of all other notebooks."
]
},
{
Expand Down Expand Up @@ -162,7 +171,7 @@
"model = gp.Model(\"standardGS\")\n",
"\n",
"# define variable of interest as a continuous \n",
"w = model.addMVar(shape=problem_size, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"w = model.addMVar(shape=problem_size, lb=0.0, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"\n",
"# set the objective function\n",
"model.setObjective(\n",
Expand All @@ -172,8 +181,8 @@
"# add sub-to-half constraints\n",
"model.addConstr(M @ w == m, name=\"sum-to-half\")\n",
"# add weight-bound constraints\n",
"model.addConstr(w >= lower_bound, name=\"lower bound\")\n",
"model.addConstr(w <= upper_bound, name=\"upper bound\")\n",
"model.addConstr(w >= lower_bound, name=\"lower-bound\")\n",
"model.addConstr(w <= upper_bound, name=\"upper-bound\")\n",
"\n",
"# solve the problem with Gurobi\n",
"model.optimize()\n",
Expand Down Expand Up @@ -231,8 +240,8 @@
"CPU model: Intel(R) Core(TM) i5-8350U CPU @ 1.70GHz, instruction set [SSE2|AVX|AVX2]\n",
"Thread count: 4 physical cores, 8 logical processors, using up to 8 threads\n",
"\n",
"Optimize a model with 23 rows, 4 columns and 25 nonzeros\n",
"Model fingerprint: 0x97227d81\n",
"Optimize a model with 22 rows, 4 columns and 24 nonzeros\n",
"Model fingerprint: 0x52c6cd50\n",
"Model has 3 quadratic objective terms\n",
"Model has 3 quadratic constraints\n",
"Coefficient statistics:\n",
Expand All @@ -242,7 +251,7 @@
" QObjective range [1e+00, 5e+00]\n",
" Bounds range [0e+00, 0e+00]\n",
" RHS range [5e-01, 1e+00]\n",
"Presolve removed 22 rows and 1 columns\n",
"Presolve removed 21 rows and 1 columns\n",
"Presolve time: 0.01s\n",
"Presolved: 9 rows, 8 columns, 14 nonzeros\n",
"Presolved model has 3 second-order cone constraints\n",
Expand All @@ -263,7 +272,7 @@
" 4 -8.12263522e-01 -8.12579837e-01 1.20e-14 2.22e-16 2.64e-05 0s\n",
" 5 -8.12499544e-01 -8.12500834e-01 8.52e-13 3.00e-15 1.07e-07 0s\n",
"\n",
"Barrier solved model in 5 iterations and 0.02 seconds (0.00 work units)\n",
"Barrier solved model in 5 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective -8.12499544e-01\n",
"\n",
"Standard: w = [0.50000 0.37500 0.12500],\n",
Expand All @@ -277,8 +286,8 @@
"model = gp.Model(\"robustGS\")\n",
"\n",
"# define variables of interest as a continuous\n",
"w = model.addMVar(shape=problem_size, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z = model.addVar(name=\"z\")\n",
"w = model.addMVar(shape=problem_size, lb=0.0, vtype=GRB.CONTINUOUS, name=\"w\")\n",
"z = model.addVar(lb=0.0, name=\"z\")\n",
"\n",
"# setup the robust objective function\n",
"model.setObjective(\n",
Expand All @@ -287,12 +296,12 @@
"\n",
"# add quadratic uncertainty constraint\n",
"model.addConstr(z**2 <= np.inner(w, omega@w), name=\"uncertainty\")\n",
"model.addConstr(z >= 0, name=\"z positive\")\n",
"# model.addConstr(z**2 <= w.transpose()@omega@w, name=\"uncertainty\")\n",
"# add sub-to-half constraints\n",
"model.addConstr(M @ w == m, name=\"sum-to-half\")\n",
"# add weight-bound constraints~\n",
"model.addConstr(w >= lower_bound, name=\"lower bound\")\n",
"model.addConstr(w <= upper_bound, name=\"upper bound\")\n",
"model.addConstr(w >= lower_bound, name=\"lower-bound\")\n",
"model.addConstr(w <= upper_bound, name=\"upper-bound\")\n",
"\n",
"# solve the problem with Gurobi\n",
"model.optimize()\n",
Expand All @@ -303,7 +312,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that even though $\\kappa = 0$, we have $z\\neq0$. This suggests that a small amount of numerical error has been introduced as a result of asking Gurobi to optimize $z$ as well as $w$. "
"Note that even though $\\kappa = 0$, we have $z\\neq0$. This suggests that a small amount of numerical error has been introduced."
]
}
],
Expand Down
Loading

0 comments on commit f21e77f

Please sign in to comment.