examples/notebooks/Solver Tolerances.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:a290494181bd86fa07d36d9005b0a896b6b9bddde6917f1aab828fbd087a5aa9"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Numerical Tolerances in CVXPY Solvers"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The problem data below is taken from  the example problem: \"Computing a sparse solution of a set of linear inequalities,\" adapted by Judson Wilson on 5/11/2014 from the origional work by Almir Mutapcic 02/28/06.  This document was created 6/4/2014 by David Stonestrom.\n",
      "\n",
      "<br> <br>\n",
      "\n",
      "The purpose of this document is to demonstrate differences in the numerical tolerances between the solvers CVXPY uses.  Each solver is applied to the same problem, and the results are presented in a table below.  The specific problem is to apply the $l1$ heuristic to find a sparse feasible point within a set of linear inequalities.  Code implementing this problem is included below the results.  \n",
      "\n",
      "<br>\n",
      "$$ \\textbf{minimize } f\\left(x\\right) = \\left|\\left|x\\right|\\right|_1$$\n",
      "\n",
      "$$\\text{subject to } Ax \\preceq b$$\n",
      "<br>\n",
      "\n",
      "\n",
      "\n",
      "\n",
      "\n"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Note that while the three objective values are very close, the number of nonzero elements at each tolerance is different for each solver.  Below is a table summarizing the results.  Fewer nonzero elements is better.  \n",
      "\n",
      "<br>\n",
      "\\begin{array}{lllll}\n",
      "\\text{Solver} & \\text{Objective} & \\text{E-8 tolerance nonzeros} & \\text{E-4 tolerance nonzeros} & \\text{Max constraint violation} \\\\ \n",
      "\\text{ECOS}   & 28.58239409 & 40 & 40 & 1.5\\text{e-}9 \\\\ \n",
      "\\text{CVXOPT} & 28.58239423 & 48 & 40 & -8.2\\text{e-}10 \\\\ \n",
      "\\text{SCS}    & 28.5795     & 50 & 48 & 4.8\\text{e-}3 \\\\\n",
      "\\end{array}\n"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import cvxpy as cvx\n",
      "import numpy\n",
      "\n",
      "# Fix random number generator so we can repeat the experiment.\n",
      "numpy.random.seed(1) #Note seed 0 leads problem data that cannot be solved.\n",
      "\n",
      "# The threshold value below which we consider an element to be zero.\n",
      "delta = 1e-8\n",
      "delta_loose = 1e-4\n",
      "# Problem dimensions (m inequalities in n-dimensional space).\n",
      "m = 100\n",
      "n = 50\n",
      "\n",
      "# Construct a feasible set of inequalities.\n",
      "# (This system is feasible for the x0 point.)\n",
      "A  = numpy.matrix(numpy.random.randn(m,n))\n",
      "x0 = numpy.matrix(numpy.random.randn(n,1))\n",
      "b  = A*x0 + numpy.random.random((m,1))\n",
      "# Create variable.\n",
      "x = cvx.Variable(n)\n",
      "\n",
      "# Create constraint.\n",
      "constraints = [A*x <= b]\n",
      "\n",
      "# Form objective.\n",
      "obj = cvx.Minimize(cvx.norm(x, 1))\n",
      "\n",
      "# Form and solve problem.\n",
      "prob = cvx.Problem(obj, constraints)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "heading",
     "level": 3,
     "metadata": {},
     "source": [
      "ECOS"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "try: \n",
      "    print prob.solve(solver=cvx.ECOS)\n",
      "    print \"status:\", prob.status\n",
      "except Exception as e:\n",
      "    print e\n",
      "\n",
      "#Number of nonzero elements in the solution (its cardinality or diversity)\n",
      "nnz_l1 = (numpy.absolute(x.value) > delta).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-8 nonzeros.'.format(n, nnz_l1)\n",
      "nnz_l1_loose = (numpy.absolute(x.value) > delta_loose).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-4 nonzeros.'.format(n, nnz_l1_loose)\n",
      "print \"Max constraint violation: \", max(A*x.value - b)[0,0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "28.5823940884\n",
        "status: optimal\n",
        "Found a feasible x in R^50 that has 40 e-8 nonzeros.\n",
        "Found a feasible x in R^50 that has 40 e-4 nonzeros.\n",
        "Max constraint violation:  1.49974610508e-09\n"
       ]
      }
     ],
     "prompt_number": 2
    },
    {
     "cell_type": "heading",
     "level": 3,
     "metadata": {},
     "source": [
      "CVXOPT"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "try: \n",
      "    print prob.solve(solver=cvx.CVXOPT)\n",
      "    print \"status:\", prob.status\n",
      "except Exception as e:\n",
      "    print e\n",
      "\n",
      "#Number of nonzero elements in the solution (its cardinality or diversity)\n",
      "nnz_l1 = (numpy.absolute(x.value) > delta).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-8 nonzeros.'.format(n, nnz_l1)\n",
      "nnz_l1_loose = (numpy.absolute(x.value) > delta_loose).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-4 nonzeros.'.format(n, nnz_l1_loose)\n",
      "print \"Max constraint violation: \", max(A*x.value - b)[0,0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "28.5823942271\n",
        "status: optimal\n",
        "Found a feasible x in R^50 that has 48 e-8 nonzeros.\n",
        "Found a feasible x in R^50 that has 40 e-4 nonzeros.\n",
        "Max constraint violation:  -8.15034262303e-10\n"
       ]
      }
     ],
     "prompt_number": 3
    },
    {
     "cell_type": "heading",
     "level": 3,
     "metadata": {},
     "source": [
      "SCS"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "try: \n",
      "    print prob.solve(solver=cvx.SCS)\n",
      "    print \"status:\", prob.status\n",
      "except Exception as e:\n",
      "    print e\n",
      "\n",
      "#Number of nonzero elements in the solution (its cardinality or diversity)\n",
      "nnz_l1 = (numpy.absolute(x.value) > delta).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-8 nonzeros.'.format(n, nnz_l1)\n",
      "nnz_l1_loose = (numpy.absolute(x.value) > delta_loose).sum()\n",
      "print 'Found a feasible x in R^{} that has {} e-4 nonzeros.'.format(n, nnz_l1_loose)\n",
      "print \"Max constraint violation: \", max(A*x.value - b)[0,0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "28.5794926837\n",
        "status: optimal\n",
        "Found a feasible x in R^50 that has 50 e-8 nonzeros.\n",
        "Found a feasible x in R^50 that has 48 e-4 nonzeros.\n",
        "Max constraint violation:  0.00477885544772\n"
       ]
      }
     ],
     "prompt_number": 4
    }
   ],
   "metadata": {}
  }
 ]
}