Issue1077 chambolle pock #1092

mehrhardt · 2017-08-22T15:57:38Z

I renamed Chambolle-Pock to PDHG / primal dual hybrid gradient in many files. I hope I have not overlooked one.

Moreover, I changed the meaning of "gamma" in PDHG. This is now replaced with gamma_primal. In addition, one can do dual acceleration with gamma_dual.

There are still 3 small issues that should be solved before any merging.

The documentation is updated, but I think there are two links that still contain the name "Chambolle-Pock" and I didnt know how to change that as I didnt know where exactly they are pointing to.
The gamma_primal part is tested in a new example of L2-TV denoising (ROF). Here I would like to compute the objective function values over the iterations to show that it gets faster with gamma_primal > 0. How would I do that? Moreover, I wanted to compare different versions (assignments to functionals f and g). Somehow they converge to different objective function values (which they should not ...).
A simple example where gamma_dual > 0 makes sense would be something with Huberized TV. Is this or something similar available in ODL?

adler-j · 2017-08-22T16:45:28Z

Try having a look at the failing tests, I'll review.

adler-j

Overall a good change I'd say. With that said the name is quite long, I'd probably prefer pdhg. You should also search the documentation for any residual "chambolle_pock" mentions (e.g. in "getting started"). Those all need to be updated.

Could also need @kohr-h input.

adler-j · 2017-08-22T16:46:59Z

doc/source/guide/primal_dual_hybrid_gradient_guide.rst

 #####################

-The `chambolle_pock_solver` was introduced in 2011 by Chambolle and Pock in the paper `A first-order primal-dual algorithm for convex problems with applications to imaging
+The `primal_dual_hybrid_gradient_solver` was studied in 2011 by Chambolle and Pock in the paper `A first-order primal-dual algorithm for convex problems with applications to imaging


I'd scrap _solver from the function name now. Too long as is.

Perhaps we could even go with pdhg?

adler-j · 2017-08-22T16:47:45Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

@@ -0,0 +1,124 @@
+"""Total variation denoising using the primal-dual hybrid gradient algorithm.


Not quite clear what this function adds? I thought we had something similar before.

The current version of this file runs 3 algorithms for the ROF (L2^2 + TV) problem. 1) the algorithm of the old example 2) a variant that most people would use for this and 3) the 1/k^2 convergence version that makes use of strong convexity.

It also shows how to compare algorithms with ODL. Thereby this add some features that other examples don't have.

Sounds good!

adler-j · 2017-08-22T16:48:26Z

examples/solvers/pdhg_denoising_ROF_primal_acceleration.py

+# Make separable sum of functionals, order must correspond to the operator K
+f = l1_norm
+
+# Data fit with non-negativity constraint


This could use a word on what is happening here.

This file is removed and the primal acceleration bit is on the algorithmic comparison

adler-j · 2017-08-22T16:49:04Z

examples/solvers/pdhg_denoising_ROF_primal_acceleration.py

+# Gradient operator
+gradient = odl.Gradient(space, method='forward')
+
+# Matrix of operators


No matrix to be seen here?

adler-j · 2017-08-22T16:49:43Z

examples/solvers/pdhg_denoising_ROF_primal_acceleration.py

+# Original image
+orig = space.element(image)
+
+# Add noise


Comments like this add very little. What you do is obvious from the code here

This is your or Holger's comment. I am happy to remove it :)

Feel free to remove old cruft :-)

adler-j · 2017-08-22T16:50:44Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py

+                             ''.format(gamma_dual_in))
+
+    if gamma_primal is not None and gamma_dual is not None:
+        raise ValueError('Only one acceleration parameter can be used')


Any reason for this? Also this should be mentioned in the function definition.

Done. Depending on which part is strongly convex, the step sizes get either multiplied or divided by theta. Both doesn't make sense.

mehrhardt · 2017-08-24T16:47:04Z

All changes are done and Travis is happy. @adler-j, do you want to have another look? Maybe @kohr-h wants to add anything?

mehrhardt · 2017-08-24T16:51:19Z

Before, I forget. The third point I raised initially about dual acceleration is still not met. I think it would be great if one looked at Huberized ROF and compares two algorithms: 1) PDHG with dual acceleration and 2) linearly convergent PDHG.

I won't be able to do that this week as this needs some non-trivial coding I suppose.

kohr-h

Looks (very) good. I have a few comments regarding documentation mostly, but otherwise this can go in for sure. Well done!

kohr-h · 2017-08-26T19:43:48Z

doc/source/release_notes.rst

@@ -220,7 +220,7 @@ Improvements
  * Major review of minor style issues. (:pull:`534`)
  * Typeset math in proximals. (:pull:`580`)

- Improved installation docs and update of Chambolle-Pock documentation. (:pull:`121`)
+- Improved installation docs and update of PDHG documentation. (:pull:`121`)


Please revert the changes in old release notes (I guess you used search-and-replace)

Agree with @kohr-h. We shouldnt edit old release notes

kohr-h · 2017-08-26T19:50:41Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+making use of the strong convexity of the problem.
+
+For further details and a description of the solution method used, see
+:ref:`PDHG` in the ODL documentation.


I'm not sure this link will work like this. (Well, currently all links are broken :-/, but that's a different story). Perhaps link to a URL in the documentation, you may need to derive the link address from the old one and the new name of the document.

Better now?

Looking good to me.

kohr-h · 2017-08-26T19:52:08Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+        self.obj_function_values_ergodic = []
+
+    def __call__(self, x):
+        # Fill in proper calls to functionals here


I guess you can remove this comment, now that the proper calls are used :-)

kohr-h · 2017-08-26T19:54:23Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

@@ -0,0 +1,221 @@
+"""Total variation denoising using PDHG.


