From e9e8b8b8cfaa5ffae46b6aac4114997f43760ff3 Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Thu, 8 Dec 2022 16:04:51 +0000 Subject: [PATCH 1/6] Pass beta and phi through constructor of g_proximal. Remove them from call in forward_backward in preparation for other (TF) g_proximal implementations. --- cpp/examples/forward_backward/inpainting.cc | 2 +- .../inpainting_credible_interval.cc | 2 +- .../forward_backward/inpainting_joint_map.cc | 2 +- cpp/sopt/forward_backward.h | 2 +- cpp/sopt/l1_g_proximal.h | 14 ++++++++++---- cpp/tests/forward_backward.cc | 2 +- cpp/tests/inpainting.cc | 2 +- 7 files changed, 16 insertions(+), 10 deletions(-) diff --git a/cpp/examples/forward_backward/inpainting.cc b/cpp/examples/forward_backward/inpainting.cc index 174d3f52a..1f9de8044 100644 --- a/cpp/examples/forward_backward/inpainting.cc +++ b/cpp/examples/forward_backward/inpainting.cc @@ -103,7 +103,7 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(false); + auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) .l1_proximal_itermax(50) diff --git a/cpp/examples/forward_backward/inpainting_credible_interval.cc b/cpp/examples/forward_backward/inpainting_credible_interval.cc index 502377be6..07e2fe4b1 100644 --- a/cpp/examples/forward_backward/inpainting_credible_interval.cc +++ b/cpp/examples/forward_backward/inpainting_credible_interval.cc @@ -105,7 +105,7 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(false); + auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) .l1_proximal_itermax(50) diff --git a/cpp/examples/forward_backward/inpainting_joint_map.cc b/cpp/examples/forward_backward/inpainting_joint_map.cc index 53b20f364..be2e91c03 100644 --- a/cpp/examples/forward_backward/inpainting_joint_map.cc +++ b/cpp/examples/forward_backward/inpainting_joint_map.cc @@ -104,7 +104,7 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(false); + auto gp = std::make_shared>(fb->beta(), fb->Phi(), false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) .l1_proximal_itermax(50) diff --git a/cpp/sopt/forward_backward.h b/cpp/sopt/forward_backward.h index 4cddc5c72..f998f96e4 100644 --- a/cpp/sopt/forward_backward.h +++ b/cpp/sopt/forward_backward.h @@ -231,7 +231,7 @@ void ForwardBackward::iteration_step(t_Vector &out, t_Vector &residual, t_Vector &z, const t_real lambda) const { p = out; f_gradient(z, residual); - g_proximal(out, gamma() * beta(), out - beta() / nu() * (Phi().adjoint() * z)); + g_proximal(out, gamma(), z); p = out + lambda * (out - p); residual = (Phi() * p) - target(); } diff --git a/cpp/sopt/l1_g_proximal.h b/cpp/sopt/l1_g_proximal.h index 8e2a51649..31a5e6b2e 100644 --- a/cpp/sopt/l1_g_proximal.h +++ b/cpp/sopt/l1_g_proximal.h @@ -35,8 +35,10 @@ class L1GProximal : public GProximal { // In the constructor we need to construct the private l1_proximal_ // object that contains the real implementation details. The tight_frame // parameter is required for internal logic in l1_proximal - L1GProximal(bool tight_frame = false) + L1GProximal(Real const &beta, t_LinearTransform const &Phi, bool tight_frame = false) : tight_frame_ (tight_frame), + beta_(beta), + Phi_(Phi), l1_proximal_() {} ~L1GProximal() {}; @@ -55,7 +57,9 @@ class L1GProximal : public GProximal { // Return g_proximal as a lambda function. Used in operator() in base class. t_Proximal proximal_function() const override { return [this](t_Vector &out, Real gamma, t_Vector const &x) { - this -> l1_proximal(out, gamma, x); + this -> l1_proximal(out, + gamma * beta_, + out - beta_ / l1_proximal_.nu() * (Phi_.adjoint() * x)); }; } @@ -78,10 +82,10 @@ class L1GProximal : public GProximal { // ~~~ #define SOPT_MACRO(VAR, TYPE) \ /** \brief Getter, forwards to l1_proximal **/ \ - TYPE const &l1_proximal_##VAR() const { return l1_proximal().VAR(); } \ + TYPE const &l1_proximal_##VAR() const { return l1_proximal_.VAR(); } \ /** \brief Setter, forwards to l1_proximal **/ \ L1GProximal &l1_proximal_##VAR(TYPE const ARG) { \ - l1_proximal().VAR(ARG); \ + l1_proximal_.VAR(ARG); \ return *this; \ } @@ -108,6 +112,8 @@ class L1GProximal : public GProximal { bool tight_frame_; proximal::L1 l1_proximal_; + Real const beta_; + t_LinearTransform const Phi_; // Helper functions for calling l1_proximal //! Calls l1 proximal operator, checking for real constraints and tight frame diff --git a/cpp/tests/forward_backward.cc b/cpp/tests/forward_backward.cc index 6f9509836..9e1fbc40b 100644 --- a/cpp/tests/forward_backward.cc +++ b/cpp/tests/forward_backward.cc @@ -83,7 +83,7 @@ TEST_CASE("Check type returned on setting variables") { CHECK(is_imaging_proximal_ref::value); // Test the types of the l1 g_proximal object separately - auto gp = std::make_shared>(false); + auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); CHECK(is_l1_g_proximal_refl1_proximal_tolerance(1e-2))>::value); CHECK(is_l1_g_proximal_refl1_proximal_nu(1))>::value); CHECK(is_l1_g_proximal_refl1_proximal_itermax(50))>::value); diff --git a/cpp/tests/inpainting.cc b/cpp/tests/inpainting.cc index 9137719b2..5088baf96 100644 --- a/cpp/tests/inpainting.cc +++ b/cpp/tests/inpainting.cc @@ -67,7 +67,7 @@ TEST_CASE("Inpainting"){ // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(false); + auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) .l1_proximal_itermax(50) From fe3eeafdedef669cae596ab1c2de8710aaa7ac78 Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Thu, 8 Dec 2022 16:35:45 +0000 Subject: [PATCH 2/6] Fix failing test --- cpp/tests/forward_backward.cc | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/cpp/tests/forward_backward.cc b/cpp/tests/forward_backward.cc index 9e1fbc40b..d28dde80b 100644 --- a/cpp/tests/forward_backward.cc +++ b/cpp/tests/forward_backward.cc @@ -25,13 +25,17 @@ sopt::t_int random_integer(sopt::t_int min, sopt::t_int max) { typedef sopt::t_real Scalar; typedef sopt::Vector t_Vector; typedef sopt::t_real t_real; +typedef sopt::LinearTransform t_LinearTransform; auto const N = 5; TEST_CASE("Forward Backward with ||x - x0||_2^2 function", "[fb]") { using namespace sopt; t_Vector const target0 = t_Vector::Random(N); - auto const g0 = [](t_Vector &out, const t_real gamma, const t_Vector &x) { - proximal::id(out, gamma, x); + t_real const beta = 0.2; + t_real const nu = 1.0; + t_LinearTransform const Phi = linear_transform_identity(); + auto const g0 = [=](t_Vector &out, const t_real gamma, const t_Vector &x) { + proximal::id(out, gamma * beta, out - beta / nu * (Phi.adjoint() * x)); }; auto const grad = [](t_Vector &out, const t_Vector &x) { out = x; }; const t_Vector x_guess = t_Vector::Random(target0.size()); @@ -45,8 +49,9 @@ TEST_CASE("Forward Backward with ||x - x0||_2^2 function", "[fb]") { auto const fb = algorithm::ForwardBackward(grad, g0, target0) .itermax(300) .gamma(0.1) - .beta(0.2) + .beta(beta) .is_converged(convergence); + auto const result = fb(std::make_tuple(x_guess, res)); CAPTURE(result.niters); CAPTURE(result.x); From d9ee4000520bc9f8e056e7de70137943db88a175 Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Wed, 14 Dec 2022 13:28:40 +0000 Subject: [PATCH 3/6] Change 3rd argument of g_proximal to similar to l1_proximal --- cpp/sopt/l2_forward_backward.h | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/cpp/sopt/l2_forward_backward.h b/cpp/sopt/l2_forward_backward.h index 2d60e092b..46564645d 100644 --- a/cpp/sopt/l2_forward_backward.h +++ b/cpp/sopt/l2_forward_backward.h @@ -259,9 +259,11 @@ typename L2ForwardBackward::Diagnostic L2ForwardBackward::operat Diagnostic result; auto const g_proximal = [this](t_Vector &out, Real gamma, t_Vector const &x) { if (this->l2_proximal_weights().size() > 1) - this->l2_proximal_weighted()(out, this->l2_proximal_weights() * gamma, x); + this->l2_proximal_weighted()(out, this->l2_proximal_weights() * gamma, + out - beta() / nu() * (Phi().adjoint() * x)); else - this->l2_proximal()(out, this->l2_proximal_weights()(0) * gamma, x); + this->l2_proximal()(out, this->l2_proximal_weights()(0) * gamma, + out - beta() / nu() * (Phi().adjoint() * x)); }; const Real sigma_factor = sigma() * sigma(); auto const f_gradient = [this, sigma_factor](t_Vector &out, t_Vector const &x) { From 674953223126c91cac7df097eacd83ae6f5091cc Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Thu, 5 Jan 2023 11:53:36 +0000 Subject: [PATCH 4/6] Set beta and Phi with setters instead of constructor --- .github/workflows/cmake.yml | 4 +-- cpp/examples/forward_backward/inpainting.cc | 6 +++-- .../inpainting_credible_interval.cc | 8 +++--- .../forward_backward/inpainting_joint_map.cc | 6 +++-- cpp/sopt/l1_g_proximal.h | 12 ++++----- cpp/sopt/l1_proximal.h | 26 ++++++++++++++++++- cpp/tests/forward_backward.cc | 4 ++- cpp/tests/inpainting.cc | 6 +++-- 8 files changed, 52 insertions(+), 20 deletions(-) diff --git a/.github/workflows/cmake.yml b/.github/workflows/cmake.yml index b33ed6de5..67c2a4507 100644 --- a/.github/workflows/cmake.yml +++ b/.github/workflows/cmake.yml @@ -87,14 +87,14 @@ jobs: uses: actions/setup-python@v4 with: python-version: '3.10' - + - name: Install Conan id: conan uses: turtlebrowser/get-conan@main - name: Conan version run: echo "${{ steps.conan.outputs.version }}" - + - name: Prepare ccache timestamp id: ccache_cache_timestamp run: echo "{date_and_time}={$(date +'%Y-%m-%d-%H;%M;%S')}" >> $GITHUB_OUTPUT diff --git a/cpp/examples/forward_backward/inpainting.cc b/cpp/examples/forward_backward/inpainting.cc index 1f9de8044..6636da358 100644 --- a/cpp/examples/forward_backward/inpainting.cc +++ b/cpp/examples/forward_backward/inpainting.cc @@ -103,14 +103,16 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); + auto gp = std::make_shared>(false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) + .l1_proximal_beta(fb.beta()) + .l1_proximal_Phi(fb.Phi()) .l1_proximal_itermax(50) .l1_proximal_positivity_constraint(true) .l1_proximal_real_constraint(true) .Psi(psi); - + // Once the properties are set, inject it into the ImagingForwardBackward object fb.g_proximal(gp); diff --git a/cpp/examples/forward_backward/inpainting_credible_interval.cc b/cpp/examples/forward_backward/inpainting_credible_interval.cc index 07e2fe4b1..7b3bdc707 100644 --- a/cpp/examples/forward_backward/inpainting_credible_interval.cc +++ b/cpp/examples/forward_backward/inpainting_credible_interval.cc @@ -105,17 +105,19 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); + auto gp = std::make_shared>(false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) + .l1_proximal_beta(fb.beta()) + .l1_proximal_Phi(fb.Phi()) .l1_proximal_itermax(50) .l1_proximal_positivity_constraint(true) .l1_proximal_real_constraint(true) .Psi(psi); - + // Once the properties are set, inject it into the ImagingForwardBackward object fb.g_proximal(gp); - + SOPT_HIGH_LOG("Starting Forward Backward"); // Alternatively, forward-backward can be called with a tuple (x, residual) as argument // Here, we default to (Φ^Ty/ν, ΦΦ^Ty/ν - y) diff --git a/cpp/examples/forward_backward/inpainting_joint_map.cc b/cpp/examples/forward_backward/inpainting_joint_map.cc index be2e91c03..af1693475 100644 --- a/cpp/examples/forward_backward/inpainting_joint_map.cc +++ b/cpp/examples/forward_backward/inpainting_joint_map.cc @@ -104,14 +104,16 @@ int main(int argc, char const **argv) { // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(fb->beta(), fb->Phi(), false); + auto gp = std::make_shared>(false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) + .l1_proximal_beta(fb->beta()) + .l1_proximal_Phi(fb->Phi()) .l1_proximal_itermax(50) .l1_proximal_positivity_constraint(true) .l1_proximal_real_constraint(true) .Psi(psi); - + // Once the properties are set, inject it into the ImagingForwardBackward object fb->g_proximal(gp); diff --git a/cpp/sopt/l1_g_proximal.h b/cpp/sopt/l1_g_proximal.h index ca5d3df6f..9847cb4a4 100644 --- a/cpp/sopt/l1_g_proximal.h +++ b/cpp/sopt/l1_g_proximal.h @@ -35,10 +35,8 @@ class L1GProximal : public GProximal { // In the constructor we need to construct the private l1_proximal_ // object that contains the real implementation details. The tight_frame // parameter is required for internal logic in l1_proximal - L1GProximal(Real const &beta, t_LinearTransform const &Phi, bool tight_frame = false) + L1GProximal(bool tight_frame = false) : tight_frame_ (tight_frame), - beta_(beta), - Phi_(Phi), l1_proximal_() {} ~L1GProximal() {}; @@ -58,8 +56,8 @@ class L1GProximal : public GProximal { t_Proximal proximal_function() const override { return [this](t_Vector &out, Real gamma, t_Vector const &x) { this -> l1_proximal(out, - gamma * beta_, - out - beta_ / l1_proximal_.nu() * (Phi_.adjoint() * x)); + gamma * l1_proximal_.beta(), + out - l1_proximal_.beta() / l1_proximal_.nu() * (l1_proximal_.Phi().adjoint() * x)); }; } @@ -101,7 +99,9 @@ class L1GProximal : public GProximal { SOPT_MACRO(positivity_constraint, bool); SOPT_MACRO(real_constraint, bool); SOPT_MACRO(nu, Real); + SOPT_MACRO(beta, Real); SOPT_MACRO(weights, t_Vector); + SOPT_MACRO(Phi, t_LinearTransform); #undef SOPT_MACRO //! Analysis operator Ψ @@ -116,8 +116,6 @@ class L1GProximal : public GProximal { bool tight_frame_; proximal::L1 l1_proximal_; - Real const beta_; - t_LinearTransform const Phi_; // Helper functions for calling l1_proximal //! Calls l1 proximal operator, checking for real constraints diff --git a/cpp/sopt/l1_proximal.h b/cpp/sopt/l1_proximal.h index b3dd7c20c..fdce3c38a 100644 --- a/cpp/sopt/l1_proximal.h +++ b/cpp/sopt/l1_proximal.h @@ -40,12 +40,16 @@ class L1TightFrame { mpi::Communicator const &adjoint_comm = mpi::Communicator()) : Psi_(linear_transform_identity()), nu_(1e0), + Phi_(linear_transform_identity()), + beta_(1e0), direct_space_comm_(direct_comm), adjoint_space_comm_(adjoint_comm), weights_(Vector::Ones(1)) {} #else L1TightFrame() - : Psi_(linear_transform_identity()), nu_(1e0), weights_(Vector::Ones(1)) {} + : Psi_(linear_transform_identity()), nu_(1e0), + Phi_(linear_transform_identity()), beta_(1e0), + weights_(Vector::Ones(1)) {} #endif #define SOPT_MACRO(NAME, TYPE) \ @@ -63,6 +67,8 @@ class L1TightFrame { SOPT_MACRO(Psi, LinearTransform>); //! Bound on the squared norm of the operator Ψ SOPT_MACRO(nu, Real); + SOPT_MACRO(Phi, LinearTransform>); + SOPT_MACRO(beta, Real); #ifdef SOPT_MPI //! Communicator for summing in direct space (input when applying Psi) SOPT_MACRO(direct_space_comm, mpi::Communicator); @@ -307,6 +313,24 @@ class L1 : protected L1TightFrame { return *this; } + //! Bounds on the squared norm of the operator Ψ + Real beta() const { return L1TightFrame::beta(); } + //! Sets the bound on the squared norm of the operator Ψ + L1 &beta(Real const &beta) { + L1TightFrame::beta(beta); + return *this; + } + + //! Linear transform applied to input prior to L1 norm + LinearTransform> const &Phi() const { return L1TightFrame::Phi(); } + //! Set Ψ and Ψ^† using a matrix + template + typename std::enable_if= 1, L1 &>::type Phi(ARGS &&... args) { + L1TightFrame::Phi(std::forward(args)...); + return *this; + } + + //! \brief Special case if Ψ ia a tight frame. //! \see L1TightFrame template diff --git a/cpp/tests/forward_backward.cc b/cpp/tests/forward_backward.cc index d28dde80b..12c2baca3 100644 --- a/cpp/tests/forward_backward.cc +++ b/cpp/tests/forward_backward.cc @@ -88,9 +88,11 @@ TEST_CASE("Check type returned on setting variables") { CHECK(is_imaging_proximal_ref::value); // Test the types of the l1 g_proximal object separately - auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); + auto gp = std::make_shared>(false); CHECK(is_l1_g_proximal_refl1_proximal_tolerance(1e-2))>::value); CHECK(is_l1_g_proximal_refl1_proximal_nu(1))>::value); + CHECK(is_l1_g_proximal_refl1_proximal_beta(1))>::value); + CHECK(is_l1_g_proximal_refl1_proximal_Phi(linear_transform_identity()))>::value); CHECK(is_l1_g_proximal_refl1_proximal_itermax(50))>::value); CHECK(is_l1_g_proximal_refl1_proximal_positivity_constraint(true))>::value); CHECK(is_l1_g_proximal_refl1_proximal_real_constraint(true))>::value); diff --git a/cpp/tests/inpainting.cc b/cpp/tests/inpainting.cc index 5088baf96..b3d09c4a5 100644 --- a/cpp/tests/inpainting.cc +++ b/cpp/tests/inpainting.cc @@ -67,14 +67,16 @@ TEST_CASE("Inpainting"){ // Create a shared pointer to an instance of the L1GProximal class // and set its properties - auto gp = std::make_shared>(fb.beta(), fb.Phi(), false); + auto gp = std::make_shared>(false); gp->l1_proximal_tolerance(1e-4) .l1_proximal_nu(1) + .l1_proximal_beta(fb.beta()) + .l1_proximal_Phi(fb.Phi()) .l1_proximal_itermax(50) .l1_proximal_positivity_constraint(true) .l1_proximal_real_constraint(true) .Psi(psi); - + // Once the properties are set, inject it into the ImagingForwardBackward object fb.g_proximal(gp); From 23c7d3db586dd69fe6554a70c7c1f88fe1e80994 Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Thu, 5 Jan 2023 12:16:15 +0000 Subject: [PATCH 5/6] Workaround for homebrew issue --- .github/workflows/cmake.yml | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.github/workflows/cmake.yml b/.github/workflows/cmake.yml index 67c2a4507..fe7b1c372 100644 --- a/.github/workflows/cmake.yml +++ b/.github/workflows/cmake.yml @@ -77,7 +77,8 @@ jobs: - name: Install Dependencies on MacOS if: ${{ contains(matrix.os, 'macos') }} run: | - brew update + # Skip brew update until https://github.com/actions/setup-python/issues/577 is fixed + # brew update brew install libtiff open-mpi libyaml ccache echo "CMAKE_PREFIX_PATH=/usr/local/opt/libomp" >> $GITHUB_ENV From d001ae164bb9abc228edb5d5760c191e0e4970c4 Mon Sep 17 00:00:00 2001 From: Tuomas Koskela Date: Fri, 6 Jan 2023 13:49:34 +0000 Subject: [PATCH 6/6] Changes suggested by code review --- cpp/sopt/l1_proximal.h | 4 +++- cpp/tests/forward_backward.cc | 6 ++++-- 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/cpp/sopt/l1_proximal.h b/cpp/sopt/l1_proximal.h index fdce3c38a..3c39d103e 100644 --- a/cpp/sopt/l1_proximal.h +++ b/cpp/sopt/l1_proximal.h @@ -66,8 +66,10 @@ class L1TightFrame { //! Linear transform applied to input prior to L1 norm SOPT_MACRO(Psi, LinearTransform>); //! Bound on the squared norm of the operator Ψ - SOPT_MACRO(nu, Real); + SOPT_MACRO(nu, Real) + //! Measurement operator, used by g_proximal SOPT_MACRO(Phi, LinearTransform>); + //! β parameter, used by g_proximal SOPT_MACRO(beta, Real); #ifdef SOPT_MPI //! Communicator for summing in direct space (input when applying Psi) diff --git a/cpp/tests/forward_backward.cc b/cpp/tests/forward_backward.cc index 12c2baca3..0fa120414 100644 --- a/cpp/tests/forward_backward.cc +++ b/cpp/tests/forward_backward.cc @@ -33,6 +33,8 @@ TEST_CASE("Forward Backward with ||x - x0||_2^2 function", "[fb]") { t_Vector const target0 = t_Vector::Random(N); t_real const beta = 0.2; t_real const nu = 1.0; + t_real const gamma = 0.1; + int const itermax = 300; t_LinearTransform const Phi = linear_transform_identity(); auto const g0 = [=](t_Vector &out, const t_real gamma, const t_Vector &x) { proximal::id(out, gamma * beta, out - beta / nu * (Phi.adjoint() * x)); @@ -47,8 +49,8 @@ TEST_CASE("Forward Backward with ||x - x0||_2^2 function", "[fb]") { CAPTURE(x_guess); CAPTURE(res); auto const fb = algorithm::ForwardBackward(grad, g0, target0) - .itermax(300) - .gamma(0.1) + .itermax(itermax) + .gamma(gamma) .beta(beta) .is_converged(convergence);