From 4e7e5e83bfa8b49a07af39f772704f3fce2d5b01 Mon Sep 17 00:00:00 2001 From: Enzo Busseti Date: Sun, 9 Jun 2024 14:45:55 +0400 Subject: [PATCH] some comments in readme --- README.md | 3 - README.rst | 96 ++++++++++++++++++++ project_euromir/__init__.py | 9 ++ project_euromir/cvxpy_solver.py | 65 +++++++++++++ project_euromir/linear_algebra.py | 57 ++++++++++++ project_euromir/loss.py | 29 ++++++ project_euromir/tests/test_linear_algebra.py | 3 + 7 files changed, 259 insertions(+), 3 deletions(-) delete mode 100644 README.md create mode 100644 README.rst create mode 100644 project_euromir/cvxpy_solver.py create mode 100644 project_euromir/linear_algebra.py create mode 100644 project_euromir/loss.py diff --git a/README.md b/README.md deleted file mode 100644 index 9dff7a4..0000000 --- a/README.md +++ /dev/null @@ -1,3 +0,0 @@ -# project_euromir - -[https://open.spotify.com/track/3ffkbz5OvPjXjOsYTsEjKu](https://open.spotify.com/track/3ffkbz5OvPjXjOsYTsEjKu) diff --git a/README.rst b/README.rst new file mode 100644 index 0000000..e528eac --- /dev/null +++ b/README.rst @@ -0,0 +1,96 @@ +Project EuroMir +=============== + +Project codename `EuroMir `_ +`(also, song) `_ +is a prototype of a conic programming solver that builds on the work done for +the older `conic program refinement +`_ research project, of which it +may inherit the name, once it's completed. + +Compared to that 2018 Stanford research project, it has a completely new +codebase (written from scratch) and it removes various unnecessary +dependencies. It also has a modified algorithm, which is guaranteed to preserve +strong convexity (unlike many similar attempts). It uses a simplified version +of the 2018 algorithm only for the final polishing. + +Algorithm (draft) +================= + +The algorithm is under development. This is the current model, see the +`scs `_ and +`conic refinement +`_ papers to +understand the notation: + +.. math:: + + \begin{array}{ll} + + \text{minimize} & \|Q u - v \|_2^2 + \| u - \Pi u \|_2^2 + \| v - \Pi^\star v \|_2^2 + + \end{array} + +The system matrix :math:`Q` from the homogeneous self-dual embedding is skew +symmetric, so at convergence it is guaranteed that :math:`u` and :math:`v` are +orthogonal. The objective function is clearly convex and has continuous +derivative. If we drop the zero cone variables (as we do), the projections +onto :math:`\mathbf{R}`, and any all-zero rows and columns of :math:`Q`, it is +also strongly convex. The conditioning and Lipschitz constant are all dependent +on the conditioning of :math:`Q`, we apply by default standard `Ruiz diagonal +pre-conditioning `_. + +We use `limited-memory BFGS +`_, as it's implemented in +the `variable-metric limited memory module of MINPACK-2 +`_, for the minimization routine. + +We then use the 2018 conic refinement algorithm (simplified, without the +normalization step), for the final refinement of the solution obtained from the +BFGS loop. + +The approach is clearly globally convergent, and it can probably be showed to +converge super-linearly. However, theoretical bounds on convergence have very +limited applicability in real mathematical engineering; smart implementation +choices and good code quality are by far more important. Only experimentation +can tell how this approach compares to mature alternative ones, such as +the various flavors of interior point and operator splitting algorithms. + + +Development +=========== + +We mirror development in Python and C (and, in the future, C with OpenMP). We +set up the build with CMake, patching ``setuptools`` to use it. All testing and +profiling is done in Python. Binding of C to Python is done using ``ctypes``: +we link at runtime. This guarantees that pre-built wheels are agnostic to the +Python version, the Numpy ABI, .... +Memory communication between Python (mostly Numpy arrays) and C is still +zero-copy and the GIL is released by ``ctypes``; there are no appreciable +overheads compared to building with ``Python.h``, since the user only calls a +``ctypes``-linked function once when solving a conic program. + +Memory +====== + +All memory is pre-allocated by the caller, which in typical usage is the CVXPY +interface. There are no ``malloc`` in our C code. We require to allocate space +for a copy of the problem data (for rescaling, that may be prevented if we +change the internal logic, adding some computational overhead), and the size of +the primal and dual variables times the L-BFGS memory, which is a small number +like 5 or 10. + +In layman terms, this is the least memory any conic solver needs, and +dramatically less than interior-point solvers. + +Usage +===== + +We will provide a CVXPY and raw Python interface as part of our packages. The +single C function the user interacts with will be also documented, for usage +from other runtime environments. + +Installation +============ + +Pre-built wheels will be available on PyPi soon. diff --git a/project_euromir/__init__.py b/project_euromir/__init__.py index 2cbac42..412737e 100644 --- a/project_euromir/__init__.py +++ b/project_euromir/__init__.py @@ -141,3 +141,12 @@ def fun(**kwargs): ('mult', 'double'), ) ) + +### +# CVXPY interface +### + +try: + from .cvxpy_solver import Solver +except ImportError: + pass diff --git a/project_euromir/cvxpy_solver.py b/project_euromir/cvxpy_solver.py new file mode 100644 index 0000000..ac25767 --- /dev/null +++ b/project_euromir/cvxpy_solver.py @@ -0,0 +1,65 @@ +# BSD 3-Clause License + +# Copyright (c) 2024-, Enzo Busseti + +# Redistribution and use in source and binary forms, with or without +# modification, are permitted provided that the following conditions are met: + +# 1. Redistributions of source code must retain the above copyright notice, this +# list of conditions and the following disclaimer. + +# 2. Redistributions in binary form must reproduce the above copyright notice, +# this list of conditions and the following disclaimer in the documentation +# and/or other materials provided with the distribution. + +# 3. Neither the name of the copyright holder nor the names of its +# contributors may be used to endorse or promote products derived from +# this software without specific prior written permission. + +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE +# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL +# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR +# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER +# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +"""CVXPY solver interface. + +Documentation: + https://www.cvxpy.org/tutorial/solvers/index.html#custom-solvers + +Relevant base class: + https://github.com/cvxpy/cvxpy/blob/master/cvxpy/reductions/solvers/conic_solvers/conic_solver.py + +Model class: + https://github.com/cvxpy/cvxpy/blob/master/cvxpy/reductions/solvers/conic_solvers/scs_conif.py +""" + +try: + from cvxpy.reductions.solvers.conic_solvers.conic_solver import ( + ConicSolver, NonNeg, Zero) +except ImportError as exc: # pragma: no cover + raise ImportError( + "Can't use CVXPY interface if CVXPY is not installed!") from exc + + +class Solver(ConicSolver): + """CVXPY solver interface. + + We follow SCS conventions for simplicity: cones ordering, in the future + exp and sdp cone definition, names in the dims dict (SCS 3), ... + """ + + MIP_CAPABLE = False + SUPPORTED_CONSTRAINTS = [Zero, NonNeg] + REQUIRES_CONSTR = False + + def name(self): + return "PROJECT_EUROMIR" + + def solve_via_data(self, *args, **kwargs): + print("Solving with a custom QP solver!") + super().solve_via_data(*args, **kwargs) diff --git a/project_euromir/linear_algebra.py b/project_euromir/linear_algebra.py new file mode 100644 index 0000000..95de1c6 --- /dev/null +++ b/project_euromir/linear_algebra.py @@ -0,0 +1,57 @@ +# BSD 3-Clause License + +# Copyright (c) 2024-, Enzo Busseti + +# Redistribution and use in source and binary forms, with or without +# modification, are permitted provided that the following conditions are met: + +# 1. Redistributions of source code must retain the above copyright notice, this +# list of conditions and the following disclaimer. + +# 2. Redistributions in binary form must reproduce the above copyright notice, +# this list of conditions and the following disclaimer in the documentation +# and/or other materials provided with the distribution. + +# 3. Neither the name of the copyright holder nor the names of its +# contributors may be used to endorse or promote products derived from +# this software without specific prior written permission. + +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE +# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL +# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR +# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER +# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +"""Basic linear algebra operations.""" + +from __future__ import annotations + +import numpy as np +import scipy as sp + + +def Q_matvec( + matrix: sp.csc_matrix, + b: np.array, + c: np.array, + input: np.array) -> np.array: + """Python implementation of HSDE's Q matvec.""" + + m = len(b) + n = len(c) + output = np.zeros(m + n + 1) + + output[:n] += matrix.T @ input[n:n+m] + output[:n] += c * input[-1] + + output[n:n+m] -= matrix @ input[:n] + output[n:n+m] += b * input[-1] + + output[-1] -= np.dot(c, input[:n]) + output[-1] -= np.dot(b, input[n:n+m]) + + return output diff --git a/project_euromir/loss.py b/project_euromir/loss.py new file mode 100644 index 0000000..e5fabea --- /dev/null +++ b/project_euromir/loss.py @@ -0,0 +1,29 @@ +# BSD 3-Clause License + +# Copyright (c) 2024-, Enzo Busseti + +# Redistribution and use in source and binary forms, with or without +# modification, are permitted provided that the following conditions are met: + +# 1. Redistributions of source code must retain the above copyright notice, this +# list of conditions and the following disclaimer. + +# 2. Redistributions in binary form must reproduce the above copyright notice, +# this list of conditions and the following disclaimer in the documentation +# and/or other materials provided with the distribution. + +# 3. Neither the name of the copyright holder nor the names of its +# contributors may be used to endorse or promote products derived from +# this software without specific prior written permission. + +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE +# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL +# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR +# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER +# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +"""Define loss function and gradient.""" diff --git a/project_euromir/tests/test_linear_algebra.py b/project_euromir/tests/test_linear_algebra.py index 925559d..43c26e5 100644 --- a/project_euromir/tests/test_linear_algebra.py +++ b/project_euromir/tests/test_linear_algebra.py @@ -34,6 +34,9 @@ import scipy.sparse as sp import project_euromir as lib +from project_euromir.linear_algebra import Q_matvec + +from .test_equilibrate import _make_Q class TestLinearAlgebra(TestCase):