I think this is a very well done and nicely documented example. Great!

kohr-h · 2017-08-26T19:54:48Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+reg_param = 0.3
+
+# l2-squared data matching
+factr = 1./reg_param * 0.5


Neither Travis nor Spyder complained about pep8 here. Can you be more specific what you think is violating pep8?

I guess he's comming at 1./reg_param missing spaces.

With that said, why not 0.5/reg_param?

Travis doesn't run the examples and thus also no PEP8 check. And indeed, Spyder doesn't complain.

autopep8 seems to do a better (more complete job) here. It's available in pip, and if you run it with the -d flag it gives you a nice diff with suggested changes. You can also run with -i to get the changes in-place.

kohr-h · 2017-08-26T20:08:37Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+# show images
+plt.figure(0)
+ax1 = plt.subplot(231)
+ax1.imshow(orig, clim=[0, 1])


Usually one needs a conversion from ODL's 'xy' axes to matplotlib's 'ij' axes. Probably the images end up correctly oriented like this, but if you do, e.g., orig.show() the image will be rotated by 90 degrees.

When wrapping a raw array representing an image in 'ij' convention, the conversion is np.rot90(img, -1), and the other way around it is np.rot90(odl_elem, 1).

It's up to you if you consider that useful here.

I would personally really suggest simply using space.element(img).show() here. Guaranteed to get everything correct

I would like to use img.show() but for this example, it is most important to create plot that are next to each other. Of course the user can drag them there but I really would not recommend that to anyone. As far as I know img.show() will always open a new figure, so it cannot be used here. Any ideas how to combine these two ideas?

Well you have a point there, and no, I don't have any idea on how to solve it. I guess this way is fine. One way would be to create a utility:

def convert_to_plt(x): # do rotations and stuff

kohr-h · 2017-08-26T20:11:29Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+obj_opt = min(obj_alg1 + obj_alg2 + obj_alg3)
+
+
+def rel_fun(x): return (np.array(x) - obj_opt)/(x[0] - obj_opt)


Function body on a new line.

kohr-h · 2017-08-26T20:13:52Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py

-    """Chambolle-Pock algorithm for non-smooth convex optimization problems.
+def pdhg(x, f, g, L, tau, sigma, niter, **kwargs):
+    """Primal-dual hybrid gradient algorithm for non-smooth convex optimization
+    problems.


This should fit on one line (perhaps simply remove "problems")

"Primal-dual hybrid gradient algorithm for convex optimization" is fine imo

kohr-h · 2017-08-26T20:19:03Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py

@@ -39,12 +40,12 @@ def chambolle_pock_solver(x, f, g, L, tau, sigma, niter, **kwargs):

    where ``L`` is an operator and ``F`` and ``G`` are functionals.


The F and G references should all be lowercase, could you please fix that while at it?

