ENH/DOC: updates in functionals and default_ops, see full commit message.

functional_basic_example:
- update comments

functional_basic_example_solver:
- update comments
- remove example when solving with Chambolle-Pock

default_ops:
- update doc in ConstantOperator
- update input checks ConstantOperator__init__

default_functionals:
- update doc in sevaral methods in several functionals

functional:
- update doc in some methods of Functional
- new imlpementation of derivative in Functional
- update doc in several other places

Author: aringh

examples/solvers/functional_basic_example.py

@@ -45,11 +45,11 @@
 import numpy as np
 import odl
 
+
 # Here we define the functional
 class MyFunctional(odl.solvers.Functional):
     """This is my functional: ||x||_2^2 + <x, y>."""
 
-    # Defining the __init__ function
     def __init__(self, domain, y):
         # This comand calls the init of Functional and sets a number of
         # parameters associated with a functional. All but domain have default
@@ -63,8 +63,7 @@ def __init__(self, domain, y):
             raise TypeError('y is not in the domain!')
         self._y = y
 
-    # Now we define a propert which returns y, so that the user can see which
-    # value is used in a particular instance of the class.
+    # Property that returns the linear term.
     @property
     def y(self):
         return self._y
@@ -76,35 +75,36 @@ def _call(self, x):
     # Next we define the gradient. Note that this is a property.
     @property
     def gradient(self):
-        # Inside this property, we define the gradient operator. This can be
-        # defined anywhere and just returned here, but in this example we will
-        # also define it here.
 
+        # The class corresponding to the gradient operator.
         class MyGradientOperator(odl.Operator):
+            """Class that implements the gradient operator of the functional
+            ``||x||_2^2 + <x,y>``.
+            """
 
-            # Define an __init__ method for this operator
             def __init__(self, functional):
                 super().__init__(domain=functional.domain,
                                  range=functional.domain)
 
                 self._functional = functional
 
-            # Define a _call method for this operator
             def _call(self, x):
                 return 2.0 * x + self._functional.y
 
         return MyGradientOperator(functional=self)
 
-    # Next we define the convex conjugate functional
+    # Next we define the convex conjugate functional.
     @property
     def conjugate_functional(self):
-        # This functional is implemented below
+        # This functional is implemented below.
         return MyFunctionalConjugate(domain=self.domain, y=self.y)
 
 
-# Here is the conjugate functional
+# Here is the conjugate functional.
 class MyFunctionalConjugate(odl.solvers.Functional):
