{ "cells": [ { "cell_type": "markdown", "id": "5257bb27", "metadata": {}, "source": [ "# Dictionary Learning\n", "\n", "Solve orthogonal dictionary learning problem taken from: Yu Bai, Qijia Jiang, and Ju Sun. \n", "[\"Subgradient descent learns orthogonal dictionaries.\"](https://arxiv.org/abs/1810.10702) arXiv preprint arXiv:1810.10702 (2018)." ] }, { "cell_type": "markdown", "id": "d5b1ab10", "metadata": {}, "source": [ "## Problem Description" ] }, { "cell_type": "markdown", "id": "28993e76", "metadata": {}, "source": [ "Given data $\\{y_i \\}_{i \\in[m]}$ generated as $y_i = A x_i$, where $A \\in R^{n \\times n}$ is a fixed unknown orthogonal matrix and each $x_i \\in R^n$ is an iid Bernoulli-Gaussian random vector with parameter $\\theta \\in (0,1)$, recover $A$. \n", "\n", "Write $Y \\doteq [y_1,...,y_m]$ and $X \\doteq [x_1,...,x_m]$. To find the column of $A$, one can perform the following optimization:\n", "\n", "$$\\min_{q \\in R^n} f(q) \\doteq \\frac{1}{m} ||q^T Y||_{1} = \\frac{1}{m} \\sum_{i=1}^m |q^T y_i|,$$\n", "$$\\text{s.t.} ||q||_2 = 1$$\n", "\n", "This problem is nonconvex due to the constraint and nonsmooth due to the objective.\n", "\n", "Based on the above statistical model, $q^T Y = q^T A X$ has the highest sparsity when $q$ is a column of $A$ (up to sign) so that $q^T A$ is 1-sparse. " ] }, { "cell_type": "markdown", "id": "08dfdd50", "metadata": {}, "source": [ "## Modules Importing\n", "Import all necessary modules and add PyGRANSO src folder to system path." ] }, { "cell_type": "code", "execution_count": 1, "id": "90ed32f9", "metadata": {}, "outputs": [], "source": [ "import time\n", "import numpy as np\n", "import torch\n", "import numpy.linalg as la\n", "from scipy.stats import norm\n", "from pygranso.pygranso import pygranso\n", "from pygranso.pygransoStruct import pygransoStruct" ] }, { "cell_type": "markdown", "id": "17a1b7fe", "metadata": {}, "source": [ "## Data Initialization \n", "Specify torch device, and generate data\n", "\n", "Use GPU for this problem. If no cuda device available, please set *device = torch.device('cpu')*" ] }, { "cell_type": "code", "execution_count": 2, "id": "8b4842e1", "metadata": {}, "outputs": [], "source": [ "device = torch.device('cuda')\n", "n = 30\n", "\n", "np.random.seed(1)\n", "m = 10*n**2 # sample complexity\n", "theta = 0.3 # sparsity level\n", "Y = norm.ppf(np.random.rand(n,m)) * (norm.ppf(np.random.rand(n,m)) <= theta) # Bernoulli-Gaussian model\n", "# All the user-provided data (vector/matrix/tensor) must be in torch tensor format. \n", "# As PyTorch tensor is single precision by default, one must explicitly set `dtype=torch.double`.\n", "# Also, please make sure the device of provided torch tensor is the same as opts.torch_device.\n", "Y = torch.from_numpy(Y).to(device=device, dtype=torch.double)" ] }, { "cell_type": "markdown", "id": "ec80716b", "metadata": {}, "source": [ "## Function Set-Up\n", "\n", "Encode the optimization variables, and objective and constraint functions.\n", "\n", "Note: please strictly follow the format of comb_fn, which will be used in the PyGRANSO main algortihm." ] }, { "cell_type": "code", "execution_count": 3, "id": "fb360e75", "metadata": {}, "outputs": [], "source": [ "# variables and corresponding dimensions.\n", "var_in = {\"q\": [n,1]}\n", "\n", "\n", "def user_fn(X_struct,Y):\n", " q = X_struct.q\n", " \n", " # objective function\n", " qtY = q.T @ Y\n", " f = 1/m * torch.norm(qtY, p = 1)\n", "\n", " # inequality constraint, matrix form\n", " ci = None\n", "\n", " # equality constraint \n", " ce = pygransoStruct()\n", " ce.c1 = q.T @ q - 1\n", "\n", " return [f,ci,ce]\n", "\n", "comb_fn = lambda X_struct : user_fn(X_struct,Y)" ] }, { "cell_type": "markdown", "id": "f0f55ace", "metadata": {}, "source": [ "## User Options\n", "Specify user-defined options for PyGRANSO" ] }, { "cell_type": "code", "execution_count": 4, "id": "f3a65b57", "metadata": {}, "outputs": [], "source": [ "opts = pygransoStruct()\n", "opts.torch_device = device\n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "\n", "opts.print_frequency = 10" ] }, { "cell_type": "markdown", "id": "8bca18c7", "metadata": {}, "source": [ "## Main Algorithm" ] }, { "cell_type": "code", "execution_count": 5, "id": "632976b3", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "\u001b[33m╔═════ QP SOLVER NOTICE ════════════════════════════════════════════════════════════════════════╗\n", "\u001b[0m\u001b[33m║ PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, ║\n", "\u001b[0m\u001b[33m║ the default is osqp. Users may provide their own wrapper for the QP solver. ║\n", "\u001b[0m\u001b[33m║ To disable this notice, set opts.quadprog_info_msg = False ║\n", "\u001b[0m\u001b[33m╚═══════════════════════════════════════════════════════════════════════════════════════════════╝\n", "\u001b[0m═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╗\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation ║ \n", "Version 1.2.0 ║ \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╣\n", "Problem specifications: ║ \n", " # of variables : 30 ║ \n", " # of inequality constraints : 0 ║ \n", " # of equality constraints : 1 ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " ║ <--- Penalty Function --> ║ ║ Total Violation ║ <--- Line Search ---> ║ <- Stationarity -> ║ \n", "Iter ║ Mu │ Value ║ Objective ║ Ineq │ Eq ║ SD │ Evals │ t ║ Grads │ Value ║ \n", "═════╬═══════════════════════════╬════════════════╬═════════════════╬═══════════════════════╬════════════════════╣\n", " 0 ║ 1.000000 │ 0.61751624522 ║ 0.61751624522 ║ - │ 0.000000 ║ - │ 1 │ 0.000000 ║ 1 │ 0.054664 ║ \n", " 10 ║ 1.000000 │ 0.60573380055 ║ 0.60513582468 ║ - │ 5.98e-04 ║ S │ 1 │ 1.000000 ║ 1 │ 0.024968 ║ \n", " 20 ║ 1.000000 │ 0.58456516016 ║ 0.58301955756 ║ - │ 0.001546 ║ S │ 1 │ 1.000000 ║ 1 │ 0.043517 ║ \n", " 30 ║ 1.000000 │ 0.50113197499 ║ 0.49475409554 ║ - │ 0.006378 ║ S │ 3 │ 0.250000 ║ 1 │ 0.121253 ║ \n", " 40 ║ 1.000000 │ 0.49278124194 ║ 0.49260444460 ║ - │ 1.77e-04 ║ S │ 4 │ 0.125000 ║ 1 │ 0.037304 ║ \n", " 50 ║ 1.000000 │ 0.49225009818 ║ 0.49217494722 ║ - │ 7.52e-05 ║ S │ 5 │ 0.062500 ║ 1 │ 0.032163 ║ \n", " 60 ║ 1.000000 │ 0.49212731751 ║ 0.49208854433 ║ - │ 3.88e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.051779 ║ \n", " 70 ║ 1.000000 │ 0.49203371691 ║ 0.49201049130 ║ - │ 2.32e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.054529 ║ \n", " 80 ║ 1.000000 │ 0.49197689465 ║ 0.49197679422 ║ - │ 1.00e-07 ║ S │ 2 │ 0.500000 ║ 1 │ 0.001300 ║ \n", " 90 ║ 1.000000 │ 0.49194701030 ║ 0.49194698105 ║ - │ 2.93e-08 ║ S │ 5 │ 0.062500 ║ 5 │ 6.92e-05 ║ \n", " 100 ║ 1.000000 │ 0.49194382838 ║ 0.49194381415 ║ - │ 1.42e-08 ║ S │ 6 │ 0.031250 ║ 10 │ 1.28e-05 ║ \n", " 110 ║ 1.000000 │ 0.49194277900 ║ 0.49194277111 ║ - │ 7.88e-09 ║ S │ 5 │ 0.062500 ║ 18 │ 6.93e-07 ║ \n", " 120 ║ 1.000000 │ 0.49194243076 ║ 0.49194242538 ║ - │ 5.38e-09 ║ S │ 6 │ 0.031250 ║ 27 │ 1.99e-07 ║ \n", " 130 ║ 1.000000 │ 0.49194218055 ║ 0.49194217869 ║ - │ 1.87e-09 ║ S │ 4 │ 0.125000 ║ 37 │ 5.60e-08 ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible ║ \n", "Optimization results: ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " F ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " B ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " MF ║ │ ║ 0.61751624522 ║ - │ 0.000000 ║ │ │ ║ │ ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "Iterations: 134 ║ \n", "Function evaluations: 525 ║ \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╝\n", "Total Wall Time: 3.2532527446746826s\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "start = time.time()\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "end = time.time()\n", "print(\"Total Wall Time: {}s\".format(end - start))\n", "print(max(abs(soln.final.x))) # should be close to 1" ] }, { "cell_type": "markdown", "id": "7219e881", "metadata": {}, "source": [ "## Various Options\n", "\n", "**(Optional)** Set fvalquit. Quit if the objective value drops below this value at a feasible \n", "iterate (that is, satisfying feasibility tolerances \n", "opts.viol_ineq_tol and opts.viol_eq_tol)\n", "\n", "In the example below, we get termination code 2 since the target objective reached at point feasible to tolerances" ] }, { "cell_type": "code", "execution_count": 6, "id": "d1b0c5bf", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "###### QP SOLVER NOTICE #########################################################################\n", "# PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, #\n", "# the default is osqp. Users may provide their own wrapper for the QP solver. #\n", "# To disable this notice, set opts.quadprog_info_msg = False #\n", "#################################################################################################\n", "==================================================================================================================\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation | \n", "Version 1.2.0 | \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang | \n", "==================================================================================================================\n", "Problem specifications: | \n", " # of variables : 30 | \n", " # of inequality constraints : 0 | \n", " # of equality constraints : 1 | \n", "==================================================================================================================\n", " | <--- Penalty Function --> | | Total Violation | <--- Line Search ---> | <- Stationarity -> | \n", "Iter | Mu | Value | Objective | Ineq | Eq | SD | Evals | t | Grads | Value | \n", "=====|===========================|================|=================|=======================|====================|\n", " 0 | 1.000000 | 0.61751624522 | 0.61751624522 | - | 0.000000 | - | 1 | 0.000000 | 1 | 0.054664 | \n", " 10 | 1.000000 | 0.60573380055 | 0.60513582468 | - | 5.98e-04 | S | 1 | 1.000000 | 1 | 0.024968 | \n", " 20 | 1.000000 | 0.58456516016 | 0.58301955756 | - | 0.001546 | S | 1 | 1.000000 | 1 | 0.043517 | \n", " 30 | 1.000000 | 0.50113197499 | 0.49475409554 | - | 0.006378 | S | 3 | 0.250000 | 1 | 0.121253 | \n", " 40 | 1.000000 | 0.49278124194 | 0.49260444460 | - | 1.77e-04 | S | 4 | 0.125000 | 1 | 0.037304 | \n", " 50 | 1.000000 | 0.49225009818 | 0.49217494722 | - | 7.52e-05 | S | 5 | 0.062500 | 1 | 0.032163 | \n", " 60 | 1.000000 | 0.49212731751 | 0.49208854433 | - | 3.88e-05 | S | 4 | 0.125000 | 1 | 0.051779 | \n", " 70 | 1.000000 | 0.49203371691 | 0.49201049130 | - | 2.32e-05 | S | 4 | 0.125000 | 1 | 0.054529 | \n", "==================================================================================================================\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible | \n", "Optimization results: | \n", "==================================================================================================================\n", " F | | | 0.49205285240 | - | 3.35e-07 | | | | | | \n", " B | | | 0.49205285240 | - | 3.35e-07 | | | | | | \n", " MF | | | 0.61751624522 | - | 0.000000 | | | | | | \n", "==================================================================================================================\n", "Iterations: 77 | \n", "Function evaluations: 256 | \n", "PyGRANSO termination code: 2 --- target objective reached at point feasible to tolerances. | \n", "==================================================================================================================\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "opts = pygransoStruct()\n", "opts.torch_device = device\n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "opts.print_frequency = 10\n", "opts.print_ascii = True\n", "\n", "\n", "opts.fvalquit = 0.4963\n", "\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "print(max(abs(soln.final.x))) # should be close to 1" ] }, { "cell_type": "markdown", "id": "425a5bfb", "metadata": {}, "source": [ "Set opt_tol. Tolerance for reaching (approximate) optimality/stationarity.\n", "See opts.ngrad, opts.evaldist, and the description of PyGRANSO's \n", "output argument soln, specifically the subsubfield .dnorm for more\n", "information.\n", "\n", "In the result below, PyGRANSO terminated when stationarity is below 1e-4" ] }, { "cell_type": "code", "execution_count": 7, "id": "150ea95d", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "###### QP SOLVER NOTICE #########################################################################\n", "# PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, #\n", "# the default is osqp. Users may provide their own wrapper for the QP solver. #\n", "# To disable this notice, set opts.quadprog_info_msg = False #\n", "#################################################################################################\n", "==================================================================================================================\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation | \n", "Version 1.2.0 | \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang | \n", "==================================================================================================================\n", "Problem specifications: | \n", " # of variables : 30 | \n", " # of inequality constraints : 0 | \n", " # of equality constraints : 1 | \n", "==================================================================================================================\n", " | <--- Penalty Function --> | | Total Violation | <--- Line Search ---> | <- Stationarity -> | \n", "Iter | Mu | Value | Objective | Ineq | Eq | SD | Evals | t | Grads | Value | \n", "=====|===========================|================|=================|=======================|====================|\n", " 0 | 1.000000 | 0.61751624522 | 0.61751624522 | - | 0.000000 | - | 1 | 0.000000 | 1 | 0.054664 | \n", " 10 | 1.000000 | 0.60573380055 | 0.60513582468 | - | 5.98e-04 | S | 1 | 1.000000 | 1 | 0.024968 | \n", " 20 | 1.000000 | 0.58456516016 | 0.58301955756 | - | 0.001546 | S | 1 | 1.000000 | 1 | 0.043517 | \n", " 30 | 1.000000 | 0.50113197499 | 0.49475409554 | - | 0.006378 | S | 3 | 0.250000 | 1 | 0.121253 | \n", " 40 | 1.000000 | 0.49278124194 | 0.49260444460 | - | 1.77e-04 | S | 4 | 0.125000 | 1 | 0.037304 | \n", " 50 | 1.000000 | 0.49225009818 | 0.49217494722 | - | 7.52e-05 | S | 5 | 0.062500 | 1 | 0.032163 | \n", " 60 | 1.000000 | 0.49212731751 | 0.49208854433 | - | 3.88e-05 | S | 4 | 0.125000 | 1 | 0.051779 | \n", " 70 | 1.000000 | 0.49203371691 | 0.49201049130 | - | 2.32e-05 | S | 4 | 0.125000 | 1 | 0.054529 | \n", " 80 | 1.000000 | 0.49197689465 | 0.49197679422 | - | 1.00e-07 | S | 2 | 0.500000 | 1 | 0.001300 | \n", " 90 | 1.000000 | 0.49194701030 | 0.49194698105 | - | 2.93e-08 | S | 5 | 0.062500 | 5 | 6.92e-05 | \n", "==================================================================================================================\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible | \n", "Optimization results: | \n", "==================================================================================================================\n", " F | | | 0.49194698105 | - | 2.93e-08 | | | | | | \n", " B | | | 0.49194698105 | - | 2.93e-08 | | | | | | \n", " MF | | | 0.61751624522 | - | 0.000000 | | | | | | \n", "==================================================================================================================\n", "Iterations: 90 | \n", "Function evaluations: 299 | \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. | \n", "==================================================================================================================\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "opts = pygransoStruct()\n", "opts.torch_device = device \n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "opts.print_frequency = 10\n", "opts.print_ascii = True\n", "\n", "opts.opt_tol = 1e-4 # default 1e-8\n", "\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "print(max(abs(soln.final.x))) # should be close to 1" ] }, { "cell_type": "markdown", "id": "e45d24d4", "metadata": {}, "source": [ "There are multiple other settings. Please uncomment to try them. Detailed description can be found by typing\n", "\n", "import pygransoOptionsAdvanced\n", "\n", "help(pygransoOptionsAdvanced)" ] }, { "cell_type": "code", "execution_count": 8, "id": "d6f3b9bd", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "\u001b[33m╔═════ QP SOLVER NOTICE ════════════════════════════════════════════════════════════════════════╗\n", "\u001b[0m\u001b[33m║ PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, ║\n", "\u001b[0m\u001b[33m║ the default is osqp. Users may provide their own wrapper for the QP solver. ║\n", "\u001b[0m\u001b[33m║ To disable this notice, set opts.quadprog_info_msg = False ║\n", "\u001b[0m\u001b[33m╚═══════════════════════════════════════════════════════════════════════════════════════════════╝\n", "\u001b[0m═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╗\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation ║ \n", "Version 1.2.0 ║ \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╣\n", "Problem specifications: ║ \n", " # of variables : 30 ║ \n", " # of inequality constraints : 0 ║ \n", " # of equality constraints : 1 ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " ║ <--- Penalty Function --> ║ ║ Total Violation ║ <--- Line Search ---> ║ <- Stationarity -> ║ \n", "Iter ║ Mu │ Value ║ Objective ║ Ineq │ Eq ║ SD │ Evals │ t ║ Grads │ Value ║ \n", "═════╬═══════════════════════════╬════════════════╬═════════════════╬═══════════════════════╬════════════════════╣\n", " 0 ║ 1.000000 │ 0.61751624522 ║ 0.61751624522 ║ - │ 0.000000 ║ - │ 1 │ 0.000000 ║ 1 │ 0.054664 ║ \n", " 10 ║ 1.000000 │ 0.60573380055 ║ 0.60513582468 ║ - │ 5.98e-04 ║ S │ 1 │ 1.000000 ║ 1 │ 0.024968 ║ \n", " 20 ║ 1.000000 │ 0.58456516016 ║ 0.58301955756 ║ - │ 0.001546 ║ S │ 1 │ 1.000000 ║ 1 │ 0.043517 ║ \n", " 30 ║ 1.000000 │ 0.50113197499 ║ 0.49475409554 ║ - │ 0.006378 ║ S │ 3 │ 0.250000 ║ 1 │ 0.121253 ║ \n", " 40 ║ 1.000000 │ 0.49278124194 ║ 0.49260444460 ║ - │ 1.77e-04 ║ S │ 4 │ 0.125000 ║ 1 │ 0.037304 ║ \n", " 50 ║ 1.000000 │ 0.49225009818 ║ 0.49217494722 ║ - │ 7.52e-05 ║ S │ 5 │ 0.062500 ║ 1 │ 0.032163 ║ \n", " 60 ║ 1.000000 │ 0.49212731751 ║ 0.49208854433 ║ - │ 3.88e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.051779 ║ \n", " 70 ║ 1.000000 │ 0.49203371691 ║ 0.49201049130 ║ - │ 2.32e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.054529 ║ \n", " 80 ║ 1.000000 │ 0.49197689465 ║ 0.49197679422 ║ - │ 1.00e-07 ║ S │ 2 │ 0.500000 ║ 1 │ 0.001300 ║ \n", " 90 ║ 1.000000 │ 0.49194701030 ║ 0.49194698105 ║ - │ 2.93e-08 ║ S │ 5 │ 0.062500 ║ 5 │ 6.92e-05 ║ \n", " 100 ║ 1.000000 │ 0.49194382838 ║ 0.49194381415 ║ - │ 1.42e-08 ║ S │ 6 │ 0.031250 ║ 10 │ 1.28e-05 ║ \n", " 110 ║ 1.000000 │ 0.49194277900 ║ 0.49194277111 ║ - │ 7.88e-09 ║ S │ 5 │ 0.062500 ║ 18 │ 6.93e-07 ║ \n", " 120 ║ 1.000000 │ 0.49194243076 ║ 0.49194242538 ║ - │ 5.38e-09 ║ S │ 6 │ 0.031250 ║ 27 │ 1.99e-07 ║ \n", " 130 ║ 1.000000 │ 0.49194218055 ║ 0.49194217869 ║ - │ 1.87e-09 ║ S │ 4 │ 0.125000 ║ 37 │ 5.60e-08 ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible ║ \n", "Optimization results: ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " F ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " B ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " MF ║ │ ║ 0.61751624522 ║ - │ 0.000000 ║ │ │ ║ │ ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "Iterations: 134 ║ \n", "Function evaluations: 525 ║ \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╝\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "opts = pygransoStruct()\n", "opts.torch_device = device \n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "opts.print_frequency = 10\n", "\n", "# Please uncomment to try different settings\n", "\n", "# Tolerance for determining when the relative decrease in the penalty\n", "# function is sufficiently small. PyGRANSO will terminate if when \n", "# the relative decrease in the penalty function is at or below this\n", "# tolerance and the current iterate is feasible to tolerances.\n", "# Generally, we don't recommend using this feature since small steps\n", "# are not necessarily indicative of being near a stationary point,\n", "# particularly for nonsmooth problems.\n", "\n", "# Termination Code 1\n", "# opts.rel_tol = 1e-2 # default 0\n", "\n", "# Tolerance for how small of a step the line search will attempt\n", "# before terminating.\n", "\n", "# Termination Code 6 or 7\n", "# opts.step_tol = 1e-6 # default 1e-12\n", "# opts.step_tol = 1e-3\n", "\n", "# Acceptable total violation tolerance of the equality constraints.\n", "# opts.viol_eq_tol = 1e-12# default 1e-6, make it smaller will make current point harder to be considered as feasible\n", "\n", "# Quit if the elapsed clock time in seconds exceeds this. unit: second\n", "# opts.maxclocktime = 1.\n", "\n", "# Number of characters wide to print values for the penalty function,\n", "# the objective function, and the total violations of the inequality \n", "# and equality constraints. \n", "# opts.print_width = 9\n", "\n", "# PyGRANSO's uses orange\n", "# printing to highlight pertinent information. However, the user\n", "# is the given option to disable it, since support cannot be\n", "# guaranteed (since it is an undocumented feature).\n", "# opts.print_use_orange = False\n", "\n", "# opts.init_step_size = 1e-2\n", "# opts.search_direction_rescaling = True\n", "\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "print(max(abs(soln.final.x))) # should be close to 1" ] }, { "cell_type": "markdown", "id": "7a825d44", "metadata": {}, "source": [ "**(For Advanced User)** Users can specify analytical gradients instead of using the auto-differentiation feature. Please check the documentation in `pygranso.py` on the required format of `combined_fn`" ] }, { "cell_type": "code", "execution_count": 9, "id": "4f3b437d", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "\u001b[33m╔═════ QP SOLVER NOTICE ════════════════════════════════════════════════════════════════════════╗\n", "\u001b[0m\u001b[33m║ PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, ║\n", "\u001b[0m\u001b[33m║ the default is osqp. Users may provide their own wrapper for the QP solver. ║\n", "\u001b[0m\u001b[33m║ To disable this notice, set opts.quadprog_info_msg = False ║\n", "\u001b[0m\u001b[33m╚═══════════════════════════════════════════════════════════════════════════════════════════════╝\n", "\u001b[0m═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╗\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation ║ \n", "Version 1.2.0 ║ \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╣\n", "Problem specifications: ║ \n", " # of variables : 30 ║ \n", " # of inequality constraints : 0 ║ \n", " # of equality constraints : 1 ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " ║ <--- Penalty Function --> ║ ║ Total Violation ║ <--- Line Search ---> ║ <- Stationarity -> ║ \n", "Iter ║ Mu │ Value ║ Objective ║ Ineq │ Eq ║ SD │ Evals │ t ║ Grads │ Value ║ \n", "═════╬═══════════════════════════╬════════════════╬═════════════════╬═══════════════════════╬════════════════════╣\n", " 0 ║ 1.000000 │ 0.61751624522 ║ 0.61751624522 ║ - │ 0.000000 ║ - │ 1 │ 0.000000 ║ 1 │ 0.054664 ║ \n", " 10 ║ 1.000000 │ 0.60573380055 ║ 0.60513582468 ║ - │ 5.98e-04 ║ S │ 1 │ 1.000000 ║ 1 │ 0.024968 ║ \n", " 20 ║ 1.000000 │ 0.58456516016 ║ 0.58301955756 ║ - │ 0.001546 ║ S │ 1 │ 1.000000 ║ 1 │ 0.043517 ║ \n", " 30 ║ 1.000000 │ 0.50113197499 ║ 0.49475409554 ║ - │ 0.006378 ║ S │ 3 │ 0.250000 ║ 1 │ 0.121253 ║ \n", " 40 ║ 1.000000 │ 0.49278124194 ║ 0.49260444460 ║ - │ 1.77e-04 ║ S │ 4 │ 0.125000 ║ 1 │ 0.037304 ║ \n", " 50 ║ 1.000000 │ 0.49225009818 ║ 0.49217494722 ║ - │ 7.52e-05 ║ S │ 5 │ 0.062500 ║ 1 │ 0.032163 ║ \n", " 60 ║ 1.000000 │ 0.49212731751 ║ 0.49208854433 ║ - │ 3.88e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.051779 ║ \n", " 70 ║ 1.000000 │ 0.49203371691 ║ 0.49201049130 ║ - │ 2.32e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.054529 ║ \n", " 80 ║ 1.000000 │ 0.49197689465 ║ 0.49197679422 ║ - │ 1.00e-07 ║ S │ 2 │ 0.500000 ║ 1 │ 0.001300 ║ \n", " 90 ║ 1.000000 │ 0.49194701030 ║ 0.49194698105 ║ - │ 2.93e-08 ║ S │ 5 │ 0.062500 ║ 5 │ 6.92e-05 ║ \n", " 100 ║ 1.000000 │ 0.49194382838 ║ 0.49194381415 ║ - │ 1.42e-08 ║ S │ 6 │ 0.031250 ║ 10 │ 1.28e-05 ║ \n", " 110 ║ 1.000000 │ 0.49194277900 ║ 0.49194277111 ║ - │ 7.88e-09 ║ S │ 5 │ 0.062500 ║ 18 │ 6.93e-07 ║ \n", " 120 ║ 1.000000 │ 0.49194243076 ║ 0.49194242538 ║ - │ 5.38e-09 ║ S │ 6 │ 0.031250 ║ 27 │ 1.99e-07 ║ \n", " 130 ║ 1.000000 │ 0.49194218055 ║ 0.49194217869 ║ - │ 1.87e-09 ║ S │ 4 │ 0.125000 ║ 37 │ 5.66e-08 ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible ║ \n", "Optimization results: ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " F ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " B ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " MF ║ │ ║ 0.61751624522 ║ - │ 0.000000 ║ │ │ ║ │ ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "Iterations: 134 ║ \n", "Function evaluations: 525 ║ \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╝\n", "Total Wall Time: 2.052833318710327s\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "# Without AD\n", "def comb_fn(X_struct):\n", " q = X_struct.q\n", "\n", " # objective function\n", " qtY = q.T @ Y\n", " f = 1/m * torch.norm(qtY, p = 1).item()\n", " f_grad = 1/m*Y@torch.sign(Y.T@q)\n", "\n", " # inequality constraint, matrix form\n", " ci = None\n", " ci_grad = None\n", "\n", " # equality constraint \n", " ce = q.T @ q - 1\n", " ce_grad = 2*q\n", "\n", " return [f,f_grad,ci,ci_grad,ce,ce_grad]\n", "\n", "opts = pygransoStruct()\n", "opts.torch_device = device\n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "opts.print_frequency = 10\n", "opts.globalAD = False # disable global auto-differentiation\n", "\n", "start = time.time()\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "end = time.time()\n", "print(\"Total Wall Time: {}s\".format(end - start))\n", "print(max(abs(soln.final.x))) # should be close to 1" ] }, { "cell_type": "markdown", "id": "5c6556b1", "metadata": {}, "source": [ "**(For Advanced User)** Alternatively, users can use the auto-differentiation feature partially." ] }, { "cell_type": "code", "execution_count": 10, "id": "50b0c894", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "\u001b[33m╔═════ QP SOLVER NOTICE ════════════════════════════════════════════════════════════════════════╗\n", "\u001b[0m\u001b[33m║ PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, ║\n", "\u001b[0m\u001b[33m║ the default is osqp. Users may provide their own wrapper for the QP solver. ║\n", "\u001b[0m\u001b[33m║ To disable this notice, set opts.quadprog_info_msg = False ║\n", "\u001b[0m\u001b[33m╚═══════════════════════════════════════════════════════════════════════════════════════════════╝\n", "\u001b[0m═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╗\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation ║ \n", "Version 1.2.0 ║ \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╣\n", "Problem specifications: ║ \n", " # of variables : 30 ║ \n", " # of inequality constraints : 0 ║ \n", " # of equality constraints : 1 ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " ║ <--- Penalty Function --> ║ ║ Total Violation ║ <--- Line Search ---> ║ <- Stationarity -> ║ \n", "Iter ║ Mu │ Value ║ Objective ║ Ineq │ Eq ║ SD │ Evals │ t ║ Grads │ Value ║ \n", "═════╬═══════════════════════════╬════════════════╬═════════════════╬═══════════════════════╬════════════════════╣\n", " 0 ║ 1.000000 │ 0.61751624522 ║ 0.61751624522 ║ - │ 0.000000 ║ - │ 1 │ 0.000000 ║ 1 │ 0.054664 ║ \n", " 10 ║ 1.000000 │ 0.60573380055 ║ 0.60513582468 ║ - │ 5.98e-04 ║ S │ 1 │ 1.000000 ║ 1 │ 0.024968 ║ \n", " 20 ║ 1.000000 │ 0.58456516016 ║ 0.58301955756 ║ - │ 0.001546 ║ S │ 1 │ 1.000000 ║ 1 │ 0.043517 ║ \n", " 30 ║ 1.000000 │ 0.50113197499 ║ 0.49475409554 ║ - │ 0.006378 ║ S │ 3 │ 0.250000 ║ 1 │ 0.121253 ║ \n", " 40 ║ 1.000000 │ 0.49278124194 ║ 0.49260444460 ║ - │ 1.77e-04 ║ S │ 4 │ 0.125000 ║ 1 │ 0.037304 ║ \n", " 50 ║ 1.000000 │ 0.49225009818 ║ 0.49217494722 ║ - │ 7.52e-05 ║ S │ 5 │ 0.062500 ║ 1 │ 0.032163 ║ \n", " 60 ║ 1.000000 │ 0.49212731751 ║ 0.49208854433 ║ - │ 3.88e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.051779 ║ \n", " 70 ║ 1.000000 │ 0.49203371691 ║ 0.49201049130 ║ - │ 2.32e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.054529 ║ \n", " 80 ║ 1.000000 │ 0.49197689465 ║ 0.49197679422 ║ - │ 1.00e-07 ║ S │ 2 │ 0.500000 ║ 1 │ 0.001300 ║ \n", " 90 ║ 1.000000 │ 0.49194701030 ║ 0.49194698105 ║ - │ 2.93e-08 ║ S │ 5 │ 0.062500 ║ 5 │ 6.92e-05 ║ \n", " 100 ║ 1.000000 │ 0.49194382838 ║ 0.49194381415 ║ - │ 1.42e-08 ║ S │ 6 │ 0.031250 ║ 10 │ 1.28e-05 ║ \n", " 110 ║ 1.000000 │ 0.49194277900 ║ 0.49194277111 ║ - │ 7.88e-09 ║ S │ 5 │ 0.062500 ║ 18 │ 6.93e-07 ║ \n", " 120 ║ 1.000000 │ 0.49194243076 ║ 0.49194242538 ║ - │ 5.38e-09 ║ S │ 6 │ 0.031250 ║ 27 │ 1.99e-07 ║ \n", " 130 ║ 1.000000 │ 0.49194218055 ║ 0.49194217869 ║ - │ 1.87e-09 ║ S │ 4 │ 0.125000 ║ 37 │ 5.66e-08 ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible ║ \n", "Optimization results: ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " F ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " B ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " MF ║ │ ║ 0.61751624522 ║ - │ 0.000000 ║ │ │ ║ │ ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "Iterations: 134 ║ \n", "Function evaluations: 525 ║ \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╝\n", "Total Wall Time: 2.0349695682525635s\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "# import the AD function\n", "from pygranso.private.tensor2vec import getCiGradVec\n", "\n", "# partial AD\n", "def comb_fn(X_struct):\n", " q = X_struct.q\n", " q.requires_grad_(True)\n", "\n", " # objective function\n", " q_tmp = q.detach().clone()\n", " qtY = q_tmp.T @ Y\n", " f = 1/m * torch.norm(qtY, p = 1).item()\n", " f_grad = 1/m*Y@torch.sign(Y.T@q_tmp)\n", "\n", " # inequality constraint, matrix form\n", " ci = None\n", " ci_grad = None\n", "\n", " # equality constraint \n", " ce = q.T @ q - 1\n", " # ce_grad = 2*q\n", " ce_grad = getCiGradVec(nvar=n,nconstr_ci_total=1,var_dim_map=var_in,X=X_struct,ci_vec_torch=ce,torch_device=device,double_precision=torch.double)\n", "\n", " # return value must be detached from the computational graph\n", " ce = ce.detach()\n", "\n", " return [f,f_grad,ci,ci_grad,ce,ce_grad]\n", "\n", "opts = pygransoStruct()\n", "opts.torch_device = device\n", "opts.maxit = 500\n", "np.random.seed(1)\n", "x0 = norm.ppf(np.random.rand(n,1))\n", "x0 /= la.norm(x0,2)\n", "opts.x0 = torch.from_numpy(x0).to(device=device, dtype=torch.double)\n", "opts.print_frequency = 10\n", "opts.globalAD = False # disable global auto-differentiation\n", "\n", "\n", "start = time.time()\n", "soln = pygranso(var_spec = var_in,combined_fn = comb_fn,user_opts = opts)\n", "end = time.time()\n", "print(\"Total Wall Time: {}s\".format(end - start))\n", "print(max(abs(soln.final.x))) # should be close to 1**(For Advanced User)** Users can specify analytical gradients instead of using the AD feature" ] }, { "cell_type": "markdown", "id": "9ca8b410", "metadata": {}, "source": [ "## Different Set-Up\n", "\n", "**(Optional)** Using torch.nn to set up the user-provided function. This setting is important in solving constrained deep learning problems." ] }, { "cell_type": "markdown", "id": "06480f09", "metadata": {}, "source": [ "### Modules Importing\n", "Import all necessary modules and add PyGRANSO src folder to system path." ] }, { "cell_type": "code", "execution_count": 11, "id": "64ac4009", "metadata": {}, "outputs": [], "source": [ "import time\n", "import numpy as np\n", "import torch\n", "import numpy.linalg as la\n", "from scipy.stats import norm\n", "from pygranso.pygranso import pygranso\n", "from pygranso.pygransoStruct import pygransoStruct\n", "\n", "from pygranso.private.getNvar import getNvarTorch\n", "import torch.nn as nn" ] }, { "cell_type": "markdown", "id": "12a4d15e", "metadata": {}, "source": [ "### Initialization\n", "Specify torch device, create torch model and generate data\n", "\n", "Use GPU for this problem. If no cuda device available, please set device = torch.device('cpu')" ] }, { "cell_type": "code", "execution_count": 12, "id": "71ab27d2", "metadata": {}, "outputs": [], "source": [ "device = torch.device('cuda')\n", "\n", "class Dict_Learning(nn.Module):\n", " \n", " def __init__(self,n):\n", " super().__init__()\n", " np.random.seed(1)\n", " q0 = norm.ppf(np.random.rand(n,1))\n", " q0 /= la.norm(q0,2)\n", " self.q = nn.Parameter( torch.from_numpy(q0) )\n", " \n", " def forward(self, Y,m):\n", " qtY = self.q.T @ Y\n", " f = 1/m * torch.norm(qtY, p = 1)\n", " return f\n", "\n", "## Data initialization\n", "n = 30\n", "np.random.seed(1)\n", "m = 10*n**2 # sample complexity\n", "theta = 0.3 # sparsity level\n", "Y = norm.ppf(np.random.rand(n,m)) * (norm.ppf(np.random.rand(n,m)) <= theta) # Bernoulli-Gaussian model\n", "# All the user-provided data (vector/matrix/tensor) must be in torch tensor format.\n", "# As PyTorch tensor is single precision by default, one must explicitly set `dtype=torch.double`.\n", "# Also, please make sure the device of provided torch tensor is the same as opts.torch_device.\n", "Y = torch.from_numpy(Y).to(device=device, dtype=torch.double)\n", "\n", "torch.manual_seed(0)\n", "\n", "model = Dict_Learning(n).to(device=device, dtype=torch.double)" ] }, { "cell_type": "markdown", "id": "9a7d01fc", "metadata": {}, "source": [ "### Function Set-Up\n", "Encode the optimization variables, and objective and constraint functions.\n", "\n", "Note: please strictly follow the format of comb_fn, which will be used in the PyGRANSO main algortihm." ] }, { "cell_type": "code", "execution_count": 13, "id": "e897ee5a", "metadata": {}, "outputs": [], "source": [ "def user_fn(model,Y,m):\n", " # objective function \n", " f = model(Y,m)\n", "\n", " q = list(model.parameters())[0]\n", "\n", " # inequality constraint\n", " ci = None\n", "\n", " # equality constraint \n", " ce = pygransoStruct()\n", " ce.c1 = q.T @ q - 1\n", "\n", " return [f,ci,ce]\n", "\n", "comb_fn = lambda model : user_fn(model,Y,m)" ] }, { "cell_type": "markdown", "id": "2db9c25d", "metadata": {}, "source": [ "### User Options\n", "Specify user-defined options for PyGRANSO" ] }, { "cell_type": "code", "execution_count": 14, "id": "265fb1bc", "metadata": {}, "outputs": [], "source": [ "opts = pygransoStruct()\n", "opts.torch_device = device\n", "opts.maxit = 500\n", "np.random.seed(1)\n", "nvar = getNvarTorch(model.parameters())\n", "opts.x0 = torch.nn.utils.parameters_to_vector(model.parameters()).detach().reshape(nvar,1)\n", "\n", "opts.print_frequency = 10" ] }, { "cell_type": "markdown", "id": "345a456d", "metadata": {}, "source": [ "### Main Algorithm" ] }, { "cell_type": "code", "execution_count": 15, "id": "478e7b1c", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", "\u001b[33m╔═════ QP SOLVER NOTICE ════════════════════════════════════════════════════════════════════════╗\n", "\u001b[0m\u001b[33m║ PyGRANSO requires a quadratic program (QP) solver that has a quadprog-compatible interface, ║\n", "\u001b[0m\u001b[33m║ the default is osqp. Users may provide their own wrapper for the QP solver. ║\n", "\u001b[0m\u001b[33m║ To disable this notice, set opts.quadprog_info_msg = False ║\n", "\u001b[0m\u001b[33m╚═══════════════════════════════════════════════════════════════════════════════════════════════╝\n", "\u001b[0m═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╗\n", "PyGRANSO: A PyTorch-enabled port of GRANSO with auto-differentiation ║ \n", "Version 1.2.0 ║ \n", "Licensed under the AGPLv3, Copyright (C) 2021-2022 Tim Mitchell and Buyun Liang ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╣\n", "Problem specifications: ║ \n", " # of variables : 30 ║ \n", " # of inequality constraints : 0 ║ \n", " # of equality constraints : 1 ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " ║ <--- Penalty Function --> ║ ║ Total Violation ║ <--- Line Search ---> ║ <- Stationarity -> ║ \n", "Iter ║ Mu │ Value ║ Objective ║ Ineq │ Eq ║ SD │ Evals │ t ║ Grads │ Value ║ \n", "═════╬═══════════════════════════╬════════════════╬═════════════════╬═══════════════════════╬════════════════════╣\n", " 0 ║ 1.000000 │ 0.61751624522 ║ 0.61751624522 ║ - │ 0.000000 ║ - │ 1 │ 0.000000 ║ 1 │ 0.054664 ║ \n", " 10 ║ 1.000000 │ 0.60573380055 ║ 0.60513582468 ║ - │ 5.98e-04 ║ S │ 1 │ 1.000000 ║ 1 │ 0.024968 ║ \n", " 20 ║ 1.000000 │ 0.58456516016 ║ 0.58301955756 ║ - │ 0.001546 ║ S │ 1 │ 1.000000 ║ 1 │ 0.043517 ║ \n", " 30 ║ 1.000000 │ 0.50113197499 ║ 0.49475409554 ║ - │ 0.006378 ║ S │ 3 │ 0.250000 ║ 1 │ 0.121253 ║ \n", " 40 ║ 1.000000 │ 0.49278124194 ║ 0.49260444460 ║ - │ 1.77e-04 ║ S │ 4 │ 0.125000 ║ 1 │ 0.037304 ║ \n", " 50 ║ 1.000000 │ 0.49225009818 ║ 0.49217494722 ║ - │ 7.52e-05 ║ S │ 5 │ 0.062500 ║ 1 │ 0.032163 ║ \n", " 60 ║ 1.000000 │ 0.49212731751 ║ 0.49208854433 ║ - │ 3.88e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.051779 ║ \n", " 70 ║ 1.000000 │ 0.49203371691 ║ 0.49201049130 ║ - │ 2.32e-05 ║ S │ 4 │ 0.125000 ║ 1 │ 0.054529 ║ \n", " 80 ║ 1.000000 │ 0.49197689465 ║ 0.49197679422 ║ - │ 1.00e-07 ║ S │ 2 │ 0.500000 ║ 1 │ 0.001300 ║ \n", " 90 ║ 1.000000 │ 0.49194701030 ║ 0.49194698105 ║ - │ 2.93e-08 ║ S │ 5 │ 0.062500 ║ 5 │ 6.92e-05 ║ \n", " 100 ║ 1.000000 │ 0.49194382838 ║ 0.49194381415 ║ - │ 1.42e-08 ║ S │ 6 │ 0.031250 ║ 10 │ 1.28e-05 ║ \n", " 110 ║ 1.000000 │ 0.49194277900 ║ 0.49194277111 ║ - │ 7.88e-09 ║ S │ 5 │ 0.062500 ║ 18 │ 6.93e-07 ║ \n", " 120 ║ 1.000000 │ 0.49194243076 ║ 0.49194242538 ║ - │ 5.38e-09 ║ S │ 6 │ 0.031250 ║ 27 │ 1.99e-07 ║ \n", " 130 ║ 1.000000 │ 0.49194218055 ║ 0.49194217869 ║ - │ 1.87e-09 ║ S │ 4 │ 0.125000 ║ 37 │ 5.60e-08 ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "F = final iterate, B = Best (to tolerance), MF = Most Feasible ║ \n", "Optimization results: ║ \n", "═════╦═══════════════════════════╦════════════════╦═════════════════╦═══════════════════════╦════════════════════╣\n", " F ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " B ║ │ ║ 0.49194215780 ║ - │ 1.30e-09 ║ │ │ ║ │ ║ \n", " MF ║ │ ║ 0.61751624522 ║ - │ 0.000000 ║ │ │ ║ │ ║ \n", "═════╩═══════════════════════════╩════════════════╩═════════════════╩═══════════════════════╩════════════════════╣\n", "Iterations: 134 ║ \n", "Function evaluations: 525 ║ \n", "PyGRANSO termination code: 0 --- converged to stationarity and feasibility tolerances. ║ \n", "═════════════════════════════════════════════════════════════════════════════════════════════════════════════════╝\n", "Total Wall Time: 2.3208694458007812s\n", "tensor([1.0000], device='cuda:0', dtype=torch.float64)\n" ] } ], "source": [ "start = time.time()\n", "soln = pygranso(var_spec= model, combined_fn = comb_fn, user_opts = opts)\n", "end = time.time()\n", "print(\"Total Wall Time: {}s\".format(end - start))\n", "print(max(abs(soln.final.x))) # should be close to 1" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.7" } }, "nbformat": 4, "nbformat_minor": 5 }