summaryrefslogtreecommitdiff
path: root/hadamard.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'hadamard.cpp')
-rw-r--r--hadamard.cpp158
1 files changed, 31 insertions, 127 deletions
diff --git a/hadamard.cpp b/hadamard.cpp
index 1afb270..8eba184 100644
--- a/hadamard.cpp
+++ b/hadamard.cpp
@@ -3,6 +3,7 @@
#include <fstream>
#include <iostream>
+#include <chrono>
class MeasureEnergy : public Measurement {
public:
@@ -33,20 +34,17 @@ public:
class MeasureTransitionRates : public ParallelMeasurement {
public:
- std::vector<std::vector<unsigned>> nAccepted;
+ std::vector<unsigned> nAccepted;
unsigned total_steps;
- MeasureTransitionRates(unsigned n) : nAccepted(n - 1) {
+ MeasureTransitionRates(unsigned n) : nAccepted(n - 1, 0) {
total_steps = 0;
- for (unsigned i = 0; i < n - 1; i++) {
- nAccepted[i].resize(i + 1);
- }
}
- void after_step(bool accepted, unsigned i, unsigned j, double, double, const MCMC&,
+ void after_step(bool accepted, unsigned i, double, double, const MCMC&,
const MCMC&) override {
if (accepted)
- nAccepted[j - 1][i]++;
+ nAccepted[i]++;
}
void after_sweep(const std::vector<MCMC>&) override { total_steps++; }
@@ -54,14 +52,13 @@ public:
int main(int argc, char* argv[]) {
unsigned n_tuning = 1e2;
- std::list<double> β0s;
- std::list<double> β1s;
- std::list<unsigned> Ns;
+ double β₀ = 0.1;
+ double β₁ = 10;
+ unsigned N = 16;
unsigned k = 2;
double ε = 0.01;
unsigned M = 10;
- unsigned N = 1e4;
unsigned m = 1e4;
int opt;
@@ -72,16 +69,16 @@ int main(int argc, char* argv[]) {
k = atoi(optarg);
break;
case 'b':
- β0s.push_back(atof(optarg));
+ β₀ = atof(optarg);
break;
case 'c':
- β1s.push_back(atof(optarg));
+ β₁ = atof(optarg);
break;
case 'e':
ε = atof(optarg);
break;
case 'n':
- Ns.push_back((unsigned)atof(optarg));
+ m = (unsigned)atof(optarg);
break;
case 't':
n_tuning = (unsigned)atof(optarg);
@@ -99,33 +96,13 @@ int main(int argc, char* argv[]) {
unsigned n = pow(2, k);
- std::list<range> rs;
- unsigned num = 0;
-
- if ((β0s.size() != β1s.size()) || (β0s.size() != Ns.size())) {
- std::cout << "You need the same number of ranges!\n";
- exit(0);
- } else {
- auto it0 = β0s.begin();
- auto it1 = β1s.begin();
- auto itN = Ns.begin();
-
- while (it0 != β0s.end()) {
- rs.push_back({*it0, *it1, *itN});
- num += *itN;
- it0++;
- it1++;
- itN++;
- }
- }
-
- std::vector<Measurement*> As(num);
+ std::vector<Measurement*> As(N);
for (Measurement*& A : As) {
A = new MeasureEnergy();
}
- MeasureTransitionRates B(num);
+ MeasureTransitionRates B(N);
- PT p(rs, n, B, As);
+ PT p(β₀, β₁, N, n, B, As);
for (MCMC& sim : p.Ms) {
sim.M = walsh(k);
@@ -135,115 +112,42 @@ int main(int argc, char* argv[]) {
std::cout << "Beginning simulation of " << n << ".\n";
std::cout << "Beginning " << n_tuning << " tuning tempering updates of " << M
<< " sweeps each.\n";
- p.tune(n_tuning, M, ε);
- std::cout << "Finished tuning, beginning " << N << " measurement tempering updates of " << M
+ std::vector<double> f = p.tune(n_tuning, M, ε, 0.05);
+ std::cout << "Finished tuning, beginning " << m << " measurement tempering updates of " << M
<< " sweeps each.\n";
- p.run(N, M);
+ p.run(m, M);
std::cout << "Finished " << n << ".\n";
- std::string rs_string = "";
-
- for (range r : rs) {
- rs_string +=
- "_" + std::to_string(r.β0) + "_" + std::to_string(r.β1) + "_" + std::to_string(r.N);
- }
+ auto tag = std::chrono::high_resolution_clock::now();
- std::string filename = "probs_" + std::to_string(n) + rs_string + ".dat";
- std::ifstream file(filename);
+ std::string filename = "hmm_" + std::to_string(n) + "_" + std::to_string(β₀) + "_" + std::to_string(β₁) + "_" + std::to_string(N) + "_" + std::to_string(tag.time_since_epoch().count()) + ".dat";
- unsigned N_old = 0;
- std::vector<std::vector<unsigned long>> data_old(B.nAccepted.size());
+ std::ofstream file(filename);
- for (unsigned i = 0; i < B.nAccepted.size(); i++) {
- data_old[i].resize(B.nAccepted[i].size());
+ for (const MCMC& M : p.Ms) {
+ file << M.β << " ";
}
-
- if (file.is_open()) {
- file >> N_old;
-
- for (unsigned i = 0; i < B.nAccepted.size(); i++) {
- for (unsigned j = 0; j < B.nAccepted[i].size(); j++) {
- double num;
- file >> num;
- data_old[i][j] = num;
- }
- }
-
- file.close();
- }
-
- std::ofstream file_out(filename);
-
- file_out << N_old + B.total_steps << "\n";
+ file << std::endl;
for (unsigned i = 0; i < B.nAccepted.size(); i++) {
- for (unsigned j = 0; j < B.nAccepted[i].size(); j++) {
- file_out << std::fixed << data_old[i][j] + B.nAccepted[i][j] << " ";
- }
- file_out << "\n";
+ file << std::fixed << B.nAccepted[i] / (double)B.total_steps << " ";
}
-
- file_out.close();
-
- std::string efilename = "energies_" + std::to_string(n) + rs_string + ".dat";
- std::ifstream efile(efilename);
-
- unsigned Ne_old = 0;
- std::vector<double> edata_old(p.As.size());
-
- if (efile.is_open()) {
- efile >> Ne_old;
-
- for (unsigned i = 0; i < p.As.size(); i++) {
- double num;
- efile >> num;
- edata_old[i] = num;
- }
-
- efile.close();
- }
-
- std::ofstream efile_out(efilename);
-
- efile_out << Ne_old + ((MeasureEnergy*)As[0])->N << "\n";
-
- for (unsigned i = 0; i < As.size(); i++) {
- efile_out << std::fixed << edata_old[i] + ((MeasureEnergy*)As[i])->totalE << " ";
- }
-
- efile_out.close();
-
- std::string ρfilename = "rhos_" + std::to_string(n) + rs_string + ".dat";
- std::ifstream ρfile(ρfilename);
-
- std::vector<std::vector<unsigned>> ρdata_old(As.size());
+ file << std::endl;
for (unsigned i = 0; i < As.size(); i++) {
- ρdata_old[i].resize(((MeasureEnergy*)As[0])->ρ_dist.size());
- }
-
- if (ρfile.is_open()) {
- for (unsigned i = 0; i < As.size(); i++) {
- for (unsigned j = 0; j < ((MeasureEnergy*)As[0])->ρ_dist.size(); j++) {
- unsigned num;
- ρfile >> num;
- ρdata_old[i][j] = num;
- }
- }
-
- ρfile.close();
+ file << std::fixed << ((MeasureEnergy*)As[i])->totalE / ((MeasureEnergy*)As[i])->N << " ";
}
- std::ofstream ρfile_out(ρfilename);
+ file << std::endl;
for (unsigned i = 0; i < As.size(); i++) {
- for (unsigned j = 0; j < ((MeasureEnergy*)As[0])->ρ_dist.size(); j++) {
- ρfile_out << std::fixed << ρdata_old[i][j] + ((MeasureEnergy*)As[i])->ρ_dist[j] << " ";
+ for (unsigned j = 0; j < ((MeasureEnergy*)As[i])->ρ_dist.size(); j++) {
+ file << std::fixed << ((MeasureEnergy*)As[i])->ρ_dist[j] << " ";
}
- ρfile_out << "\n";
+ file << std::endl;
}
- ρfile_out.close();
+ file.close();
return 0;
}