added LP, PO examples from the paper

examples/LP_cone.py

@@ -0,0 +1,23 @@
+import cvxpy as cp
+import numpy as np
+np.random.seed(0)
+
+from newton_admm import newton_admm, problems
+from plot import plotter
+
+# Benchmark LP
+name = "LP"
+m,n = 300,600
+beta, prob, data = problems.lp(m,n)
+prob.solve(solver="ECOS", abstol=1e-14, reltol=1e-14, feastol=1e-14, 
+           verbose=True, max_iters=50)
+betastar = np.array(beta.value).flatten()
+fstar = prob.value
+
+beta_from_x = data['beta_from_x']
+newton_out = newton_admm(data, data['dims'], benchmark=(beta_from_x, betastar), 
+                         maxiters=400, verbose=1)
+scs_out = newton_admm(data, data['dims'], benchmark=(beta_from_x, betastar), 
+                         admm_maxiters=20000, maxiters=0, verbose=1) 
+
+plotter(newton_out, scs_out, fstar, name, xmax=49, ymin=1e-12, ymax=1e2)

examples/PO_cone.py

@@ -0,0 +1,29 @@
+import cvxpy as cp
+import numpy as np
+np.random.seed(0)
+
+from newton_admm import newton_admm, problems
+from plot import plotter
+
+# Benchmark LP
+name = "PO"
+p = 2500
+beta, prob, data = problems.portfolio_opt(p)
+beta_from_x = data['beta_from_x']
+# ECOS cannot solve this to sufficient precision!
+# prob.solve(solver="ECOS", abstol=1e-14, reltol=1e-14, feastol=1e-14, 
+#            verbose=True)
+# prob.solve(solver="SCS", eps=1e-14)
+# betastar = np.array(beta.value).flatten()
+# fstar = prob.value
+
+baseline_out = newton_admm(data, data['dims'], maxiters=100, res_tol = 1e-14, verbose=10)
+betastar = beta_from_x(baseline_out['x'])
+fstar = baseline_out['info']['fstar']
+
+newton_out = newton_admm(data, data['dims'], benchmark=(beta_from_x, betastar), 
+                         maxiters=100, ridge=1e-4, verbose=1)
+scs_out = newton_admm(data, data['dims'], benchmark=(beta_from_x, betastar), 
+                         admm_maxiters=2000, maxiters=0, verbose=100) 
+
+plotter(newton_out, scs_out, fstar, name, xmax=100, ymin=1e-12, ymax=1e2) 

examples/plot.py

@@ -0,0 +1,52 @@
+
+import numpy as np
+
+import matplotlib
+matplotlib.use('svg')
+matplotlib.rcParams['pdf.fonttype'] = 42
+matplotlib.rcParams['ps.fonttype'] = 42
+import matplotlib.pyplot as plt
+
+
+def plotter(newton_out, scs_out, fstar, name, ours="Newton-ADMM", 
+            xmax=None, ymin=1e-12, ymax=1e-2, lw=3, fontsize=14, 
+            loc= "upper right"): 
+    newton_out = newton_out['info']['benchmark']
+    scs_out = scs_out['info']['benchmark']
+
+    solve_times_SCS_ret = np.cumsum(scs_out['time'])
+    solve_times_newt_ret = np.cumsum(newton_out['time'])
+    errors_SCS_ret = scs_out['error']
+    errors_newt_ret = newton_out['error']
+    subopts_SCS = np.abs(fstar - np.array(scs_out['fval']))
+    subopts_newt = np.abs(fstar - np.array(newton_out['fval']))
+
+    # Make suboptimality plot.
+    fig,ax = plt.subplots()
+    ax.semilogy(solve_times_SCS_ret, subopts_SCS, lw=lw, color="blue", label="SCS")
+    ax.semilogy(solve_times_newt_ret, subopts_newt, lw=lw, color="crimson", label=ours)
+
+    ax.set_xlabel("Seconds", fontsize=fontsize)
+    ax.set_ylabel("Suboptimality", fontsize=fontsize)
+    plt.tick_params(axis='both', which='major', labelsize=fontsize)
+    ax.set_xlim([0, xmax])
+    ax.set_ylim([ymin, ymax])
+    ax.legend(loc=loc, fontsize=fontsize)
+    fig.savefig("subopts_" + name + ".pdf", bbox_inches="tight")
+
+    # Make solution error plot.
+    lw = 3
+    fontsize = 14
+    loc = "upper right"
+
+    fig,ax = plt.subplots()
+    ax.semilogy(solve_times_SCS_ret, errors_SCS_ret, lw=lw, color="blue", label="SCS")
+    ax.semilogy(solve_times_newt_ret, errors_newt_ret, lw=lw, color="crimson", label=ours)
+
+    ax.set_xlabel("Seconds", fontsize=fontsize)
+    ax.set_ylabel("Estimation error", fontsize=fontsize)
+    plt.tick_params(axis='both', which='major', labelsize=fontsize)
+    ax.set_xlim([0, xmax])
+    ax.set_ylim([ymin, ymax])
+    ax.legend(loc=loc, fontsize=fontsize)
+    fig.savefig("estims_" + name + ".pdf", bbox_inches="tight")