# Project Euler 548
# Gozinta Chains — sum n <= 10^16 with g(n)=n.
extern {
function calloc(n: i64, size: i64) -> ptr<void>
function free(p: ptr<void>) -> void
}
const LIMIT: i64 = 10000000000000000
const MAX_EXP_SUM: i64 = 54
const MAX_N: i64 = 110
function gcd(a0: i64, b0: i64) -> i64 {
let mut a: i64 = a0
let mut b: i64 = b0
while b != 0 {
let t: i64 = a % b
a = b
b = t
}
return a
}
function init_comb(comb: ptr<i64>, binom: ptr<i64>) -> void {
let mut n: i64 = 0
while n <= MAX_N {
comb[n * MAX_N + 0] = 1
let mut k: i64 = 1
while k <= n {
if k < n {
comb[n * MAX_N + k] = comb[(n - 1) * MAX_N + (k - 1)] + comb[(n - 1) * MAX_N + k]
} else {
comb[n * MAX_N + k] = 1
}
k = k + 1
}
n = n + 1
}
n = 0
while n <= MAX_EXP_SUM {
let mut k: i64 = 0
while k <= n {
binom[n * (MAX_EXP_SUM + 1) + k] = comb[n * MAX_N + k]
k = k + 1
}
n = n + 1
}
}
function gozinta_from_sig(sig: ptr<i64>, sig_len: i64, comb: ptr<i64>, binom: ptr<i64>, lim: i64) -> i64 {
let mut s: i64 = 0
let mut i: i64 = 0
while i < sig_len { s = s + sig[i]; i = i + 1 }
if s == 0 { return 1 }
let P: ptr<i64> = calloc(s + 2, 8)
if P == null { return 0 }
let mut t: i64 = 1
while t <= s {
let mut prod: i64 = 1
let tt: i64 = t - 1
i = 0
while i < sig_len {
let a: i64 = sig[i]
prod = prod * comb[(a + t - 1) * MAX_N + tt]
i = i + 1
}
P[t] = prod
t = t + 1
}
let mut total: i64 = 0
let mut m: i64 = 1
while m <= s {
let mut A_m: i64 = 0
t = 1
while t <= m {
let term: i64 = binom[m * (MAX_EXP_SUM + 1) + t] * P[t]
if ((m - t) & 1) != 0 { A_m = A_m - term } else { A_m = A_m + term }
t = t + 1
}
total = total + A_m
if lim > 0 && total > lim { free(P); return lim + 1 }
m = m + 1
}
free(P)
return total
}
function sig_equals_g(sig: ptr<i64>, sig_len: i64, comb: ptr<i64>, binom: ptr<i64>) -> bool {
let g: i64 = gozinta_from_sig(sig, sig_len, comb, binom, LIMIT)
if g > LIMIT { return false }
let gs: ptr<i64> = calloc(16, 8)
if gs == null { return false }
let mut n: i64 = g
let mut glen: i64 = 0
if n == 1 {
free(gs)
if sig_len == 0 { return true }
return false
}
let mut p: i64 = 2
while p * p <= n {
if n % p == 0 {
let mut e: i64 = 0
while n % p == 0 {
n = n / p
e = e + 1
}
let mut j: i64 = glen
while j > 0 && gs[j - 1] < e {
gs[j] = gs[j - 1]
j = j - 1
}
gs[j] = e
glen = glen + 1
}
p = p + 1
}
if n > 1 {
let mut j: i64 = glen
while j > 0 && gs[j - 1] < 1 {
gs[j] = gs[j - 1]
j = j - 1
}
gs[j] = 1
glen = glen + 1
}
if glen != sig_len { free(gs); return false }
let mut i: i64 = 0
while i < sig_len {
if sig[i] != gs[i] { free(gs); return false }
i = i + 1
}
free(gs)
return true
}
function enumerate_sigs(first_primes: ptr<i64>, sig: ptr<i64>, sig_len: i64,
pos: i64, prev_e: i64, sum_used: i64, prod: i64,
comb: ptr<i64>, binom: ptr<i64>, total: ptr<i64>) -> void {
if sig_len > 0 {
if sig_equals_g(sig, sig_len, comb, binom) {
let g: i64 = gozinta_from_sig(sig, sig_len, comb, binom, LIMIT)
if g <= LIMIT { total[0] = total[0] + g }
}
}
if sum_used == MAX_EXP_SUM || pos >= 60 { return }
let p: i64 = first_primes[pos]
let max_e: i64 = prev_e
let rem: i64 = MAX_EXP_SUM - sum_used
let mut me: i64 = max_e
if rem < me { me = rem }
let mut power: i64 = p
let mut e: i64 = 1
while e <= me {
let new_prod: i64 = prod * power
if new_prod > LIMIT { break }
sig[sig_len] = e
enumerate_sigs(first_primes, sig, sig_len + 1, pos + 1, e, sum_used + e, new_prod, comb, binom, total)
if power > LIMIT / p { break }
power = power * p
e = e + 1
}
}
function main() -> i32 {
let comb: ptr<i64> = calloc((MAX_N + 1) * (MAX_N + 1), 8)
let binom: ptr<i64> = calloc((MAX_EXP_SUM + 1) * (MAX_EXP_SUM + 2), 8)
let first_primes: ptr<i64> = calloc(60, 8)
let sig: ptr<i64> = calloc(16, 8)
let total: ptr<i64> = calloc(1, 8)
if comb == null || binom == null || first_primes == null || sig == null || total == null {
return 1
}
init_comb(comb, binom)
total[0] = 1
let mut np: i64 = 0
let mut p: i64 = 2
while np < 60 {
if p == 2 || p == 3 || (p > 3 && (p % 2) != 0 && (p % 3) != 0) {
let mut ok: i32 = 1
if p > 3 {
let mut d: i64 = 5
while d * d <= p {
if p % d == 0 || p % (d + 2) == 0 { ok = 0; break }
d = d + 6
}
}
if ok == 1 {
first_primes[np] = p
np = np + 1
}
}
p = p + 1
}
enumerate_sigs(first_primes, sig, 0, 0, MAX_EXP_SUM, 0, 1, comb, binom, total)
printf("%lld\n", total[0])
free(comb); free(binom); free(first_primes); free(sig); free(total)
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 gcd_i64_i64(int64_t a0, int64_t b0);
void init_comb_ptr_i64_ptr_i64(int64_t* comb, int64_t* binom);
int64_t gozinta_from_sig_ptr_i64_i64_ptr_i64_ptr_i64_i64(int64_t* sig, int64_t sig_len, int64_t* comb, int64_t* binom, int64_t lim);
bool sig_equals_g_ptr_i64_i64_ptr_i64_ptr_i64(int64_t* sig, int64_t sig_len, int64_t* comb, int64_t* binom);
void enumerate_sigs_ptr_i64_ptr_i64_i64_i64_i64_i64_i64_ptr_i64_ptr_i64_ptr_i64(int64_t* first_primes, int64_t* sig, int64_t sig_len, int64_t pos, int64_t prev_e, int64_t sum_used, int64_t prod, int64_t* comb, int64_t* binom, int64_t* total);
int32_t main(void);
static const int64_t LIMIT = 10000000000000000;
static const int64_t MAX_EXP_SUM = 54;
static const int64_t MAX_N = 110;
int64_t gcd_i64_i64(int64_t a0, int64_t b0) {
int64_t a = a0;
int64_t b = b0;
while (b != 0) {
int64_t t = FLOW_CHECKED_MOD((a), (b));
a = b;
b = t;
}
return a;
}
void init_comb_ptr_i64_ptr_i64(int64_t* comb, int64_t* binom) {
int64_t n = 0;
while (n <= MAX_N) {
comb[((n * MAX_N) + 0)] = 1;
int64_t k = 1;
while (k <= n) {
if (k < n) {
comb[((n * MAX_N) + k)] = (comb[(((n - 1) * MAX_N) + (k - 1))] + comb[(((n - 1) * MAX_N) + k)]);
} else {
comb[((n * MAX_N) + k)] = 1;
}
k = (k + 1);
}
n = (n + 1);
}
n = 0;
while (n <= MAX_EXP_SUM) {
int64_t k = 0;
while (k <= n) {
binom[((n * (MAX_EXP_SUM + 1)) + k)] = comb[((n * MAX_N) + k)];
k = (k + 1);
}
n = (n + 1);
}
}
int64_t gozinta_from_sig_ptr_i64_i64_ptr_i64_ptr_i64_i64(int64_t* sig, int64_t sig_len, int64_t* comb, int64_t* binom, int64_t lim) {
int64_t s = 0;
int64_t i = 0;
while (i < sig_len) {
s = (s + sig[i]);
i = (i + 1);
}
if (s == 0) {
return 1;
}
int64_t* P = (int64_t*)(calloc((s + 2), 8));
if (P == NULL) {
return 0;
}
int64_t t = 1;
while (t <= s) {
int64_t prod = 1;
int64_t tt = (t - 1);
i = 0;
while (i < sig_len) {
int64_t a = sig[i];
prod = (prod * comb[((((a + t) - 1) * MAX_N) + tt)]);
i = (i + 1);
}
P[t] = prod;
t = (t + 1);
}
int64_t total = 0;
int64_t m = 1;
while (m <= s) {
int64_t A_m = 0;
t = 1;
while (t <= m) {
int64_t term = (binom[((m * (MAX_EXP_SUM + 1)) + t)] * P[t]);
if (((m - t) & 1) != 0) {
A_m = (A_m - term);
} else {
A_m = (A_m + term);
}
t = (t + 1);
}
total = (total + A_m);
if ((lim > 0 && total > lim)) {
free(P);
return (lim + 1);
}
m = (m + 1);
}
free(P);
return total;
}
bool sig_equals_g_ptr_i64_i64_ptr_i64_ptr_i64(int64_t* sig, int64_t sig_len, int64_t* comb, int64_t* binom) {
int64_t g = gozinta_from_sig_ptr_i64_i64_ptr_i64_ptr_i64_i64(sig, sig_len, comb, binom, LIMIT);
if (g > LIMIT) {
return 0;
}
int64_t* gs = (int64_t*)(calloc(16, 8));
if (gs == NULL) {
return 0;
}
int64_t n = g;
int64_t glen = 0;
if (n == 1) {
free(gs);
if (sig_len == 0) {
return 1;
}
return 0;
}
int64_t p = 2;
while ((p * p) <= n) {
if (FLOW_CHECKED_MOD((n), (p)) == 0) {
int64_t e = 0;
while (FLOW_CHECKED_MOD((n), (p)) == 0) {
n = FLOW_CHECKED_DIV((n), (p));
e = (e + 1);
}
int64_t j = glen;
while ((j > 0 && gs[(j - 1)] < e)) {
gs[j] = gs[(j - 1)];
j = (j - 1);
}
gs[j] = e;
glen = (glen + 1);
}
p = (p + 1);
}
if (n > 1) {
int64_t j = glen;
while ((j > 0 && gs[(j - 1)] < 1)) {
gs[j] = gs[(j - 1)];
j = (j - 1);
}
gs[j] = 1;
glen = (glen + 1);
}
if (glen != sig_len) {
free(gs);
return 0;
}
int64_t i = 0;
while (i < sig_len) {
if (sig[i] != gs[i]) {
free(gs);
return 0;
}
i = (i + 1);
}
free(gs);
return 1;
}
void enumerate_sigs_ptr_i64_ptr_i64_i64_i64_i64_i64_i64_ptr_i64_ptr_i64_ptr_i64(int64_t* first_primes, int64_t* sig, int64_t sig_len, int64_t pos, int64_t prev_e, int64_t sum_used, int64_t prod, int64_t* comb, int64_t* binom, int64_t* total) {
if (sig_len > 0) {
if (sig_equals_g_ptr_i64_i64_ptr_i64_ptr_i64(sig, sig_len, comb, binom)) {
int64_t g = gozinta_from_sig_ptr_i64_i64_ptr_i64_ptr_i64_i64(sig, sig_len, comb, binom, LIMIT);
if (g <= LIMIT) {
total[0] = (total[0] + g);
}
}
}
if ((sum_used == MAX_EXP_SUM || pos >= 60)) {
return;
}
int64_t p = first_primes[pos];
int64_t max_e = prev_e;
int64_t rem = (MAX_EXP_SUM - sum_used);
int64_t me = max_e;
if (rem < me) {
me = rem;
}
int64_t power = p;
int64_t e = 1;
while (e <= me) {
int64_t new_prod = (prod * power);
if (new_prod > LIMIT) {
break;
}
sig[sig_len] = e;
enumerate_sigs_ptr_i64_ptr_i64_i64_i64_i64_i64_i64_ptr_i64_ptr_i64_ptr_i64(first_primes, sig, (sig_len + 1), (pos + 1), e, (sum_used + e), new_prod, comb, binom, total);
if (power > FLOW_CHECKED_DIV((LIMIT), (p))) {
break;
}
power = (power * p);
e = (e + 1);
}
}
int32_t main(void) {
int64_t* comb = (int64_t*)(calloc(((MAX_N + 1) * (MAX_N + 1)), 8));
int64_t* binom = (int64_t*)(calloc(((MAX_EXP_SUM + 1) * (MAX_EXP_SUM + 2)), 8));
int64_t* first_primes = (int64_t*)(calloc(60, 8));
int64_t* sig = (int64_t*)(calloc(16, 8));
int64_t* total = (int64_t*)(calloc(1, 8));
if (((((comb == NULL || binom == NULL) || first_primes == NULL) || sig == NULL) || total == NULL)) {
return 1;
}
init_comb_ptr_i64_ptr_i64(comb, binom);
total[0] = 1;
int64_t np = 0;
int64_t p = 2;
while (np < 60) {
if (((p == 2 || p == 3) || ((p > 3 && FLOW_CHECKED_MOD((p), (2)) != 0) && FLOW_CHECKED_MOD((p), (3)) != 0))) {
int32_t ok = 1;
if (p > 3) {
int64_t d = 5;
while ((d * d) <= p) {
if ((FLOW_CHECKED_MOD((p), (d)) == 0 || FLOW_CHECKED_MOD((p), ((d + 2))) == 0)) {
ok = 0;
break;
}
d = (d + 6);
}
}
if (ok == 1) {
first_primes[np] = p;
np = (np + 1);
}
}
p = (p + 1);
}
enumerate_sigs_ptr_i64_ptr_i64_i64_i64_i64_i64_i64_ptr_i64_ptr_i64_ptr_i64(first_primes, sig, 0, 0, MAX_EXP_SUM, 0, 1, comb, binom, total);
printf("%lld\n", total[0]);
free(comb);
free(binom);
free(first_primes);
free(sig);
free(total);
return 0;
}