kohr-h · 2017-08-26T20:21:07Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py

-        ``F^*`` to be uniformly convex.
+        with ``tau`` and ``sigma`` as initial values. Requires ``G`` to be
+        strongly convex. ``gamma_primal`` is the strong convexity constant of
+        ``G``. Acceleration can either be done on the primal part or the dual


I don't quite understand this. It seems impossible that gamma_primal is both a derived property of G and a user choice.
Is the strong convexity constant maybe an upper bound for the choice of acceleration parameter?

True. In practice, one wants to have gamma as large as possible. Also the use of "strong convexity constant" is often not unique. Some authors refer to it as the largest constant such that property X holds, others as any constant such that X holds.

adler-j

Some minor comments, overall very nice!

adler-j · 2017-08-28T13:13:15Z

doc/source/release_notes.rst

@@ -220,7 +220,7 @@ Improvements
  * Major review of minor style issues. (:pull:`534`)
  * Typeset math in proximals. (:pull:`580`)

- Improved installation docs and update of Chambolle-Pock documentation. (:pull:`121`)
+- Improved installation docs and update of PDHG documentation. (:pull:`121`)


Agree with @kohr-h. We shouldnt edit old release notes

adler-j · 2017-08-28T13:14:56Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+
+# Add noise and convert to space element
+image += np.random.normal(0, 0.1, shape)
+noisy = space.element(image)


noisy = orig + 0.1 * odl.phantom.white_noise(space)

adler-j · 2017-08-28T13:16:03Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+# Run algorithms 2 and 3
+x_alg2 = x_start.copy()
+callback.reset()
+callback.callbacks[1].reset()


This should not be needed? Given the line above

adler-j · 2017-08-28T13:17:15Z

examples/solvers/pdhg_denoising_ROF_algorithm_comparison.py

+# show images
+plt.figure(0)
+ax1 = plt.subplot(231)
+ax1.imshow(orig, clim=[0, 1])


I would personally really suggest simply using space.element(img).show() here. Guaranteed to get everything correct

adler-j · 2017-08-28T13:18:33Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py

-    """Chambolle-Pock algorithm for non-smooth convex optimization problems.
+def pdhg(x, f, g, L, tau, sigma, niter, **kwargs):
+    """Primal-dual hybrid gradient algorithm for non-smooth convex optimization
+    problems.


"Primal-dual hybrid gradient algorithm for convex optimization" is fine imo

mehrhardt · 2017-08-30T08:35:37Z

Did the corrections. I would propose to merge it for now and we add the dual acceleration example at some other point. For instance when implementing Stochastic PDHG.

kohr-h · 2017-08-30T08:39:37Z

Did the corrections. I would propose to merge it for now and we add the dual acceleration example at some other point. For instance when implementing Stochastic PDHG.

I agree. Could you open a new issue for this so we have it on the agenda?

mehrhardt · 2017-08-30T08:44:28Z

Done, see #1101

kohr-h · 2017-08-30T10:12:27Z

I ran the examples and found some issues, PR coming up.

kohr-h · 2017-08-30T10:15:14Z

See mehrhardt#1

There's one issue left, and I don't know the source: The TGV tomography example craps out after 1 iteration. If you set niter = 1 the result is still OK (not NaN), but the next iteration makes everything NaN.
Any idea where this comes from? I haven't checked against the old version of the example, so I can't say if this PR introduces the issue.

mehrhardt · 2017-08-30T10:33:11Z

Interestingly, if I run the example the nan values sometimes occur and sometimes not. @kohr-h, can you confirm this or do you always get nan?

kohr-h · 2017-08-30T12:10:45Z

Interestingly, if I run the example the nan values sometimes occur and sometimes not. @kohr-h, can you confirm this or do you always get nan?

I'll check again, but if I remember correctly I ran it about 5 times (also varying tau and sigma`) and got NaN in the second iteration every time.

kohr-h · 2017-08-30T13:05:41Z

For me it seems to reliably fail. The reason are apparently the derivatives: after 1 iteration I get this (the white parts are NaN values):

kohr-h · 2017-08-30T13:08:16Z

But actually, on the master branch we have the same issue (I just checked), so this PR is not the cause.

So IMO it's fine that we just merge this and make a new issue for the TGV example.

mehrhardt · 2017-09-19T10:30:02Z