-    """Calculations give that this funtional has the analytic expression
+    """Conjugate functional to ``||x||_2^2 + <x,y>``.
+
+    Calculations give that this funtional has the analytic expression
     f^*(x) = ||x||^2/2 - ||x-y||^2/4 + ||y||^2/2 - <x,y>.
     """
     def __init__(self, domain, y):

examples/solvers/functional_basic_example_solver.py

@@ -21,10 +21,8 @@
 using the default functionals. The problem we will solve is to minimize
 1/2 * ||x - g||_2^2 + lam*||x||_1, for some vector g and some constant lam,
 subject to that all components in x are greater than or equal to 0. The
-theoretical optimal solution to this problem is x = (g - lam)_+, where ( )+
-denotes the positive part of the element, i.e., (z_i)_+ = z_i if z_i >= 0, and
-0 otherwise.
-"""
+theoretical optimal solution to this problem is x = (g - lam)_+, where ( )_+
+denotes the positive part of the element, i.e., (z_i)_+ = max(z_i, 0)."""
 
 import numpy as np
 import odl
@@ -33,12 +31,12 @@
 n = 10
 space = odl.rn(n)
 
-# Create parameters. First half of g are ones, second half are minus ones.
+# Create parameters.
 g = space.element(np.hstack((np.ones(n/2), -np.ones(n/2))))
 lam = 0.5
 
 # Note that with the values above, the optimal solution is given by a vector
-# with fires half of the elements equal to 0.5, the second half equal to 0.
+# with first half of the elements equal to 0.5, the second half equal to 0.
 
 # Create the L1-norm functional and multiplyit with the constant lam.
 lam_l1_func = lam * odl.solvers.L1Norm(space)
@@ -62,53 +60,15 @@
 sigma = 0.5
 
 # Starting point, and also updated inplace in the solver
-x_fbpd = space.zero()
+x = space.element(np.random.randn(n))
+print('Initial guess: x = {}'.format(x.asarray()))
 
 # Optional: pass callback objects to solver
 callback = (odl.solvers.CallbackPrintIteration())
 
 # Run the algorithm
-odl.solvers.forward_backward_pd(x=x_fbpd, prox_f=prox_f, prox_cc_g=[prox_cc_g],
+odl.solvers.forward_backward_pd(x=x, prox_f=prox_f, prox_cc_g=[prox_cc_g],
                                 L=[L], grad_h=grad_h, tau=tau, sigma=[sigma],
                                 niter=niter, callback=callback)
 
-print(x_fbpd.asarray())
-
-
-# The problem can also be solved using, e.g., the Chambolle-Pock algorithm.
-# Here we create the necessary proximal factories of the conjugate functionals
-# (see the Chambolle-Pock algorithm and examples on this for more information).
-prox_cc_l2 = trans_l2_func.conjugate_functional.proximal
-prox_cc_l1 = lam_l1_func.conjugate_functional.proximal
-
-# Combined the proximals for use in the solver
-proximal_dual = odl.solvers.combine_proximals(prox_cc_l2, prox_cc_l1)
-
-# Create the matrix of operators for the Chambolle-Pock solver
-op = odl.BroadcastOperator(odl.IdentityOperator(space),
-                           odl.IdentityOperator(space))
-
-# Create the proximal operator for the constraint
-proximal_primal = odl.solvers.proximal_nonnegativity(op.domain)
-
-# The operator norm is sqrt(2), since only identity operators are used
-op_norm = np.sqrt(2)
-
-# Some solver parameters
-niter = 50  # Number of iterations
-tau = 1.0 / op_norm  # Step size for the primal variable
-sigma = 1.0 / op_norm  # Step size for the dual variable
-
-# Optional: pass callback objects to solver
-callback = (odl.solvers.CallbackPrintIteration())
-
-# Starting point, and also updated inplace in the solver
-x_cp = op.domain.zero()
-
-# Run the algorithm
-odl.solvers.chambolle_pock_solver(op=op, x=x_cp, tau=tau, sigma=sigma,
-                                  proximal_primal=proximal_primal,
-                                  proximal_dual=proximal_dual, niter=niter,
-                                  callback=callback)
-
-print(x_cp.asarray())
+print('Solution guess: x = {}'.format(x.asarray()))

odl/operator/default_ops.py

@@ -826,24 +826,22 @@ class ConstantOperator(Operator):
     """
 
     def __init__(self, vector, domain=None, range=None):
-        """Initialize an instance.
+        """Initialize a new instance.
 
         Parameters
         ----------
         vector : `LinearSpaceVector`
             The vector constant to be returned
         domain : `LinearSpace`, optional
-            default is vector.space
-            The domain of the operator.
+            Domain of the operator. Default: ``vector.space``
         range : `LinearSpace`, optional
-            default : vector.space
-            The range of the operator.
+            Range of the operator. Default: ``vector.space``
         """
-        if domain is None or range is None:
-            if not isinstance(vector, LinearSpaceVector):
-                raise TypeError('If either domain or range is unspecified, '
-                                '`vector` {!r} has to be a LinearSpaceVector '
-                                'instance'.format(vector))
+        if ((domain is None or range is None) and
+                not isinstance(vector, LinearSpaceVector)):
+                raise TypeError('If either domain or range is unspecified '
+                                '`vector` must be LinearSpaceVector, got {!r}.'
+                                ''.format(vector))
 
         if domain is None:
             domain = vector.space