|
15 | 15 | # You should have received a copy of the GNU General Public License
|
16 | 16 | # along with ODL. If not, see <http://www.gnu.org/licenses/>.
|
17 | 17 |
|
18 |
| -"""Basic examples on how to use the fucntional class together with solvers. |
| 18 | +"""Basic example on how to use the functional class together with solvers. |
19 | 19 |
|
20 | 20 | This file shows an example of how to set up and solve an optimization problem
|
21 | 21 | using the default functionals. The problem we will solve is to minimize
|
22 | 22 | 1/2 * ||x - g||_2^2 + lam*||x||_1, for some vector g and some constant lam,
|
23 | 23 | subject to that all components in x are greater than or equal to 0. The
|
24 | 24 | theoretical optimal solution to this problem is x = (g - lam)_+, where ( )+
|
25 |
| -denotes the positive part of (i.e., (z)_+ = z if z >= 0, 0 otherwise). |
| 25 | +denotes the positive part of the element, i.e., (z_i)_+ = z_i if z_i >= 0, and |
| 26 | +0 otherwise. |
26 | 27 | """
|
27 | 28 |
|
28 | 29 | import numpy as np
|
|
46 | 47 | l2_func = 1.0 / 2.0 * odl.solvers.L2NormSquare(space)
|
47 | 48 | trans_l2_func = l2_func.translate(g)
|
48 | 49 |
|
49 |
| -# The problem will be solved using the Chambolle-Pock algorithm. Here we create |
50 |
| -# the necessary proximal factories of the conjugate functionals (see the |
51 |
| -# Chambolle-Pock algorithm and examples on this for more information). |
| 50 | +# The problem will be solved using the forward-backward primal-dual algorithm. |
| 51 | +# In this setting we let f = nonnegativity contraint, g = l1-norm, and h = |
| 52 | +# the squared l2-norm. Here we create necessary proximal and gradient operators |
| 53 | +# from the functionals. |
| 54 | +prox_f = odl.solvers.proximal_nonnegativity(space) |
| 55 | +prox_cc_g = lam_l1_func.conjugate_functional.proximal |
| 56 | +L = odl.IdentityOperator(space) |
| 57 | +grad_h = trans_l2_func.gradient |
| 58 | + |
| 59 | +# Some solver parameters |
| 60 | +niter = 50 |
| 61 | +tau = 0.5 |
| 62 | +sigma = 0.5 |
| 63 | + |
| 64 | +# Starting point, and also updated inplace in the solver |
| 65 | +x_fbpd = space.zero() |
| 66 | + |
| 67 | +# Optional: pass callback objects to solver |
| 68 | +callback = (odl.solvers.CallbackPrintIteration()) |
| 69 | + |
| 70 | +# Run the algorithm |
| 71 | +odl.solvers.forward_backward_pd(x=x_fbpd, prox_f=prox_f, prox_cc_g=[prox_cc_g], |
| 72 | + L=[L], grad_h=grad_h, tau=tau, sigma=[sigma], |
| 73 | + niter=niter, callback=callback) |
| 74 | + |
| 75 | +print(x_fbpd.asarray()) |
| 76 | + |
| 77 | + |
| 78 | +# The problem can also be solved using, e.g., the Chambolle-Pock algorithm. |
| 79 | +# Here we create the necessary proximal factories of the conjugate functionals |
| 80 | +# (see the Chambolle-Pock algorithm and examples on this for more information). |
52 | 81 | prox_cc_l2 = trans_l2_func.conjugate_functional.proximal
|
53 | 82 | prox_cc_l1 = lam_l1_func.conjugate_functional.proximal
|
54 | 83 |
|
|
62 | 91 | # Create the proximal operator for the constraint
|
63 | 92 | proximal_primal = odl.solvers.proximal_nonnegativity(op.domain)
|
64 | 93 |
|
65 |
| -# The operator norm is 1, since only identity operators are used |
66 |
| -op_norm = 1 |
| 94 | +# The operator norm is sqrt(2), since only identity operators are used |
| 95 | +op_norm = np.sqrt(2) |
67 | 96 |
|
68 | 97 | # Some solver parameters
|
69 |
| -niter = 100 # Number of iterations |
| 98 | +niter = 50 # Number of iterations |
70 | 99 | tau = 1.0 / op_norm # Step size for the primal variable
|
71 | 100 | sigma = 1.0 / op_norm # Step size for the dual variable
|
72 | 101 |
|
73 | 102 | # Optional: pass callback objects to solver
|
74 | 103 | callback = (odl.solvers.CallbackPrintIteration())
|
75 | 104 |
|
76 | 105 | # Starting point, and also updated inplace in the solver
|
77 |
| -x = op.domain.zero() |
| 106 | +x_cp = op.domain.zero() |
78 | 107 |
|
79 | 108 | # Run the algorithm
|
80 |
| -odl.solvers.chambolle_pock_solver( |
81 |
| - op, x, tau=tau, sigma=sigma, proximal_primal=proximal_primal, |
82 |
| - proximal_dual=proximal_dual, niter=niter, callback=callback) |
| 109 | +odl.solvers.chambolle_pock_solver(op=op, x=x_cp, tau=tau, sigma=sigma, |
| 110 | + proximal_primal=proximal_primal, |
| 111 | + proximal_dual=proximal_dual, niter=niter, |
| 112 | + callback=callback) |
83 | 113 |
|
84 |
| -print(x.asarray()) |
| 114 | +print(x_cp.asarray()) |
0 commit comments