So, what needs to be done here? There seems to be a conflict in two files. How do we resolve this? I am happy to do it if you tell me how.

adler-j · 2017-09-19T10:31:25Z

Basically what you need to do is:

git checkout master
git pull
git checkout Issue1077_ChambollePock
git rebase master

Then you'll get some conflicts that you need to solve, and once that is done we can merge this.

adler-j · 2017-09-19T11:05:25Z

So something went very wrong in the rebase here :O

mehrhardt · 2017-09-19T11:11:13Z

Can you be more specific? What is it what you don't like?

mehrhardt · 2017-09-19T11:24:52Z

What I did was pretty much what you said:

git checkout upstream/master
git pull
git checkout Issue1077_ChambollePock
git rebase upstream/master

I resolved the conflicts by choosing my files and deleting some others. Then continued with the rebase and pushed the final result.

kohr-h · 2017-09-19T12:06:59Z

Seems to be somehow rebased on the wrong commit. There are old commits apparently added as new ones, which can be seen from the numerous items in history on this page.
You can undo the whole thing by running git reflog, look for the last thing before the git rebase, copy the hash of that commit and git checkout <last_good_commit_hash>. Then you can delete the botched branch git branch -D Issue1077_ChambollePock and git checkout -b Issue1077_ChambollePock to make a good version again.
Then start from that commit again.

It might be that the first part

git checkout upstream/master
git pull

didn't do what you expected it to do. It's good that you rebase on upstream/master, not your local one, since that's also a source of error. My recommended workflow is like this:

git fetch upstream  # updates the remote branches
git rebase upstream/master

Simpler, less error-prone.

mehrhardt · 2017-09-19T13:22:52Z

I did a bit different strategy:

	$ git reset --hard 053f4b251fb408d1944acf97c107b1b0fe6cb78d
	$ git rebase upstream/master

which seems to have done exactly what was needed. Is that OK? Otherwise I can always undo this since I have not pushed these changes. While my branch seems to be exactly what I wanted, I am a bit puzzled that I did not receive conflicts.

I would now push my changes to the remote branch and thereby getting rid of the noise above. Is that fine?

kohr-h · 2017-09-19T13:31:23Z

I did a bit different strategy:

$ git reset --hard 053f4b2
$ git rebase upstream/master

which seems to have done exactly what was needed. Is that OK? Otherwise I can always undo this since I have not pushed these changes. While my branch seems to be exactly what I wanted, I am a bit puzzled that I did not receive conflicts.

I would now push my changes to the remote branch and thereby getting rid of the noise above. Is that fine?

If that works for you, all good. It's possible that you didn't get conflicts because they could be auto-merged. If you check the messages during rebase and see something like auto-merge then that's it. Usually this happens when the same file has been edited upstream, but in a different place.

Just go ahead and push, there's nothing that can't be fixed 😉

mehrhardt · 2017-09-19T13:41:20Z

Yes, there were auto-merges. I double checked the changes of the current version of this pull request and it looks good to me. I could not find anything that was missing or incorrect.

kohr-h

One minor thing that slipped, after that we can merge

kohr-h · 2017-09-19T13:45:58Z

odl/solvers/nonsmooth/primal_dual_hybrid_gradient.py


-    where :math:`\|K\|` is the operator norm of :math:`K`.
+    where :math:`\|K\|` is the operator norm of :math:`L`.


Should be L in the first one, too

mehrhardt · 2017-09-19T13:52:56Z

Done. I found even another "K" that should have been an "L". :)

kohr-h · 2017-09-19T13:58:46Z

Done. I found even another "K" that should have been an "L". :)

Perfect. Just that nasty PEP8 check blew up :-)

kohr-h · 2017-09-19T14:11:12Z

Ok, I'll merge like this, no need to rebase on master again. Thanks @matthiasje!

adler-j · 2017-09-19T14:15:44Z

Now this is breaking changes for some downstream stuff (e.g. all my machine learning) that needs to be fixed.

We need to have a look at that.

kohr-h · 2017-09-19T14:44:06Z

Maybe you should use a pinned git version of ODL in the requirements.txt? 😛

adler-j reviewed Aug 22, 2017

View reviewed changes

kohr-h approved these changes Aug 26, 2017

View reviewed changes

