Skip to content

Rename nu and gamma in padmm examples #377

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 2 commits into from
Jan 13, 2025
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
9 changes: 5 additions & 4 deletions cpp/example/padmm_mpi_random_coverage.cc
Original file line number Diff line number Diff line change
Expand Up @@ -86,19 +86,20 @@ std::shared_ptr<sopt::algorithm::ImagingProximalADMM<t_complex>> padmm_factory(
auto const epsilon = utilities::calculate_l2_radius(uv_data.vis.size(), sigma);
#endif
PURIFY_LOW_LOG("SARA Size = {}, Rank = {}", sara.size(), comm.rank());
const t_real gamma =
const t_real regulariser_strength =
utilities::step_size(uv_data.vis, measurements,
std::make_shared<sopt::LinearTransform<Vector<t_complex>> const>(Psi),
sara.size()) *
1e-3;
PURIFY_LOW_LOG("Epsilon {}, Rank = {}", epsilon, comm.rank());
PURIFY_LOW_LOG("Gamma {}, SARA Size = {}, Rank = {}", gamma, sara.size(), comm.rank());
PURIFY_LOW_LOG("Regulariser_Strength {}, SARA Size = {}, Rank = {}", regulariser_strength,
sara.size(), comm.rank());

// shared pointer because the convergence function need access to some data that we would rather
// not reproduce. E.g. padmm definition is self-referential.
auto padmm = std::make_shared<sopt::algorithm::ImagingProximalADMM<t_complex>>(uv_data.vis);
padmm->itermax(50)
.gamma(comm.all_reduce(gamma, MPI_MAX))
.regulariser_strength(comm.all_reduce(regulariser_strength, MPI_MAX))
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon)
#if PURIFY_PADMM_ALGORITHM == 2
Expand All @@ -115,7 +116,7 @@ std::shared_ptr<sopt::algorithm::ImagingProximalADMM<t_complex>> padmm_factory(
.l1_proximal_real_constraint(true)
.residual_tolerance(epsilon)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(*measurements);
sopt::ScalarRelativeVariation<t_complex> conv(padmm->relative_variation(),
Expand Down
21 changes: 12 additions & 9 deletions cpp/example/padmm_mpi_real_data.cc
Original file line number Diff line number Diff line change
Expand Up @@ -60,19 +60,19 @@ std::shared_ptr<sopt::algorithm::ImagingProximalADMM<t_complex>> padmm_factory(
#elif PURIFY_PADMM_ALGORITHM == 3 || PURIFY_PADMM_ALGORITHM == 1
auto const epsilon = 3 * std::sqrt(2 * uv_data.size()) * sigma;
#endif
const t_real gamma =
const t_real regulariser_strength =
utilities::step_size(uv_data.vis, measurements,
std::make_shared<sopt::LinearTransform<Vector<t_complex>> const>(Psi),
sara.size()) *
1e-3;
PURIFY_MEDIUM_LOG("Epsilon {}", epsilon);
PURIFY_MEDIUM_LOG("Gamma {}", gamma);
PURIFY_MEDIUM_LOG("Regulariser_Strength {}", regulariser_strength);

// shared pointer because the convergence function need access to some data that we would rather
// not reproduce. E.g. padmm definition is self-referential.
auto padmm = std::make_shared<sopt::algorithm::ImagingProximalADMM<t_complex>>(uv_data.vis);
padmm->itermax(50)
.gamma(comm.all_reduce<t_real>(gamma, MPI_MAX))
.regulariser_strength(comm.all_reduce<t_real>(regulariser_strength, MPI_MAX))
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon)
#if PURIFY_PADMM_ALGORITHM == 2
Expand All @@ -89,7 +89,7 @@ std::shared_ptr<sopt::algorithm::ImagingProximalADMM<t_complex>> padmm_factory(
.l1_proximal_real_constraint(true)
.residual_tolerance(epsilon)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(*measurements);
sopt::ScalarRelativeVariation<t_complex> conv(padmm->relative_variation(),
Expand Down Expand Up @@ -127,13 +127,16 @@ std::shared_ptr<sopt::algorithm::ImagingProximalADMM<t_complex>> padmm_factory(
const auto algo_update = [uv_data, imsizex, imsizey, padmm_weak, iter,
comm](const Vector<t_complex> &x) -> bool {
auto padmm = padmm_weak.lock();
if (comm.is_root()) PURIFY_MEDIUM_LOG("Step size γ {}", padmm->gamma());
if (comm.is_root()) PURIFY_MEDIUM_LOG("Step size γ {}", padmm->regulariser_strength());
*iter = *iter + 1;
Vector<t_complex> const alpha = padmm->Psi().adjoint() * x;
const t_real new_gamma = comm.all_reduce(alpha.real().cwiseAbs().maxCoeff(), MPI_MAX) * 1e-3;
if (comm.is_root()) PURIFY_MEDIUM_LOG("Step size γ update {}", new_gamma);
padmm->gamma(((std::abs(padmm->gamma() - new_gamma) > 0.2) and *iter < 200) ? new_gamma
: padmm->gamma());
const t_real new_regulariser_strength =
comm.all_reduce(alpha.real().cwiseAbs().maxCoeff(), MPI_MAX) * 1e-3;
if (comm.is_root()) PURIFY_MEDIUM_LOG("Step size γ update {}", new_regulariser_strength);
padmm->regulariser_strength(
((std::abs(padmm->regulariser_strength() - new_regulariser_strength) > 0.2) and *iter < 200)
? new_regulariser_strength
: padmm->regulariser_strength());
// updating parameter

Vector<t_complex> const residual = padmm->Phi().adjoint() * (uv_data.vis - padmm->Phi() * x);
Expand Down
11 changes: 6 additions & 5 deletions cpp/example/padmm_random_coverage.cc
Original file line number Diff line number Diff line change
Expand Up @@ -90,10 +90,11 @@ void padmm(const std::string &name, const Image<t_complex> &M31, const std::stri
auto const padmm =
sopt::algorithm::ImagingProximalADMM<t_complex>(uv_data.vis)
.itermax(500)
.gamma((Psi.adjoint() * (measurements_transform->adjoint() * uv_data.vis).eval())
.cwiseAbs()
.maxCoeff() *
1e-3)
.regulariser_strength(
(Psi.adjoint() * (measurements_transform->adjoint() * uv_data.vis).eval())
.cwiseAbs()
.maxCoeff() *
1e-3)
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon)
.tight_frame(false)
Expand All @@ -107,7 +108,7 @@ void padmm(const std::string &name, const Image<t_complex> &M31, const std::stri
#ifdef PURIFY_CImg
.is_converged(show_image)
#endif
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(*measurements_transform);

Expand Down
19 changes: 11 additions & 8 deletions cpp/example/padmm_real_data.cc
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,8 @@ void padmm(const std::string &name, const t_uint &imsizex, const t_uint &imsizey
pfitsio::write2d(Image<t_real>::Map(dimage.data(), imsizey, imsizex), dirty_image_fits);
pfitsio::write2d(Image<t_real>::Map(psf.data(), imsizey, imsizex), psf_image_fits);
auto const epsilon = 3 * std::sqrt(2 * uv_data.size()) * sigma;
auto const gamma = (measurements_transform->adjoint() * uv_data.vis).real().maxCoeff() * 1e-3;
auto const regulariser_strength =
(measurements_transform->adjoint() * uv_data.vis).real().maxCoeff() * 1e-3;
PURIFY_HIGH_LOG("Using epsilon of {}", epsilon);
#ifdef PURIFY_CImg
auto const canvas = std::make_shared<CDisplay>(
Expand All @@ -79,7 +80,7 @@ void padmm(const std::string &name, const t_uint &imsizex, const t_uint &imsizey
#endif
auto padmm = std::make_shared<sopt::algorithm::ImagingProximalADMM<t_complex>>(uv_data.vis);
padmm->itermax(500)
.gamma(gamma)
.regulariser_strength(regulariser_strength)
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon)
.tight_frame(false)
Expand All @@ -90,7 +91,7 @@ void padmm(const std::string &name, const t_uint &imsizex, const t_uint &imsizey
.l1_proximal_real_constraint(true)
.residual_convergence(epsilon)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(*measurements_transform);

Expand All @@ -101,14 +102,16 @@ void padmm(const std::string &name, const t_uint &imsizex, const t_uint &imsizey
const auto algo_update = [uv_data, imsizex, imsizey, padmm_weak,
iter](const Vector<t_complex> &x) -> bool {
auto padmm = padmm_weak.lock();
PURIFY_MEDIUM_LOG("Step size γ {}", padmm->gamma());
PURIFY_MEDIUM_LOG("Step size γ {}", padmm->regulariser_strength());
*iter = *iter + 1;
Vector<t_complex> const alpha = padmm->Psi().adjoint() * x;
// updating parameter
const t_real new_gamma = alpha.real().cwiseAbs().maxCoeff() * 1e-3;
PURIFY_MEDIUM_LOG("Step size γ update {}", new_gamma);
padmm->gamma(((std::abs(padmm->gamma() - new_gamma) > 0.2) and *iter < 200) ? new_gamma
: padmm->gamma());
const t_real new_regulariser_strength = alpha.real().cwiseAbs().maxCoeff() * 1e-3;
PURIFY_MEDIUM_LOG("Step size γ update {}", new_regulariser_strength);
padmm->regulariser_strength(
((std::abs(padmm->regulariser_strength() - new_regulariser_strength) > 0.2) and *iter < 200)
? new_regulariser_strength
: padmm->regulariser_strength());

Vector<t_complex> const residual = padmm->Phi().adjoint() * (uv_data.vis - padmm->Phi() * x);

Expand Down
8 changes: 4 additions & 4 deletions cpp/example/padmm_reweighted_simulation.cc
Original file line number Diff line number Diff line change
Expand Up @@ -83,15 +83,15 @@ int main(int nargs, char const **args) {
Vector<t_complex> initial_estimate = Vector<t_complex>::Zero(dimage.size());

auto const epsilon = utilities::calculate_l2_radius(uv_data.vis.size(), sigma);
auto const purify_gamma =
auto const purify_regulariser_strength =
(Psi.adjoint() * (measurements_transform.adjoint() * uv_data.vis).eval()).real().maxCoeff() *
1e-3;

PURIFY_HIGH_LOG("Starting sopt!");
PURIFY_MEDIUM_LOG("Epsilon {}", epsilon);
PURIFY_MEDIUM_LOG("Gamma {}", purify_gamma);
PURIFY_MEDIUM_LOG("Regulariser_Strength {}", purify_regulariser_strength);
auto const padmm = sopt::algorithm::ImagingProximalADMM<t_complex>(uv_data.vis)
.gamma(purify_gamma)
.regulariser_strength(purify_regulariser_strength)
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon * 1.001)
.tight_frame(false)
Expand All @@ -102,7 +102,7 @@ int main(int nargs, char const **args) {
.l1_proximal_real_constraint(true)
.residual_convergence(epsilon * 1.001)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(measurements_transform);
// Timing reconstruction
Expand Down
10 changes: 5 additions & 5 deletions cpp/example/padmm_simulation.cc
Original file line number Diff line number Diff line change
Expand Up @@ -94,7 +94,7 @@ int main(int nargs, char const **args) {
Vector<t_complex> initial_estimate = Vector<t_complex>::Zero(dimage.size());

auto const epsilon = utilities::calculate_l2_radius(uv_data.vis.size(), sigma);
auto const purify_gamma =
auto const purify_regulariser_strength =
(Psi.adjoint() * (measurements_transform.adjoint() * uv_data.vis).eval()).real().maxCoeff() *
1e-3;
t_int iters = 0;
Expand All @@ -104,9 +104,9 @@ int main(int nargs, char const **args) {
};
PURIFY_HIGH_LOG("Starting sopt!");
PURIFY_MEDIUM_LOG("Epsilon {}", epsilon);
PURIFY_MEDIUM_LOG("Gamma {}", purify_gamma);
PURIFY_MEDIUM_LOG("Regulariser_Strength {}", purify_regulariser_strength);
auto const padmm = sopt::algorithm::ImagingProximalADMM<t_complex>(uv_data.vis)
.gamma(purify_gamma)
.regulariser_strength(purify_regulariser_strength)
.relative_variation(1e-3)
.l2ball_proximal_epsilon(epsilon * 1.001)
.tight_frame(false)
Expand All @@ -117,12 +117,12 @@ int main(int nargs, char const **args) {
.l1_proximal_real_constraint(true)
.residual_convergence(epsilon * 1.001)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.itermax(100)
.is_converged(convergence_function)
.Phi(measurements_transform);
// Timing reconstructionu
// Timing reconstruction

std::clock_t c_start = std::clock();
auto const diagnostic = padmm();
Expand Down
10 changes: 5 additions & 5 deletions cpp/example/sara_padmm_random_coverage.cc
Original file line number Diff line number Diff line change
Expand Up @@ -74,19 +74,19 @@ int main(int, char **) {
pfitsio::write2d(Image<t_real>::Map(dimage.data(), M31.rows(), M31.cols()), dirty_image_fits);

auto const epsilon = utilities::calculate_l2_radius(uv_data.vis.size(), sigma);
auto const purify_gamma =
auto const purify_regulariser_strength =
(Psi.adjoint() * (measurements_transform.adjoint() * uv_data.vis).eval())
.cwiseAbs()
.maxCoeff() *
beta;

// auto purify_gamma = 3 * utilities::median((Psi.adjoint() * (measurements_transform.adjoint() *
// (uv_data.vis - y0))).real().cwiseAbs())/0.6745;
// auto purify_regulariser_strength = 3 * utilities::median((Psi.adjoint() *
// (measurements_transform.adjoint() * (uv_data.vis - y0))).real().cwiseAbs())/0.6745;

SOPT_INFO("Using epsilon of {} \n", epsilon);
auto const padmm = sopt::algorithm::ImagingProximalADMM<t_complex>(uv_data.vis)
.itermax(1000)
.gamma(purify_gamma)
.regulariser_strength(purify_regulariser_strength)
.relative_variation(1e-6)
.l2ball_proximal_epsilon(epsilon)
.tight_frame(false)
Expand All @@ -97,7 +97,7 @@ int main(int, char **) {
.l1_proximal_real_constraint(true)
.residual_convergence(epsilon * 1.001)
.lagrange_update_scale(0.9)
.nu(1e0)
.sq_op_norm(1e0)
.Psi(Psi)
.Phi(measurements_transform);

Expand Down