Skip to content
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.

Sign up for GitHub

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

Gaussian HMC - integration within cooling gaussians and complexity improvements #301

Closed
wants to merge 5 commits into from
Closed

Gaussian HMC - integration within cooling gaussians and complexity improvements #301

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
cca44a0
Gaussian HMC in Cooling Gaussian
vgnecula May 27, 2024
07b98bc
Trigonometric Trick for imporved per-step complexity of Gaussian HMC
vgnecula May 27, 2024
0426762
Cooling Gaussians Integration Examples
vgnecula May 28, 2024
ede9bf0
Integration of Gaussian HMC Length Step as a passed parameter in Cool…
vgnecula May 30, 2024
4f76d52
Gaussian HMC -> length parameter integration in cooling gaussians - d…
vgnecula Jun 4, 2024
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
34 changes: 17 additions & 17 deletions include/random_walks/gaussian_hamiltonian_monte_carlo_exact_walk.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -224,22 +224,23 @@ private :
}
}

inline void update_position(Point &p, Point &v, NT const& T, NT const& omega)
{
NT C, Phi;
for (size_t i = 0; i < p.dimension(); i++)
{
C = std::sqrt(p[i] * p[i] + (v[i] * v[i]) / (omega * omega));
Phi = std::atan((-v[i]) / (p[i] * omega));
if (v[i] < 0.0 && Phi < 0.0) {
Phi += M_PI;
} else if (v[i] > 0.0 && Phi > 0.0) {
Phi -= M_PI;
}
p.set_coord(i, C * std::cos(omega * T + Phi));
v.set_coord(i, -C * omega * std::sin(omega * T + Phi));
}
inline void update_position(Point& p, Point& v, const NT& T, const NT& omega) {
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think I prefer the constant reference of type T to be written as T const& is there any reason for rewriting as const T&?


NT next_p, next_v;

NT sinVal = std::sin(omega * T);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this is the trick, thanks, great!

NT cosVal = std::cos(omega * T);

NT factor1 = sinVal / omega;
NT factor2 = -omega * sinVal;

for (size_t i = 0; i < p.dimension(); i++) {
next_p = cosVal * p[i] + v[i] * factor1;
next_v = factor2 * p[i] + cosVal * v[i];

p.set_coord(i, next_p);
v.set_coord(i, next_v);
}
}

inline double get_max_distance(std::vector<Point> &pointset, Point const& q, double &rad)
Expand Down Expand Up @@ -271,5 +272,4 @@ private :
};


#endif // RANDOM_WALKS_GAUSSIAN_HMC_WALK_HPP

#endif // RANDOM_WALKS_GAUSSIAN_HMC_WALK_HPP
10 changes: 4 additions & 6 deletions include/volume/volume_cooling_gaussians.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,9 +19,7 @@
#include <chrono>

#include "cartesian_geom/cartesian_kernel.h"
#include "random_walks/gaussian_helpers.hpp"
#include "random_walks/gaussian_ball_walk.hpp"
#include "random_walks/gaussian_cdhr_walk.hpp"
#include "random_walks/random_walks.hpp"
#include "sampling/random_point_generators.hpp"
#include "volume/math_helpers.hpp"

Expand Down Expand Up @@ -457,7 +455,7 @@ double volume_cooling_gaussians(Polytope& Pin,

template
<
typename WalkTypePolicy = GaussianCDHRWalk,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please don't change the default random walk in volume_cooling_gaussians<>.
You can create a new example in folder examples that it will call this function templated with GaussianHamiltonianMonteCarloExactWalk. You can write 3-4 examples that compute volumes of polytopes from 10 to 50 dimensions. You can use the volesti polytope generators.

typename WalkTypePolicy = GaussianHamiltonianMonteCarloExactWalk,
typename RandomNumberGenerator = BoostRandomNumberGenerator<boost::mt11213b, double>,
typename Polytope
>
Expand All @@ -472,7 +470,7 @@ double volume_cooling_gaussians(Polytope &Pin,

template
<
typename WalkTypePolicy = GaussianCDHRWalk,
typename WalkTypePolicy = GaussianHamiltonianMonteCarloExactWalk,
typename RandomNumberGenerator = BoostRandomNumberGenerator<boost::mt11213b, double>,
typename Polytope
>
Expand All @@ -487,4 +485,4 @@ double volume_cooling_gaussians(Polytope &Pin,
return volume_cooling_gaussians<WalkTypePolicy>(Pin, rng, error, walk_length);
}

#endif // VOLUME_COOLING_GAUSSIANS_HPP
#endif // VOLUME_COOLING_GAUSSIANS_HPP
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

please keep the empty ending line

Loading