adler-j reviewed Aug 28, 2017

View reviewed changes

mehrhardt added 9 commits September 19, 2017 14:12

STY: rename files from chambolle_pock* to pdhg*

53c2876

STY: rename pdhg.py to primal_dual_hybrid_gradient.py

2d912a6

rename Chambolle-Pock to PDHG inside files

1f268c5

TYPO: remove typo in filename

90c660c

BUG: correct acceleration in PDHG

f85f441

TYPO: in PDHG description

a5f02ae

STY: finish renaming of examples

6d7c2c3

STY: Update guide

ed5848e

Add two new examples

ba7de18

mehrhardt and others added 13 commits September 19, 2017 14:12

remove 'solver' from PDHG name

06ca665

remove white spaces

0497d05

correct indentation

8697df7

revert release notes

9e73228

STY: update pdhg example based on comments

c034aaa

STY: update pdhg based on comments

05b033c

STY: minor style corrections in pdhg example

2bdfd02

BUG: fix colorbar ticks when values are nan

07aebbe

DOC: fix bad PDHG convergence criterion and a doc link

88758e5

STY: pep8 fixes in examples

f9f7bec

MAINT: regularize more in deconv example with added noise

941d145

MAINT: plot images in grayscale in ROF comparison example

0db2977

MAINT: run more iterations in TGV denoising example

121131c

mehrhardt force-pushed the Issue1077_ChambollePock branch from 8057f78 to 121131c Compare September 19, 2017 13:35

kohr-h approved these changes Sep 19, 2017

View reviewed changes

TYPO: change name of operator from K to L

e5b5162

TYPO: remove empty line

89ea55b

kohr-h merged commit cc909ef into odlgroup:master Sep 19, 2017

mehrhardt deleted the Issue1077_ChambollePock branch September 19, 2017 14:13

This was referenced Sep 19, 2017

Change chambolle_pock_solver to pdhg chongchenmath/test#1

Open

Change chambolle_pock_solver to pdhg adler-j/learned_primal_dual#1

Closed

Change chambolle_pock_solver to pdhg adler-j/learned_gradient_tomography#1

Closed

		@@ -0,0 +1,124 @@
		"""Total variation denoising using the primal-dual hybrid gradient algorithm.

		obj_opt = min(obj_alg1 + obj_alg2 + obj_alg3)


		def rel_fun(x): return (np.array(x) - obj_opt)/(x[0] - obj_opt)

		@@ -39,12 +40,12 @@ def chambolle_pock_solver(x, f, g, L, tau, sigma, niter, **kwargs):

		where ``L`` is an operator and ``F`` and ``G`` are functionals.


		where :math:`\\|K\\|` is the operator norm of :math:`K`.
		where :math:`\\|K\\|` is the operator norm of :math:`L`.

Issue1077 chambolle pock #1092

Issue1077 chambolle pock #1092

Conversation

mehrhardt commented Aug 22, 2017

adler-j commented Aug 22, 2017

adler-j left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mehrhardt commented Aug 24, 2017

mehrhardt commented Aug 24, 2017

kohr-h left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

kohr-h Aug 26, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

adler-j Aug 30, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

adler-j left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mehrhardt commented Aug 30, 2017

kohr-h commented Aug 30, 2017

mehrhardt commented Aug 30, 2017

kohr-h commented Aug 30, 2017

kohr-h commented Aug 30, 2017

mehrhardt commented Aug 30, 2017

kohr-h commented Aug 30, 2017

kohr-h commented Aug 30, 2017

kohr-h commented Aug 30, 2017 • edited Loading

mehrhardt commented Sep 19, 2017

adler-j commented Sep 19, 2017

adler-j commented Sep 19, 2017

mehrhardt commented Sep 19, 2017

mehrhardt commented Sep 19, 2017 • edited Loading

kohr-h commented Sep 19, 2017

mehrhardt commented Sep 19, 2017

kohr-h commented Sep 19, 2017

mehrhardt commented Sep 19, 2017

kohr-h left a comment

Choose a reason for hiding this comment

kohr-h Aug 26, 2017 •

edited

Loading

adler-j Aug 30, 2017 •

edited

Loading

kohr-h commented Aug 30, 2017 •

edited

Loading

mehrhardt commented Sep 19, 2017 •

edited

Loading