Files
@ 82182a2f068a
Branch filter:
Location: AENC/resampling_chain/montecarlo/resampler_segment.cpp - annotation
82182a2f068a
2.2 KiB
text/x-c++src
Initial commit with Monte Carlo sampler
82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a 82182a2f068a | #include <iostream>
#include <random>
#include "resampler.hpp"
std::mt19937 mt(std::random_device{}());
class SegmentResampler : public ResamplerBase<SegmentResampler, std::mt19937> {
public:
SegmentResampler(float p_, size_t n) : ResamplerBase(p_, n, mt) {
setVertex(0, BAD);
}
// No periodic boundary
void resampleNeighbors(int i) {
resampleVertex(i);
if (i != 0)
resampleVertex(i - 1);
if (i + 1 != (int)getN())
resampleVertex(i + 1);
}
};
void doRun(float p, size_t n) {
std::cerr << "Running p = " << p << " and n = " << n << std::endl;
constexpr int samples = 500;
int terminateCount = 0;
int hitnCount = 0;
for (int j = 0; j < samples; ++j) {
SegmentResampler sampler(p, n);
for (int t = 0; t != 500000; t++) {
sampler.doMove();
if (sampler.numBads() == 0) {
terminateCount++;
break;
}
if (sampler.getState()[n - 1] == BAD) {
hitnCount++;
break;
}
}
}
if (terminateCount + hitnCount != samples) {
std::cerr << (samples - terminateCount - hitnCount) << '/' << samples
<< " runs with timelimit!\n";
}
std::cout << '{' << p << ',' << n << ',';
std::cout << hitnCount << '/' << hitnCount + terminateCount << '}';
std::cout << std::flush;
}
int main() {
// Measure probability of hitting vertex N starting from single-bad
// Interesting range is p in [0.60,0.75]
std::cout << "(* Probability of hitting site n when starting from single "
"BAD at 0 on segment [0,n].\n";
std::cout << "Data-format: {p, n, probability} *)\n";
std::cout << "{\n";
for (size_t n = 600; ; n += 200) {
doRun(0.10f, n);
std::cout << ',';
for (int i = 0; i <= 100; ++i) {
float p = 0.60f + 0.001f * i;
doRun(p, n);
std::cout << ',';
}
doRun(0.75f, n);
std::cout << ',';
doRun(0.80f, n);
std::cout << ',';
doRun(0.90f, n);
if (n == 1000)
break;
else
std::cout << ',';
}
std::cout << '}' << std::endl;
return 0;
}
|