Merge pull request #6 from fabian-sp/f-refactor

Refactor

README.md

@@ -7,9 +7,11 @@ This repository contains a `Python` implementation of the SQP-GS (*Sequential Qu
 
 The algorithm can solve problems of the form
 
+```
     min  f(x)
     s.t. g(x) <= 0
          h(x) = 0
+```
 
 where `f`, `g` and `h` are locally Lipschitz functions. Hence, the algorithm can solve problems with nonconvex and nonsmooth objective and constraints. For details, we refer to the original paper.
 
@@ -19,14 +21,17 @@ The code was tested for a 2-dim nonsmooth version of the Rosenbrock function, co
 
 To reproduce this experiment, see the file `example_rosenbrock.py`.
 
-![SQP-GS trajectories for a 2-dim example](rosenbrock.png "SQP-GS trajectories for a 2-dim example")
+![SQP-GS trajectories for a 2-dim example](data/img/rosenbrock.png "SQP-GS trajectories for a 2-dim example")
 
 
 ## Implementation details
 The solver can be called via 
 
-    SQP_GS(f, gI, gE)
-
+```python
+    from ncopt.sqpgs import SQPGS
+    problem = SQPGS(f, gI, gE)
+    problem.solve()
+```
 It has three main arguments, called `f`, `gI` and `gE`. `f` is the objective. `gI` and `gE` are lists of inequality and equality constraint functions. Each element of `gI` and `gE` as well as the objective `f` needs to be an instance of a class which contains the following properties. The constraint functions are allowed to have multi-dimensional output.
 
 ### Attributes
@@ -45,6 +50,7 @@ Moreover, we implemented a class for a constraint coming from a Pytorch neural n
 
 
 ## References
-* [1] F. E. Curtis and M. L. Overton, A sequential quadratic programming algorithm for nonconvex, nonsmooth constrained optimization, SIAM Journal on Optimization, 22 (2012), pp. 474–500, https://doi.org/10.1137/090780201.
+[1] Frank E. Curtis and Michael L. Overton, A sequential quadratic programming algorithm for nonconvex, nonsmooth constrained optimization, 
+SIAM Journal on Optimization 2012 22:2, 474-500, https://doi.org/10.1137/090780201.
 
-* [2] F. E. Curtis, MATLAB implementation of SLQP-GS, https://coral.ise.lehigh.edu/frankecurtis/software/.
+[2] Frank E. Curtis and Michael L. Overton, MATLAB implementation of SLQP-GS, https://coral.ise.lehigh.edu/frankecurtis/software/.

data/checkpoints/README.md

@@ -0,0 +1 @@
+# Folder for storing Pytorch checkpoints

data/img/README.md

@@ -0,0 +1 @@
+# Folder for storing images

example_rosenbrock.py

@@ -1,52 +1,84 @@
 """
-author: Fabian Schaipp
+@author: Fabian Schaipp
 """
 import numpy as np
 import matplotlib.pyplot as plt
 
-from ncopt.sqpgs import SQP_GS
+from ncopt.sqpgs import SQPGS
 from ncopt.funs import f_rosenbrock, g_max, g_linear
-#from ncopt.torch_obj import Net
 
-#%%
+#%% Setup
+
 f = f_rosenbrock()
 g = g_max()
 
-A = np.eye(2); b = np.ones(2)*5; g1 = g_linear(A, b)
-#D = Net(model)
+# could add this equality constraint
+# A = np.eye(2); b = np.ones(2)*5; g1 = g_linear(A, b)
 
 # inequality constraints (list of functions)
 gI = [g]
-# equality constraints (list of scalar functions)
+# equality constraints (list of functions)
 gE = []
 
-xstar = np.array([1/np.sqrt(2), 0.5])
+xstar = np.array([1/np.sqrt(2), 0.5])       # solution
 
-#%%
-X, Y = np.meshgrid(np.linspace(-2,2,100), np.linspace(-2,2,100))
+#%% How to use the solver
+
+problem = SQPGS(f, gI, gE, x0=None, tol=1e-6, max_iter=100, verbose=True)
+x = problem.solve()
+
+print(f"Distance to solution:", f'{np.linalg.norm(x-xstar):.9f}')
+
+#%% Solve from multiple starting points and plot
+
+_x, _y = np.linspace(-2,2,100), np.linspace(-2,2,100)
+X, Y = np.meshgrid(_x, _y)
 Z = np.zeros_like(X)
 
-for j in np.arange(100):
-    for i in np.arange(100):
+for j in np.arange(X.shape[0]):
+    for i in np.arange(X.shape[1]):
         Z[i,j] = f.eval(np.array([X[i,j], Y[i,j]]))
 
+#%%
+from matplotlib.lines import Line2D
+
+fig, ax = plt.subplots(figsize=(5,4))
 
-fig, ax = plt.subplots()
-ax.contourf(X,Y,Z, levels = 20)
-ax.scatter(xstar[0], xstar[1], marker = "*", s = 200, c = "gold", alpha = 1, zorder = 200)   
+np.random.seed(123)
 
+# Plot contour and solution
+ax.contourf(X, Y, Z, cmap='gist_heat', levels=20, alpha=0.7, 
+            antialiased=True, lw=0, zorder=0)
+# ax.contourf(X, Y, Z, colors='lightgrey', levels=20, alpha=0.7, 
+#             antialiased=True, lw=0, zorder=0)
+# ax.contour(X, Y, Z, cmap='gist_heat', levels=8, alpha=0.7, 
+#             antialiased=True, linewidths=4, zorder=0)
+ax.scatter(xstar[0], xstar[1], marker="*", s=200, c="gold", zorder=1)   
 
 for i in range(20):
-    x0 = np.random.randn(2)# np.zeros(2)
-    x_k, x_hist, SP = SQP_GS(f, gI, gE, x0, tol = 1e-6, max_iter = 100, verbose = False)
+    x0 = np.random.randn(2)
+    problem = SQPGS(f, gI, gE, x0, tol=1e-6, max_iter=100, verbose=False,
+                    store_history=True)
+    x_k = problem.solve()
     print(x_k)
-    ax.plot(x_hist[:,0], x_hist[:,1], c = "silver", lw = 0.7, ls = '--', alpha = 0.5)
-    ax.scatter(x_k[0], x_k[1], marker = "+", s = 50, c = "k", alpha = 1, zorder = 210)
+
+    x_hist = problem.x_hist
+    ax.plot(x_hist[:,0], x_hist[:,1], c="silver", lw=1, ls='--', alpha=0.5, zorder=2)
+    ax.scatter(x_k[0], x_k[1], marker="+", s=50, c="k", zorder=3)
+    ax.scatter(x0[0], x0[1], marker="o", s=30, c="silver", zorder=3)
 
 ax.set_xlim(-2,2)
 ax.set_ylim(-2,2)
 
+fig.suptitle("Trajectory for multiple starting points")
 
+legend_elements = [Line2D([0], [0], marker='*', lw=0, color='gold', label='Solution',
+                          markersize=15),
+                   Line2D([0], [0], marker='o', lw=0, color='silver', label='Starting point',
+                          markersize=8),
+                   Line2D([0], [0], marker='+', lw=0, color='k', label='Final iterate',
+                          markersize=8),]
+ax.legend(handles=legend_elements, ncol=3, fontsize=8)
 
-
-
+fig.tight_layout()
+fig.savefig('data/img/rosenbrock.png')