summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--langevin.cpp19
1 files changed, 7 insertions, 12 deletions
diff --git a/langevin.cpp b/langevin.cpp
index e4b6c8d..9da4ef6 100644
--- a/langevin.cpp
+++ b/langevin.cpp
@@ -12,16 +12,14 @@
using Rng = randutils::random_generator<pcg32>;
-Vector initializeVector(unsigned N, double a, Rng& r) {
+template <class Distribution, class Generator>
+Vector randomVector(unsigned N, Distribution d, Generator& r) {
Vector z(N);
- complex_normal_distribution<> dist(0, a, 0);
for (unsigned i = 0; i < N; i++) {
- z(i) = dist(r.engine());
+ z(i) = d(r);
}
- z *= sqrt(N) / sqrt(z.dot(z)); // Normalize.
-
return z;
}
@@ -38,10 +36,7 @@ Vector langevin(const Tensor& J, const Vector& z0, double T, double γ0,
while (!quit(W, steps)) {
double γ = pow(r.variate<double, std::normal_distribution>(0, γ0), 2);
- Vector η(z.size());
- for (unsigned i = 0; i < z.size(); i++) {
- η(i) = d(r.engine());
- }
+ Vector η = randomVector(z.size(), d, r.engine());
Vector zNext = z - γ * dW + η;
zNext *= sqrt(zNext.size()) / sqrt(zNext.dot(zNext));
@@ -105,10 +100,10 @@ int main(int argc, char* argv[]) {
double σ = sqrt(factorial(p) / (2.0 * pow(N, p - 1)));
Rng r;
- complex_normal_distribution<> d(0, σ, κ);
- Tensor J = generateCouplings<Scalar, PSPIN_P>(N, d, r.engine());
- Vector z = initializeVector(N, 100, r);
+ Tensor J = generateCouplings<Scalar, PSPIN_P>(N, complex_normal_distribution<>(0, σ, κ), r.engine());
+ Vector z = randomVector(N, complex_normal_distribution<>(0, 1, 0), r.engine());
+ z *= sqrt(N) / sqrt(z.dot(z)); // Normalize.
std::function<bool(double, unsigned)> findSaddle = [δ](double W, unsigned) {
std::cout << W << std::endl;