Last active
October 8, 2017 05:27
-
-
Save haotian-liu/0f3f2bdab700285de13d3b9cb4167146 to your computer and use it in GitHub Desktop.
Optimized algorithm.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include <iostream> | |
#include "gmp.h" | |
using namespace std; | |
uint64_t SimpCalc(int Y, int X) { | |
auto f = new uint64_t[X + 1]; | |
uint64_t sum; | |
memset(f, 0, sizeof(uint64_t) * (X + 1)); | |
f[1] = 2; | |
for (int i = 2; i <= Y; i++) { | |
sum = 0; | |
for (int j = X; j >= 2; j--) { | |
f[j] = f[j - 1]; | |
sum += f[j]; | |
} | |
f[1] = sum; | |
} | |
sum = 0; | |
for (int i = 1; i < X; i++) { | |
sum += f[i]; | |
} | |
return sum; | |
} | |
void GMPCalc(int Y, int X, mpz_t sum) { | |
auto f = new mpz_t[X + 1]; | |
for (int i = 0; i <= X; i++) { mpz_init(f[i]); } | |
mpz_set_ui(f[1], 2); | |
for (int i = 2; i <= Y; i++) { | |
mpz_set_ui(sum, 0); | |
for (int j = X; j >= 2; j--) { | |
mpz_set(f[j], f[j - 1]); | |
mpz_add(sum, sum, f[j]); | |
} | |
mpz_set(f[1], sum); | |
} | |
mpz_set_ui(sum, 0); | |
for (int i = 1; i < X; i++) { | |
mpz_add(sum, sum, f[i]); | |
} | |
for (int i = 0; i <= X; i++) { mpz_clear(f[i]); } | |
delete[] f; | |
} | |
double OptCalc(int X, int Y) { | |
auto f = new double[X + 1]; | |
int shifts = 0; | |
double sum; | |
double Threshold = (uint64_t)UINT32_MAX * UINT16_MAX; | |
memset(f, 0, sizeof(double) * (X + 1)); | |
f[1] = 2; | |
for (int i = 2; i <= Y; i++) { | |
sum = 0; | |
for (int j = X; j >= 2; j--) { | |
f[j] = f[j - 1]; | |
sum += f[j]; | |
} | |
f[1] = sum; | |
if (sum > Threshold) { | |
shifts += 32; | |
for (int j=1; j<X; j++) { | |
f[j] /= UINT32_MAX; | |
} | |
} | |
} | |
sum = 0; | |
for (int i = 1; i < X; i++) { | |
sum += f[i]; | |
} | |
uint64_t all = 1; | |
all <<= Y - shifts; | |
double ans = static_cast<double>(sum) / all; | |
if (ans >= 1) return 0; | |
return 1 - ans; | |
} | |
double Calc(int X, int Y) { | |
double p; | |
if (Y < 63) { | |
uint64_t all = 1, ans; | |
ans = SimpCalc(Y, X); | |
all <<= Y; | |
ans = all - ans; | |
p = static_cast<double>(ans) / all; | |
} else { | |
mpz_t ans, all, one; | |
mpf_t all_f, ans_f, p_mpf; | |
mpz_inits(ans, all, one, nullptr); | |
mpf_inits(all_f, ans_f, p_mpf, nullptr); | |
GMPCalc(Y, X, ans); | |
mpz_set_ui(one, 1); | |
mpz_mul_2exp(all, one, Y); | |
mpz_sub(ans, all, ans); | |
mpf_set_z(all_f, all); | |
mpf_set_z(ans_f, ans); | |
mpf_div(p_mpf, ans_f, all_f); | |
p = mpf_get_d(p_mpf); | |
mpz_clears(ans, all, one, nullptr); | |
mpf_clears(all_f, ans_f, p_mpf, nullptr); | |
} | |
return p; | |
} | |
int main() { | |
int Y, X; | |
cin >> X >> Y; | |
cout << OptCalc(X, Y) << endl; | |
// cout << Calc(X, Y) << endl; | |
return 0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment