Skip to content

Commit

Permalink
Updated Julian's $n = 2$ example
Browse files Browse the repository at this point in the history
  • Loading branch information
@Foggalong
Foggalong committed Apr 23, 2024
1 parent a44afc2 commit a894426
Showing 1 changed file with 76 additions and 69 deletions.
145 changes: 76 additions & 69 deletions Example2/n2-test.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 21,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -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 All @@ -45,7 +45,7 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 22,
"metadata": {},
"outputs": [],
"source": [
Expand Down Expand Up @@ -97,27 +97,25 @@
},
{
"cell_type": "code",
"execution_count": 3,
"execution_count": 23,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Set parameter Username\n",
"Academic license - for non-commercial use only - expires 2025-02-26\n",
"Gurobi Optimizer version 11.0.0 build v11.0.0rc2 (linux64 - \"Ubuntu 22.04.4 LTS\")\n",
"\n",
"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 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 +131,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",
"Barrier solved model in 7 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective 1.75000000e+00\n",
"\n",
"Standard: [0.25000 0.75000]\n"
"Standard: [0.00000001 0.99999999]\n"
]
}
],
Expand All @@ -154,7 +152,8 @@
" \"A02.txt\",\n",
" \"EBV02.txt\",\n",
" \"S02.txt\",\n",
" 2)\n",
" 2\n",
")\n",
"\n",
"lam = 0.5\n",
"\n",
Expand All @@ -164,8 +163,9 @@
"\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",
" # 0.5*w_std@(sigma@w_std) - lam*w_std.transpose()@mubar,\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 All @@ -184,7 +184,7 @@
},
{
"cell_type": "code",
"execution_count": 4,
"execution_count": 24,
"metadata": {},
"outputs": [
{
Expand All @@ -197,43 +197,42 @@
"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",
"Model fingerprint: 0x41322026\n",
"Model has 2 quadratic objective terms\n",
"Model has 2 quadratic constraints\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: 7 rows, 8 columns, 12 nonzeros\n",
"Presolved model has 3 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.700e+01\n",
" Factor NZ : 2.800e+01\n",
" Factor Ops : 1.400e+02 (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 -2.28881626e-02 1.22656250e+00 9.39e-01 1.43e-01 3.26e-01 0s\n",
" 1 6.91883948e-01 1.06699842e+00 1.98e-01 1.57e-07 5.84e-02 0s\n",
" 2 6.44324078e-01 7.21385315e-01 2.17e-07 1.73e-13 7.71e-03 0s\n",
" 3 6.67762956e-01 6.69433547e-01 6.21e-09 6.66e-16 1.67e-04 0s\n",
" 4 6.68355745e-01 6.68371340e-01 1.53e-13 2.00e-15 1.56e-06 0s\n",
" 5 6.68366724e-01 6.68367382e-01 3.57e-11 6.05e-13 6.58e-08 0s\n",
"\n",
"Barrier solved model in 6 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective -6.73611064e-01\n",
"Barrier solved model in 5 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective 6.68366724e-01\n",
"\n",
"Standard: [0.25000 0.75000]\n",
"Robust: [0.41660 0.58340]\n"
"Standard: [0.00000001 0.99999999]\n",
"Robust: [0.85714285 0.14285715]\n"
]
}
],
Expand All @@ -249,13 +248,14 @@
"z_rbs = model_rbs.addVar(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",
" # 0.5*w_rbs.transpose()@(sigma@w_rbs) - lam*w_rbs.transpose()@mubar - kappa*z_rbs,\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**2 >= np.inner(w_rbs, omega@w_rbs), name=\"uncertainty\")\n",
"model_rbs.addConstr(z_rbs >= 0, name=\"z positive\")\n",
"\n",
"# solve problem with Gurobi, print alongside standard solution for comparison\n",
Expand All @@ -268,12 +268,12 @@
"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 eight decimal places."
]
},
{
"cell_type": "code",
"execution_count": 5,
"execution_count": 25,
"metadata": {},
"outputs": [
{
Expand All @@ -286,49 +286,49 @@
"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",
"Model fingerprint: 0xe31c34b4\n",
"Model has 2 quadratic objective terms\n",
"Model has 2 quadratic constraints\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: 7 rows, 8 columns, 12 nonzeros\n",
"Presolved model has 3 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.700e+01\n",
" Factor NZ : 2.800e+01\n",
" Factor Ops : 1.400e+02 (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.04e+00 6.81e-01 5.08e-01 0s\n",
" 1 1.54200502e+00 1.52822126e+00 4.89e-02 1.42e-01 7.24e-02 0s\n",
" 2 1.64909037e+00 1.63852570e+00 5.38e-08 3.21e-02 1.66e-02 0s\n",
" 3 1.74359514e+00 1.75057927e+00 3.23e-10 8.27e-05 7.48e-04 0s\n",
" 4 1.74999141e+00 1.75002072e+00 1.01e-12 9.47e-08 2.99e-06 0s\n",
" 5 1.74999468e+00 1.75000778e+00 3.11e-12 3.00e-08 1.33e-06 0s\n",
" 6 1.74999882e+00 1.75000092e+00 1.06e-11 1.95e-09 2.11e-07 0s\n",
"\n",
"Barrier solved model in 6 iterations and 0.03 seconds (0.00 work units)\n",
"Optimal objective -5.62499955e-01\n",
"Optimal objective 1.74999882e+00\n",
"\n",
"Standard: [0.2500000031 0.7499999969]\n",
"Robust: [0.2499847190 0.7500152810]\n"
"Standard: [0.00000001 0.99999999]\n",
"Robust: [0.00000002 0.99999998]\n"
]
}
],
"source": [
"# other than kappa, reuse problem variables from the last cell\n",
"kappa = 0.0\n",
"kappa = 0\n",
"\n",
"# initialise robust genetic selection model\n",
"model_rbs = gp.Model(\"n02robust_k0\")\n",
Expand All @@ -338,21 +338,28 @@
"z_rbs = model_rbs.addVar(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",
" # 0.5*w_rbs.transpose()@(sigma@w_rbs) - lam*w_rbs.transpose()@mubar - kappa*z_rbs,\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**2 >= np.inner(w_rbs, omega@w_rbs), name=\"uncertainty\")\n",
"model_rbs.addConstr(z_rbs >= 0, name=\"z positive\")\n",
"\n",
"# solve problem with Gurobi, print alongside standard solution for comparison\n",
"model_rbs.optimize()\n",
"np.set_printoptions(formatter={'float_kind':\"{:.10f}\".format})\n",
"print(f\"Standard: {w_std.X}\")\n",
"print(f\"Robust: {w_rbs.X}\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
Expand Down

0 comments on commit a894426

Please sign in to comment.