# Project Euler 809: Rational Recurrence Relation
# f(22/7) mod 10^15 via Ackermann-Peter function and tetration.
# Pure Flow port of the native C solver. Uses i128 for CRT intermediates.
extern {
function calloc(n: i64, size: i64) -> ptr<void>
function free(p: ptr<void>) -> void
}
# Iterative extended GCD: returns g, sets x[0],y[0] with a*x + b*y = g.
function egcd(a0: i64, b0: i64, x: ptr<i64>, y: ptr<i64>) -> i64 {
let mut old_r: i64 = a0
let mut r: i64 = b0
let mut old_s: i64 = 1
let mut s: i64 = 0
let mut old_t: i64 = 0
let mut t: i64 = 1
while r != 0 {
let q: i64 = old_r / r
let tr: i64 = old_r - q * r
old_r = r
r = tr
let ts: i64 = old_s - q * s
old_s = s
s = ts
let tt: i64 = old_t - q * t
old_t = t
t = tt
}
x[0] = old_s
y[0] = old_t
return old_r
}
function inv_mod(a: i64, m: i64) -> i64 {
let xp: ptr<i64> = calloc(1, 8) as ptr<i64>
let yp: ptr<i64> = calloc(1, 8) as ptr<i64>
egcd(a % m, m, xp, yp)
let mut x: i64 = xp[0]
free(xp)
free(yp)
if x < 0 {
x = x + m
}
return x % m
}
# CRT: x ≡ r1 (mod m1), x ≡ r2 (mod m2), m1,m2 coprime. Result < m1*m2.
function crt(r1: i64, m1: i64, r2: i64, m2: i64) -> i64 {
let mut diff: i128 = ((r2 - r1) % m2) as i128
if diff < 0 {
diff = diff + (m2 as i128)
}
let k: i128 = diff * (inv_mod(m1 % m2, m2) as i128) % (m2 as i128)
let result: i128 = (r1 as i128) + (m1 as i128) * k
let prod: i128 = (m1 as i128) * (m2 as i128)
return (result % prod) as i64
}
# phi of n = 2^a * 5^b (rem assumed 1).
function phi_2_5(n: i64) -> i64 {
let mut x: i64 = n
let mut a: i64 = 0
while x % 2 == 0 {
a = a + 1
x = x / 2
}
let mut b: i64 = 0
while x % 5 == 0 {
b = b + 1
x = x / 5
}
let mut res: i64 = n
if a > 0 {
res = res / 2
}
if b > 0 {
res = (res / 5) * 4
}
return res
}
function totient_chain_len(n0: i64) -> i64 {
let mut n: i64 = n0
let mut steps: i64 = 0
while n != 1 {
n = phi_2_5(n)
steps = steps + 1
}
return steps
}
# min(2 tetrated height times, cap) for base 2.
function tetration_cap(height: i64, cap: i64) -> i64 {
if cap <= 0 {
return 0
}
let v: i64 = 2
if height <= 1 {
if v < cap {
return v
}
return cap
}
let mut i: i64 = 2
let mut vv: i64 = v
while i <= height {
if vv >= 60 {
return cap
}
vv = 1 << vv
if vv >= cap {
return cap
}
i = i + 1
}
return vv
}
# 2 tetrated height times (mod 2^a).
function tetration_mod_pow2(height: i64, a: i64) -> i64 {
if a <= 0 {
return 0
}
let modv: i64 = 1 << a
if height == 1 {
return 2 % modv
}
let exp: i64 = tetration_cap(height - 1, a)
if exp >= a {
return 0
}
return (1 << exp) % modv
}
function powmod_128(base: i64, exp0: i64, modv: i64) -> i64 {
let mut r: i128 = 1 % (modv as i128)
let mut b: i128 = base % (modv as i128)
let mut exp: i64 = exp0
while exp != 0 {
if exp % 2 == 1 {
r = r * b % (modv as i128)
}
b = b * b % (modv as i128)
exp = exp / 2
}
return r as i64
}
# 2 tetrated height times (mod modv) for modv of the form 2^a * 5^b.
function tetration_mod(height: i64, modv: i64) -> i64 {
if modv == 1 {
return 0
}
let mut x: i64 = modv
let mut a: i64 = 0
while x % 2 == 0 {
a = a + 1
x = x / 2
}
let mut b: i64 = 0
while x % 5 == 0 {
b = b + 1
x = x / 5
}
if b == 0 {
return tetration_mod_pow2(height, a)
}
if a == 0 {
if height == 1 {
return 2 % modv
}
let exp: i64 = tetration_mod(height - 1, phi_2_5(modv))
return powmod_128(2, exp, modv)
}
let m2: i64 = 1 << a
let mut m5: i64 = 1
let mut i: i64 = 0
while i < b {
m5 = m5 * 5
i = i + 1
}
let r2: i64 = tetration_mod_pow2(height, a)
let r5: i64 = tetration_mod(height, m5)
return crt(r2, m2, r5, m5) % (m2 * m5)
}
function stable_tetration_mod(modv: i64) -> i64 {
let height: i64 = totient_chain_len(modv) + 1
return tetration_mod(height, modv)
}
function main() -> i32 {
let mod2: i64 = 1 << 15
let mut mod5: i64 = 1
let mut i: i64 = 0
while i < 15 {
mod5 = mod5 * 5
i = i + 1
}
let mut r2: i64 = (-3) % mod2
if r2 < 0 {
r2 = r2 + mod2
}
let tower_mod5: i64 = stable_tetration_mod(mod5)
let mut r5: i64 = (tower_mod5 - 3) % mod5
if r5 < 0 {
r5 = r5 + mod5
}
let result: i64 = crt(r2, mod2, r5, mod5)
let mod10_15: i64 = mod2 * mod5
printf("%lld\n", result % mod10_15)
return 0
}
Generated C
#include <stdint.h>
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
/* Flow runtime helpers */
typedef struct flow_temp_node { struct flow_temp_node* next; } flow_temp_node;
static flow_temp_node* flow_temp_head = NULL;
static int flow_temp_atexit_set = 0;
__attribute__((unused)) static void flow_temp_free_all(void) {
while (flow_temp_head) {
flow_temp_node* n = flow_temp_head;
flow_temp_head = n->next;
free(n);
}
}
__attribute__((unused)) static void* flow_temp_alloc(size_t nbytes) {
flow_temp_node* node = (flow_temp_node*)malloc(sizeof(flow_temp_node) + nbytes);
if (!node) return NULL;
node->next = flow_temp_head;
flow_temp_head = node;
if (!flow_temp_atexit_set) {
flow_temp_atexit_set = 1;
atexit(flow_temp_free_all);
}
return (void*)(node + 1);
}
#ifndef FLOW_DIAG
#define FLOW_DIAG(msg) fprintf(stderr, "%s", (msg))
#endif
#ifndef FLOW_LOG
#define FLOW_LOG(fmt, ...) printf(fmt, __VA_ARGS__)
#endif
#ifndef FLOW_LOG_EMPTY
#define FLOW_LOG_EMPTY(fmt) printf(fmt)
#endif
static char* flow_strcat(const char* a, const char* b) {
size_t la = strlen(a ? a : ""), lb = strlen(b ? b : "");
char* r = (char*)flow_temp_alloc(la + lb + 1);
if (!r) return NULL;
if (la) memcpy(r, a, la);
if (lb) memcpy(r + la, b, lb);
r[la + lb] = '\0';
return r;
}
#define __flow_in_arr(arr, val) __extension__ ({ \
int _found = 0; \
size_t _n = sizeof(arr)/sizeof((arr)[0]); \
for (size_t _i = 0; _i < _n; _i++) { \
if ((arr)[_i] == (val)) { _found = 1; break; } \
} _found; })
/* Unified fault handler (MISRA #279) — override with -DFLOW_FAULT_HANDLER=fn */
#ifndef FLOW_FAULT_HANDLER
__attribute__((unused)) static inline void flow_fault_handler(const char* msg) {
fprintf(stderr, "flow: %s\n", msg ? msg : "fault");
abort();
#if defined(__GNUC__) || defined(__clang__)
__builtin_unreachable();
#endif
}
#else
#define flow_fault_handler FLOW_FAULT_HANDLER
#endif
#define flow_div_by_zero_handler() flow_fault_handler("division by zero")
#define flow_shift_ub_handler() flow_fault_handler("invalid shift (amount out of range or left-shift of negative)")
#ifndef FLOW_CHECKED_DIV
#define FLOW_CHECKED_DIV(L, R) (((R) != 0) ? ((L) / (R)) : (flow_div_by_zero_handler(), (L) * 0))
#endif
#ifndef FLOW_CHECKED_MOD
#define FLOW_CHECKED_MOD(L, R) (((R) != 0) ? ((L) % (R)) : (flow_div_by_zero_handler(), (L) * 0))
#endif
#ifndef FLOW_CHECKED_SHL
#define FLOW_CHECKED_SHL(L, R) ((((R) >= 0) && ((unsigned long long)(R) < (sizeof(L) * 8ull)) && ((L) >= 0)) ? ((L) << (R)) : (flow_shift_ub_handler(), (L) * 0))
#endif
#ifndef FLOW_CHECKED_SHR
#define FLOW_CHECKED_SHR(L, R) ((((R) >= 0) && ((unsigned long long)(R) < (sizeof(L) * 8ull))) ? ((L) >> (R)) : (flow_shift_ub_handler(), (L) * 0))
#endif
#include <math.h>
void* _ui_state = NULL;
static inline float i32_to_f32(int32_t v) { return (float)v; }
/* Host stub for @gpu kernels (device codegen replaces this). */
static inline int32_t gpu_thread_id(void) { return 0; }
int64_t egcd_i64_i64_ptr_i64_ptr_i64(int64_t a0, int64_t b0, int64_t* x, int64_t* y);
int64_t inv_mod_i64_i64(int64_t a, int64_t m);
int64_t crt_i64_i64_i64_i64(int64_t r1, int64_t m1, int64_t r2, int64_t m2);
int64_t phi_2_5_i64(int64_t n);
int64_t totient_chain_len_i64(int64_t n0);
int64_t tetration_cap_i64_i64(int64_t height, int64_t cap);
int64_t tetration_mod_pow2_i64_i64(int64_t height, int64_t a);
int64_t powmod_128_i64_i64_i64(int64_t base, int64_t exp0, int64_t modv);
int64_t tetration_mod_i64_i64(int64_t height, int64_t modv);
int64_t stable_tetration_mod_i64(int64_t modv);
int32_t main(void);
int64_t egcd_i64_i64_ptr_i64_ptr_i64(int64_t a0, int64_t b0, int64_t* x, int64_t* y) {
int64_t old_r = a0;
int64_t r = b0;
int64_t old_s = 1;
int64_t s = 0;
int64_t old_t = 0;
int64_t t = 1;
while (r != 0) {
int64_t q = FLOW_CHECKED_DIV((old_r), (r));
int64_t tr = (old_r - (q * r));
old_r = r;
r = tr;
int64_t ts = (old_s - (q * s));
old_s = s;
s = ts;
int64_t tt = (old_t - (q * t));
old_t = t;
t = tt;
}
x[0] = old_s;
y[0] = old_t;
return old_r;
}
int64_t inv_mod_i64_i64(int64_t a, int64_t m) {
int64_t* xp = (int64_t*)(((int64_t*)(calloc(1, 8))));
int64_t* yp = (int64_t*)(((int64_t*)(calloc(1, 8))));
egcd_i64_i64_ptr_i64_ptr_i64(FLOW_CHECKED_MOD((a), (m)), m, xp, yp);
int64_t x = xp[0];
free(xp);
free(yp);
if (x < 0) {
x = (x + m);
}
return FLOW_CHECKED_MOD((x), (m));
}
int64_t crt_i64_i64_i64_i64(int64_t r1, int64_t m1, int64_t r2, int64_t m2) {
__int128 diff = ((__int128)(FLOW_CHECKED_MOD(((r2 - r1)), (m2))));
if (diff < 0) {
diff = (diff + ((__int128)(m2)));
}
__int128 k = FLOW_CHECKED_MOD(((diff * ((__int128)(inv_mod_i64_i64(FLOW_CHECKED_MOD((m1), (m2)), m2))))), (((__int128)(m2))));
__int128 result = (((__int128)(r1)) + (((__int128)(m1)) * k));
__int128 prod = (((__int128)(m1)) * ((__int128)(m2)));
return ((int64_t)(FLOW_CHECKED_MOD((result), (prod))));
}
int64_t phi_2_5_i64(int64_t n) {
int64_t x = n;
int64_t a = 0;
while (FLOW_CHECKED_MOD((x), (2)) == 0) {
a = (a + 1);
x = FLOW_CHECKED_DIV((x), (2));
}
int64_t b = 0;
while (FLOW_CHECKED_MOD((x), (5)) == 0) {
b = (b + 1);
x = FLOW_CHECKED_DIV((x), (5));
}
int64_t res = n;
if (a > 0) {
res = FLOW_CHECKED_DIV((res), (2));
}
if (b > 0) {
res = (FLOW_CHECKED_DIV((res), (5)) * 4);
}
return res;
}
int64_t totient_chain_len_i64(int64_t n0) {
int64_t n = n0;
int64_t steps = 0;
while (n != 1) {
n = phi_2_5_i64(n);
steps = (steps + 1);
}
return steps;
}
int64_t tetration_cap_i64_i64(int64_t height, int64_t cap) {
if (cap <= 0) {
return 0;
}
int64_t v = 2;
if (height <= 1) {
if (v < cap) {
return v;
}
return cap;
}
int64_t i = 2;
int64_t vv = v;
while (i <= height) {
if (vv >= 60) {
return cap;
}
vv = FLOW_CHECKED_SHL((1), (vv));
if (vv >= cap) {
return cap;
}
i = (i + 1);
}
return vv;
}
int64_t tetration_mod_pow2_i64_i64(int64_t height, int64_t a) {
if (a <= 0) {
return 0;
}
int64_t modv = FLOW_CHECKED_SHL((1), (a));
if (height == 1) {
return FLOW_CHECKED_MOD((2), (modv));
}
int64_t exp = tetration_cap_i64_i64((height - 1), a);
if (exp >= a) {
return 0;
}
return FLOW_CHECKED_MOD((FLOW_CHECKED_SHL((1), (exp))), (modv));
}
int64_t powmod_128_i64_i64_i64(int64_t base, int64_t exp0, int64_t modv) {
__int128 r = FLOW_CHECKED_MOD((1), (((__int128)(modv))));
__int128 b = FLOW_CHECKED_MOD((base), (((__int128)(modv))));
int64_t exp = exp0;
while (exp != 0) {
if (FLOW_CHECKED_MOD((exp), (2)) == 1) {
r = FLOW_CHECKED_MOD(((r * b)), (((__int128)(modv))));
}
b = FLOW_CHECKED_MOD(((b * b)), (((__int128)(modv))));
exp = FLOW_CHECKED_DIV((exp), (2));
}
return ((int64_t)(r));
}
int64_t tetration_mod_i64_i64(int64_t height, int64_t modv) {
if (modv == 1) {
return 0;
}
int64_t x = modv;
int64_t a = 0;
while (FLOW_CHECKED_MOD((x), (2)) == 0) {
a = (a + 1);
x = FLOW_CHECKED_DIV((x), (2));
}
int64_t b = 0;
while (FLOW_CHECKED_MOD((x), (5)) == 0) {
b = (b + 1);
x = FLOW_CHECKED_DIV((x), (5));
}
if (b == 0) {
return tetration_mod_pow2_i64_i64(height, a);
}
if (a == 0) {
if (height == 1) {
return FLOW_CHECKED_MOD((2), (modv));
}
int64_t exp = tetration_mod_i64_i64((height - 1), phi_2_5_i64(modv));
return powmod_128_i64_i64_i64(2, exp, modv);
}
int64_t m2 = FLOW_CHECKED_SHL((1), (a));
int64_t m5 = 1;
int64_t i = 0;
while (i < b) {
m5 = (m5 * 5);
i = (i + 1);
}
int64_t r2 = tetration_mod_pow2_i64_i64(height, a);
int64_t r5 = tetration_mod_i64_i64(height, m5);
return FLOW_CHECKED_MOD((crt_i64_i64_i64_i64(r2, m2, r5, m5)), ((m2 * m5)));
}
int64_t stable_tetration_mod_i64(int64_t modv) {
int64_t height = (totient_chain_len_i64(modv) + 1);
return tetration_mod_i64_i64(height, modv);
}
int32_t main(void) {
int64_t mod2 = FLOW_CHECKED_SHL((1), (15));
int64_t mod5 = 1;
int64_t i = 0;
while (i < 15) {
mod5 = (mod5 * 5);
i = (i + 1);
}
int64_t r2 = FLOW_CHECKED_MOD(((-3)), (mod2));
if (r2 < 0) {
r2 = (r2 + mod2);
}
int64_t tower_mod5 = stable_tetration_mod_i64(mod5);
int64_t r5 = FLOW_CHECKED_MOD(((tower_mod5 - 3)), (mod5));
if (r5 < 0) {
r5 = (r5 + mod5);
}
int64_t result = crt_i64_i64_i64_i64(r2, mod2, r5, mod5);
int64_t mod10_15 = (mod2 * mod5);
printf("%lld\n", FLOW_CHECKED_MOD((result), (mod10_15)));
return 0;
}