From 1b38d23e4a9911e8f7377396c1c1fd51e95ea71b Mon Sep 17 00:00:00 2001 From: Jonas Adler Date: Thu, 14 Jul 2016 14:41:38 +0200 Subject: [PATCH 1/5] Add as_cvx_operator and example --- README.md | 1 + doc/source/refs.rst | 3 + examples/solvers/proximal_lang_poisson.py | 67 ++++++++++++++ examples/solvers/proximal_lang_tomography.py | 93 ++++++++++++++++++++ odl/operator/oputils.py | 48 +++++++++- setup.py | 3 +- 6 files changed, 213 insertions(+), 2 deletions(-) create mode 100644 examples/solvers/proximal_lang_poisson.py create mode 100644 examples/solvers/proximal_lang_tomography.py diff --git a/README.md b/README.md index dc802b83ffe..fc1e97c5c25 100644 --- a/README.md +++ b/README.md @@ -56,6 +56,7 @@ where FOO, BAR, etc are the dependencies. Others require more complicated instal | [FFTW](https://github.com/pyFFTW/pyFFTW) | Accelerated [FourierTransform](http://odl.readthedocs.io/generated/odl.trafos.fourier.FourierTransform.html) | fftw | | [PyWavelets](https://github.com/PyWavelets/pywt) | Computation of the [WaveletTransform](http://odl.readthedocs.io/generated/odl.trafos.wavelet.WaveletTransform.html) | pywavelets | | [matplotlib](http://matplotlib.org/) | Visualization through the [show](http://odl.readthedocs.io/generated/odl.discr.lp_discr.DiscreteLpVector.show.html) command | show | +| [proximal](http://github.com/comp-imaging/ProxImaL) | Solution of some convex optimization problems | proximal | | [pytest](http://pytest.org/latest/) | Unit tests | testing | diff --git a/doc/source/refs.rst b/doc/source/refs.rst index ed57397b4cb..248a4c79dba 100644 --- a/doc/source/refs.rst +++ b/doc/source/refs.rst @@ -42,6 +42,9 @@ References .. [GNS2009] Griva, I, Nash, S G, and Sofer, A. *Linear and nonlinear optimization*. Siam, 2009. +.. [Hei+2016] Heide F. et. al. *ProxImaL: Efficient Image Optimization using + Proximal Algorithms*. ACM Transactions on Graphics (TOG), 2016. + .. [Kva1991] Kvaalen, E. *A faster Broyden method*. BIT Numerical Mathematics 31 (1991), pp 369--372. diff --git a/examples/solvers/proximal_lang_poisson.py b/examples/solvers/proximal_lang_poisson.py new file mode 100644 index 00000000000..bd29c7e3c2b --- /dev/null +++ b/examples/solvers/proximal_lang_poisson.py @@ -0,0 +1,67 @@ +# Copyright 2014-2016 The ODL development group +# +# This file is part of ODL. +# +# ODL is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# ODL is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with ODL. If not, see . + +"""Poissons problem using the proximal solver. + +Solves the optimization problem + + min_x 10 ||laplacian(x) - g||_2^2 + || |grad(x)| ||_1 + +Where ``laplacian`` is the spatial laplacian, ``grad`` the spatial +gradient and ``g`` is given noisy data. +""" + +# Imports for common Python 2/3 codebase +from __future__ import print_function, division, absolute_import +from future import standard_library +standard_library.install_aliases() + +import numpy as np +import odl +import proximal + +# Create space, a square from [0, 0] to [100, 100] with (100 x 100) points +space = odl.uniform_discr([0, 0], [100, 100], [100, 100]) + +# Create odl operator for laplacian +laplacian = odl.Laplacian(space) + +# Create right hand side +phantom = odl.phantom.shepp_logan(space, modified=True) +phantom.show('phantom') +rhs = laplacian(phantom) +rhs += odl.phantom.white_noise(space) * np.std(rhs) * 0.1 +rhs.show('rhs') + +# Convert laplacian to cvx operator +cvx_laplacian = odl.as_proximal_lang_operator(laplacian) + +# Convert to array +rhs_arr = rhs.asarray() + +# Set up optimization problem +x = proximal.Variable(space.shape) +funcs = [10 * proximal.sum_squares(cvx_laplacian(x) - rhs_arr), + proximal.norm1(proximal.grad(x))] + +# Solve the problem +prob = proximal.Problem(funcs) +prob.solve(verbose=True) + +# Convert back to odl and display result +result_odl = space.element(x.value) +result_odl.show('result') diff --git a/examples/solvers/proximal_lang_tomography.py b/examples/solvers/proximal_lang_tomography.py new file mode 100644 index 00000000000..f14d91a600a --- /dev/null +++ b/examples/solvers/proximal_lang_tomography.py @@ -0,0 +1,93 @@ +# Copyright 2014-2016 The ODL development group +# +# This file is part of ODL. +# +# ODL is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# ODL is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with ODL. If not, see . + +"""Total variation tomography using the proximal solver. + +Solves the optimization problem + + min_{0 <= x <= 1} ||A(x) - g||_2^2 + 0.2 || |grad(x)| ||_1 + +Where ``A`` is a parallel beam forward projector, ``grad`` the spatial +gradient and ``g`` is given noisy data. +""" + +# Imports for common Python 2/3 codebase +from __future__ import print_function, division, absolute_import +from future import standard_library +standard_library.install_aliases() + +import numpy as np +import odl +import proximal + + +# --- Set up the forward operator (ray transform) --- # + + +# Discrete reconstruction space: discretized functions on the rectangle +# [-20, 20]^2 with 300 samples per dimension. +reco_space = odl.uniform_discr( + min_corner=[-20, -20], max_corner=[20, 20], nsamples=[300, 300], + dtype='float32') + +# Make a parallel beam geometry with flat detector +# Angles: uniformly spaced, n = 360, min = 0, max = 2 * pi +angle_partition = odl.uniform_partition(0, 2 * np.pi, 360) +# Detector: uniformly sampled, n = 558, min = -30, max = 30 +detector_partition = odl.uniform_partition(-30, 30, 558) +geometry = odl.tomo.Parallel2dGeometry(angle_partition, detector_partition) + +# The implementation of the ray transform to use, options: +# 'scikit' Requires scikit-image (can be installed by +# running ``pip install scikit-image``). +# 'astra_cpu', 'astra_cuda' Require astra tomography to be installed. +# Astra is much faster than scikit. Webpage: +# https://github.com/astra-toolbox/astra-toolbox +impl = 'astra_cuda' + +# Ray transform aka forward projection. +ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl=impl) + +# Convert ray transform to proximal language operator +cvx_ray_trafo = odl.as_proximal_lang_operator(ray_trafo) + +# Create sinogram of forward projected phantom with noise +phantom = odl.phantom.shepp_logan(reco_space, modified=True) +phantom.show('phantom') +data = ray_trafo(phantom) +data += odl.phantom.white_noise(ray_trafo.range) * np.mean(data) * 0.1 +data.show('noisy data') + +# Convert to array +rhs_arr = data.asarray() + +# Set up optimization problem +# Note that proximal is not aware of the problem scaling in ODL, hence the norm +# does not match the norm in the space. +x = proximal.Variable(reco_space.shape) +funcs = [proximal.sum_squares(cvx_ray_trafo(x) - rhs_arr), + 0.2 * proximal.norm1(proximal.grad(x)), + proximal.nonneg(x), + proximal.nonneg(1 - x)] + +# Solve the problem +prob = proximal.Problem(funcs) +prob.solve(verbose=True) + +# Convert back to odl and display result +result_odl = reco_space.element(x.value) +result_odl.show('result') diff --git a/odl/operator/oputils.py b/odl/operator/oputils.py index a286676701e..b8c9a5b06d8 100644 --- a/odl/operator/oputils.py +++ b/odl/operator/oputils.py @@ -28,7 +28,8 @@ from odl.space.base_ntuples import FnBase from odl.space import ProductSpace -__all__ = ('matrix_representation', 'power_method_opnorm', 'as_scipy_operator') +__all__ = ('matrix_representation', 'power_method_opnorm', 'as_scipy_operator', + 'as_proximal_lang_operator') def matrix_representation(op): @@ -268,6 +269,51 @@ def rmatvec(v): rmatvec=rmatvec, dtype=dtype) + +def as_proximal_lang_operator(op, norm_bound=None): + """Wrap ``op`` as a ``proximal.BlackBox``. + + This is intended to be used with the proximal language solvers. + + Parameters + ---------- + op : `Operator` + A linear operator that should be wrapped. + norm_bound : float, optional + An upper bound on the spectral norm of the operator. Note that this is + the norm as defined by proximal, and hence use the unweighted spaces. + + Returns + ------- + ``proximal.BlackBox`` : proximal_lang_operator + The wrapped operator. + + Notes + ----- + If the data representation of ``op``'s domain and range is of type + `NumpyFn` this incurs no significant overhead. If the data type is `CudaFn` + or other nonlocal type, the overhead is significant. + + References + ---------- + For documentation on the proximal language (ProxImaL) see [Hei+2016]_. + """ + + # TODO: use out parameter once "as editable array" is added + + def forward(inp, out): + out[:] = op(inp).asarray() + + def adjoint(inp, out): + out[:] = op.adjoint(inp).asarray() + + import proximal + return proximal.LinOpFactory(input_shape=op.domain.shape, + output_shape=op.range.shape, + forward=forward, + adjoint=adjoint, + norm_bound=norm_bound) + if __name__ == '__main__': # pylint: disable=wrong-import-position from odl.util.testutils import run_doctests diff --git a/setup.py b/setup.py index cce869142da..aa387ed19c1 100644 --- a/setup.py +++ b/setup.py @@ -138,7 +138,8 @@ def run_tests(self): 'fftw': 'pyfftw', 'pywavelets': 'Pywavelets', 'scikit' : 'scikit-image', - 'all': test_requires + ['matplotlib', 'pyfftw', 'Pywavelets', 'scikit-image'] + 'proximal' : 'proximal', + 'all': test_requires + ['matplotlib', 'pyfftw', 'Pywavelets', 'scikit-image', 'proximal'] }, cmdclass={'test': PyTest}, From 10c5185b3074006254a1aa754cf6853a8bf2617f Mon Sep 17 00:00:00 2001 From: Jonas Adler Date: Tue, 16 Aug 2016 17:54:08 +0200 Subject: [PATCH 2/5] DOC: improve doc of proximal interface --- doc/source/guide/in_depth/index.rst | 1 + .../guide/in_depth/proximal_lang_guide.rst | 52 +++++++++++++++++++ examples/solvers/proximal_lang_poisson.py | 6 +-- examples/solvers/proximal_lang_tomography.py | 4 +- 4 files changed, 58 insertions(+), 5 deletions(-) create mode 100644 doc/source/guide/in_depth/proximal_lang_guide.rst diff --git a/doc/source/guide/in_depth/index.rst b/doc/source/guide/in_depth/index.rst index d5d421888df..3792e2f0d17 100644 --- a/doc/source/guide/in_depth/index.rst +++ b/doc/source/guide/in_depth/index.rst @@ -11,4 +11,5 @@ This is a more in depth guide to the different parts of ODL. linearspace_guide vectorization_guide numpy_guide + proximal_lang_guide chambolle_pock_guide diff --git a/doc/source/guide/in_depth/proximal_lang_guide.rst b/doc/source/guide/in_depth/proximal_lang_guide.rst new file mode 100644 index 00000000000..bad8aefb0d7 --- /dev/null +++ b/doc/source/guide/in_depth/proximal_lang_guide.rst @@ -0,0 +1,52 @@ +.. _proximal_lang_in_depth: + +####################### +Using ODL with ProxImaL +####################### + +`Proximal +`_ is a Python-embedded modeling language for image optimization problems and can be used with ODL to solve typical inverse problems phrased as optimization problems. The package is especially suited for non-differentiable problems such as total variance denoising. + +Here is a minimal example of solving Poisson's equation equation on an interval with a TV type regularizer (:math:`\min_x \ 10||-\Delta x - rhs||_2^2 + ||\nabla x||_1`):: + + >>> space = odl.uniform_discr(0, 1, 5) + >>> op = -odl.Laplacian(space) + >>> proximal_lang_op = odl.as_proximal_lang_operator(op) + >>> rhs = space.element(lambda x: (x>0.4) & (x<0.6)) # indicator function on [0.4, 0.6] + >>> x = proximal.Variable(space.shape) + >>> prob = proximal.Problem([10 * proximal.sum_squares(x - rhs.asarray()), + >>> proximal.norm1(proximal.grad(x))]) + >>> prob.solve() + >>> x.value + array([ 0.02352054, 0.02647946, 0.9 , 0.02647946, 0.02352054]) + +Notable differences between ODL and ProxImaL +============================================ + +It may be tempting to try to convert an arbitrary problem from ODL into ProxImaL, but some differences exist. + +Norms +----- +Norms in ODL are scaled according to the underlying function space. Hence a sequence of statements converging discretizations give rise to a converging norm:: + + >>> for n in range(2, 10000): + ... X = odl.uniform_discr(0, 1, n) + ... print(X.element(lambda x: x).norm()) + 0.559016994375 + 0.576628129734 + 0.577343052266 + 0.577350268468 + >>> 1 / np.sqrt(3) # exact result + 0.577350269189 + +this is not the case in proximal, where the norm depends on the number of discretization points. Hence a scaling that is correct for a problem in ODL needs not be correct in proximal. This also changes the definition of things like the operator norm. + +This also has the added effect of changing the definition of derived features, like the spectral norm of operators. + +Spaces +------ +ODL can represent some complicated spaces, like :math:`\mathbb{R}^3 \times \mathbb{C}^2` through the `ProductSpace` class:: + + >>> space = odl.ProductSpace(odl.rn(3), odl.cn(2)) + +This can then be used in solvers and other structures. \ No newline at end of file diff --git a/examples/solvers/proximal_lang_poisson.py b/examples/solvers/proximal_lang_poisson.py index bd29c7e3c2b..99a761bb653 100644 --- a/examples/solvers/proximal_lang_poisson.py +++ b/examples/solvers/proximal_lang_poisson.py @@ -47,15 +47,15 @@ rhs += odl.phantom.white_noise(space) * np.std(rhs) * 0.1 rhs.show('rhs') -# Convert laplacian to cvx operator -cvx_laplacian = odl.as_proximal_lang_operator(laplacian) +# Convert laplacian to ProxImaL operator +proximal_lang_laplacian = odl.as_proximal_lang_operator(laplacian) # Convert to array rhs_arr = rhs.asarray() # Set up optimization problem x = proximal.Variable(space.shape) -funcs = [10 * proximal.sum_squares(cvx_laplacian(x) - rhs_arr), +funcs = [10 * proximal.sum_squares(proximal_lang_laplacian(x) - rhs_arr), proximal.norm1(proximal.grad(x))] # Solve the problem diff --git a/examples/solvers/proximal_lang_tomography.py b/examples/solvers/proximal_lang_tomography.py index f14d91a600a..dad2d586cb1 100644 --- a/examples/solvers/proximal_lang_tomography.py +++ b/examples/solvers/proximal_lang_tomography.py @@ -63,7 +63,7 @@ ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl=impl) # Convert ray transform to proximal language operator -cvx_ray_trafo = odl.as_proximal_lang_operator(ray_trafo) +proximal_lang_ray_trafo = odl.as_proximal_lang_operator(ray_trafo) # Create sinogram of forward projected phantom with noise phantom = odl.phantom.shepp_logan(reco_space, modified=True) @@ -79,7 +79,7 @@ # Note that proximal is not aware of the problem scaling in ODL, hence the norm # does not match the norm in the space. x = proximal.Variable(reco_space.shape) -funcs = [proximal.sum_squares(cvx_ray_trafo(x) - rhs_arr), +funcs = [proximal.sum_squares(proximal_lang_ray_trafo(x) - rhs_arr), 0.2 * proximal.norm1(proximal.grad(x)), proximal.nonneg(x), proximal.nonneg(1 - x)] From 3ad710f9f2060baf98e46b64dd4ee0efccddb4f3 Mon Sep 17 00:00:00 2001 From: Jonas Adler Date: Wed, 17 Aug 2016 10:00:36 +0200 Subject: [PATCH 3/5] DOC: doc updates to proximal language interface --- doc/source/refs.rst | 2 +- examples/solvers/proximal_lang_poisson.py | 20 +++++++----------- examples/solvers/proximal_lang_tomography.py | 22 ++++++++------------ odl/operator/operator.py | 2 +- odl/operator/oputils.py | 10 +++++---- 5 files changed, 25 insertions(+), 31 deletions(-) diff --git a/doc/source/refs.rst b/doc/source/refs.rst index 248a4c79dba..260de325481 100644 --- a/doc/source/refs.rst +++ b/doc/source/refs.rst @@ -42,7 +42,7 @@ References .. [GNS2009] Griva, I, Nash, S G, and Sofer, A. *Linear and nonlinear optimization*. Siam, 2009. -.. [Hei+2016] Heide F. et. al. *ProxImaL: Efficient Image Optimization using +.. [Hei+2016] Heide, F et al. *ProxImaL: Efficient Image Optimization using Proximal Algorithms*. ACM Transactions on Graphics (TOG), 2016. .. [Kva1991] Kvaalen, E. *A faster Broyden method*. BIT Numerical diff --git a/examples/solvers/proximal_lang_poisson.py b/examples/solvers/proximal_lang_poisson.py index 99a761bb653..f100953aee9 100644 --- a/examples/solvers/proximal_lang_poisson.py +++ b/examples/solvers/proximal_lang_poisson.py @@ -15,34 +15,30 @@ # You should have received a copy of the GNU General Public License # along with ODL. If not, see . -"""Poissons problem using the proximal solver. +"""Poisson's problem using the ProxImaL solver. Solves the optimization problem min_x 10 ||laplacian(x) - g||_2^2 + || |grad(x)| ||_1 -Where ``laplacian`` is the spatial laplacian, ``grad`` the spatial +Where ``laplacian`` is the spatial Laplacian, ``grad`` the spatial gradient and ``g`` is given noisy data. """ -# Imports for common Python 2/3 codebase -from __future__ import print_function, division, absolute_import -from future import standard_library -standard_library.install_aliases() - import numpy as np import odl import proximal -# Create space, a square from [0, 0] to [100, 100] with (100 x 100) points +# Create space defined on a square from [0, 0] to [100, 100] with (100 x 100) +# points space = odl.uniform_discr([0, 0], [100, 100], [100, 100]) -# Create odl operator for laplacian +# Create ODL operator for the Laplacian laplacian = odl.Laplacian(space) # Create right hand side phantom = odl.phantom.shepp_logan(space, modified=True) -phantom.show('phantom') +phantom.show('original image') rhs = laplacian(phantom) rhs += odl.phantom.white_noise(space) * np.std(rhs) * 0.1 rhs.show('rhs') @@ -58,10 +54,10 @@ funcs = [10 * proximal.sum_squares(proximal_lang_laplacian(x) - rhs_arr), proximal.norm1(proximal.grad(x))] -# Solve the problem +# Solve the problem using ProxImaL prob = proximal.Problem(funcs) prob.solve(verbose=True) # Convert back to odl and display result result_odl = space.element(x.value) -result_odl.show('result') +result_odl.show('result from ProxImaL') diff --git a/examples/solvers/proximal_lang_tomography.py b/examples/solvers/proximal_lang_tomography.py index dad2d586cb1..0b90f69fc26 100644 --- a/examples/solvers/proximal_lang_tomography.py +++ b/examples/solvers/proximal_lang_tomography.py @@ -15,7 +15,7 @@ # You should have received a copy of the GNU General Public License # along with ODL. If not, see . -"""Total variation tomography using the proximal solver. +"""Tomography with TV regularization using the ProxImaL solver. Solves the optimization problem @@ -25,11 +25,6 @@ gradient and ``g`` is given noisy data. """ -# Imports for common Python 2/3 codebase -from __future__ import print_function, division, absolute_import -from future import standard_library -standard_library.install_aliases() - import numpy as np import odl import proximal @@ -54,12 +49,12 @@ # The implementation of the ray transform to use, options: # 'scikit' Requires scikit-image (can be installed by # running ``pip install scikit-image``). -# 'astra_cpu', 'astra_cuda' Require astra tomography to be installed. +# 'astra_cpu', 'astra_cuda' Requires astra tomography to be installed. # Astra is much faster than scikit. Webpage: # https://github.com/astra-toolbox/astra-toolbox impl = 'astra_cuda' -# Ray transform aka forward projection. +# Initialize the ray transform (forward projection). ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl=impl) # Convert ray transform to proximal language operator @@ -72,22 +67,23 @@ data += odl.phantom.white_noise(ray_trafo.range) * np.mean(data) * 0.1 data.show('noisy data') -# Convert to array +# Convert to array for ProxImaL rhs_arr = data.asarray() # Set up optimization problem -# Note that proximal is not aware of the problem scaling in ODL, hence the norm -# does not match the norm in the space. +# Note that proximal is not aware of the underlying space and only works with +# matrices. Hence the norm in proximal does not match the norm in the ODL space +# exactly. x = proximal.Variable(reco_space.shape) funcs = [proximal.sum_squares(proximal_lang_ray_trafo(x) - rhs_arr), 0.2 * proximal.norm1(proximal.grad(x)), proximal.nonneg(x), proximal.nonneg(1 - x)] -# Solve the problem +# Solve the problem using ProxImaL prob = proximal.Problem(funcs) prob.solve(verbose=True) # Convert back to odl and display result result_odl = reco_space.element(x.value) -result_odl.show('result') +result_odl.show('ProxImaL result') diff --git a/odl/operator/operator.py b/odl/operator/operator.py index 19e134dbc09..a05057a8e0a 100644 --- a/odl/operator/operator.py +++ b/odl/operator/operator.py @@ -538,7 +538,7 @@ def _call(self, x, out=None, **kwargs): pattern. See the `documentation - `_ + `_ for more info on in-place vs. out-of-place evaluation. Parameters diff --git a/odl/operator/oputils.py b/odl/operator/oputils.py index b8c9a5b06d8..2f717c0d25a 100644 --- a/odl/operator/oputils.py +++ b/odl/operator/oputils.py @@ -239,7 +239,7 @@ def as_scipy_operator(op): ----- If the data representation of ``op``'s domain and range is of type `NumpyFn` this incurs no significant overhead. If the data type is `CudaFn` - or other nonlocal type, the overhead is significant. + or some other nonlocal type, the overhead is significant. """ if not op.is_linear: raise ValueError('`op` needs to be linear') @@ -273,12 +273,14 @@ def rmatvec(v): def as_proximal_lang_operator(op, norm_bound=None): """Wrap ``op`` as a ``proximal.BlackBox``. - This is intended to be used with the proximal language solvers. + This is intended to be used with the `proximal language solvers. + `_ Parameters ---------- op : `Operator` - A linear operator that should be wrapped. + Linear operator to be wrapped. Its domain and range must implement + ``shape``, and elements in these need to implement ``asarray``. norm_bound : float, optional An upper bound on the spectral norm of the operator. Note that this is the norm as defined by proximal, and hence use the unweighted spaces. @@ -292,7 +294,7 @@ def as_proximal_lang_operator(op, norm_bound=None): ----- If the data representation of ``op``'s domain and range is of type `NumpyFn` this incurs no significant overhead. If the data type is `CudaFn` - or other nonlocal type, the overhead is significant. + or some other nonlocal type, the overhead is significant. References ---------- From 06e2d3f005ade0ac8dad7370a72820cb926f8c01 Mon Sep 17 00:00:00 2001 From: Jonas Adler Date: Wed, 17 Aug 2016 10:14:42 +0200 Subject: [PATCH 4/5] DOC: add release info for proximal interface --- doc/source/prs.inc | 1 + doc/source/release_notes.rst | 7 +++---- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/doc/source/prs.inc b/doc/source/prs.inc index 4f1b6afc016..bae88c4a593 100644 --- a/doc/source/prs.inc +++ b/doc/source/prs.inc @@ -58,6 +58,7 @@ .. _PR 491: https://github.com/odlgroup/odl/pull/491 .. _PR 492: https://github.com/odlgroup/odl/pull/492 .. _PR 493: https://github.com/odlgroup/odl/pull/493 +.. _PR 494: https://github.com/odlgroup/odl/pull/494 .. _PR 495: https://github.com/odlgroup/odl/pull/495 .. _PR 500: https://github.com/odlgroup/odl/pull/500 .. _PR 502: https://github.com/odlgroup/odl/pull/502 diff --git a/doc/source/release_notes.rst b/doc/source/release_notes.rst index 89772bf0d6a..715fa82ded6 100644 --- a/doc/source/release_notes.rst +++ b/doc/source/release_notes.rst @@ -11,8 +11,7 @@ Next release New features ------------ -- First implementation of a deformation model based on linearized deformations using displacement - vector fields. (`PR 488`) +- Add option to interface with ProxImaL solvers using ODL operators. (`PR 494`) ODL 0.3.1 Release Notes (2016-08-15) ==================================== @@ -23,7 +22,7 @@ splitting. New features ------------ -- New solvers based on the Douglas-Rachford and forward-backward splitting schemes. (`PR 478`_, +- New solvers based on the Douglas-Rachford and forward-backward splitting schemes. (`PR 478`_, `PR 480`_) - ``NormOperator`` and ``DistOperator`` added. (`PR 487`_) - Single-element ``NtuplesBase`` vectors can now be converted to ``float``, ``complex`` etc. @@ -38,7 +37,7 @@ Improvements (`PR 489`_, `PR 482`_, `PR 491`_) - Clearer separation between attributes that are intended as part of the subclassing API and those that are not. (`PR 471`_) -- Chambolle-Pock solver accepts also non-linear operators and has better documentation now. +- Chambolle-Pock solver accepts also non-linear operators and has better documentation now. (`PR 490`_) - Clean-up of imports. (`PR 492`_) - All solvers now check that the given start value ``x`` is in ``op.domain``. (`PR 502`_) From bd7da98c143575ef30805b8c992c808212c892ee Mon Sep 17 00:00:00 2001 From: Jonas Adler Date: Wed, 17 Aug 2016 11:00:33 +0200 Subject: [PATCH 5/5] DOC: rename proximal -> ProxImaL where appropriate --- doc/source/guide/glossary.rst | 15 +++++++++------ doc/source/guide/in_depth/proximal_lang_guide.rst | 4 ++-- odl/operator/oputils.py | 6 +++--- 3 files changed, 14 insertions(+), 11 deletions(-) diff --git a/doc/source/guide/glossary.rst b/doc/source/guide/glossary.rst index d6e4b679b07..7f84fd3cf9b 100644 --- a/doc/source/guide/glossary.rst +++ b/doc/source/guide/glossary.rst @@ -76,17 +76,20 @@ Glossary proximal Given a closed proper convex function :math:`f`, the proximal operator is defined by - + .. math:: - + \text{prox}_f(v) = \arg\min_x f(x) + (1/2)||x - v||_2^2 - + + proximal is also occationally used instead of ProxImaL, then refering to the proximal + modelling language for the solution of convex optimization problems. + proximal factory - A proximal factory associated with a function :math:`f` is a `callable`, + A proximal factory associated with a function :math:`f` is a `callable`, which returns the proximal of the scaled function :math:`\sigma f` when called with a - scalar :math:`\sigma`. This is used due to the fact that optimization methods often use + scalar :math:`\sigma`. This is used due to the fact that optimization methods often use :math:`\text{prox}_{\sigma f}` for varying :math:`\sigma`. - + range Set of elements to which an operator maps, i.e. in which the result of an operator evaluation lies. diff --git a/doc/source/guide/in_depth/proximal_lang_guide.rst b/doc/source/guide/in_depth/proximal_lang_guide.rst index bad8aefb0d7..3c34002eeff 100644 --- a/doc/source/guide/in_depth/proximal_lang_guide.rst +++ b/doc/source/guide/in_depth/proximal_lang_guide.rst @@ -39,7 +39,7 @@ Norms in ODL are scaled according to the underlying function space. Hence a sequ >>> 1 / np.sqrt(3) # exact result 0.577350269189 -this is not the case in proximal, where the norm depends on the number of discretization points. Hence a scaling that is correct for a problem in ODL needs not be correct in proximal. This also changes the definition of things like the operator norm. +this is not the case in ProxImaL, where the norm depends on the number of discretization points. Hence a scaling that is correct for a problem in ODL needs not be correct in proximal. This also changes the definition of things like the operator norm. This also has the added effect of changing the definition of derived features, like the spectral norm of operators. @@ -49,4 +49,4 @@ ODL can represent some complicated spaces, like :math:`\mathbb{R}^3 \times \math >>> space = odl.ProductSpace(odl.rn(3), odl.cn(2)) -This can then be used in solvers and other structures. \ No newline at end of file +This can then be used in solvers and other structures. ProxImaL currently lacks an equivalent structure. \ No newline at end of file diff --git a/odl/operator/oputils.py b/odl/operator/oputils.py index 2f717c0d25a..cfa4610938b 100644 --- a/odl/operator/oputils.py +++ b/odl/operator/oputils.py @@ -273,8 +273,8 @@ def rmatvec(v): def as_proximal_lang_operator(op, norm_bound=None): """Wrap ``op`` as a ``proximal.BlackBox``. - This is intended to be used with the `proximal language solvers. - `_ + This is intended to be used with the `ProxImaL language solvers. + `_ Parameters ---------- @@ -283,7 +283,7 @@ def as_proximal_lang_operator(op, norm_bound=None): ``shape``, and elements in these need to implement ``asarray``. norm_bound : float, optional An upper bound on the spectral norm of the operator. Note that this is - the norm as defined by proximal, and hence use the unweighted spaces. + the norm as defined by ProxImaL, and hence use the unweighted spaces. Returns -------