# Project Euler 643: 2-Friendly
# f(n) = sum_{t>=1} (S(floor(n/2^t)) - 1) mod 1e9+7
# where S(m) = sum_{k=1..m} phi(k).
# Uses Du Jiao sieve with pre-enumerated floor(N/i) values.
import euler.nt { isqrt }
extern {
function calloc(n: i64, size: i64) -> ptr<void>
function free(p: ptr<void>) -> void
}
const MOD: i64 = 1000000007
const LIMIT: i64 = 10000000
const INV2: i64 = 500000004
# Binary search for key in descending-sorted array, return index or -1
function bsearch_desc(keys: ptr<i64>, n: i64, key: i64) -> i64 {
let mut lo: i64 = 0
let mut hi: i64 = n - 1
while lo <= hi {
let mid: i64 = (lo + hi) / 2
if keys[mid] == key { return mid }
if keys[mid] > key { lo = mid + 1 }
else { hi = mid - 1 }
}
return -1
}
function main() -> i32 {
let N: i64 = 100000000000 # 10^11
# Step 1: Precompute phi prefix sums up to LIMIT
let phi: ptr<i64> = calloc(LIMIT + 1, 8)
if phi == null { return 1 }
phi[1] = 1
let is_comp: ptr<i8> = calloc(LIMIT + 1, 1)
let primes_arr: ptr<i64> = calloc(LIMIT / 10 + 100, 8)
let mut npc: i64 = 0
for i in 2..(LIMIT + 1) {
if is_comp[i] == 0 {
primes_arr[npc] = i
npc = npc + 1
phi[i] = i - 1
}
for j in 0..npc {
let p: i64 = primes_arr[j]
let ip: i64 = i * p
if ip > LIMIT { break }
is_comp[ip] = 1
if i % p == 0 {
phi[ip] = phi[i] * p
break
} else {
phi[ip] = phi[i] * (p - 1)
}
}
}
# Prefix sum of phi mod MOD
let pref: ptr<i64> = calloc(LIMIT + 1, 8)
let mut s: i64 = 0
for i in 1..(LIMIT + 1) {
s = (s + phi[i]) % MOD
pref[i] = s
}
free(is_comp)
free(primes_arr)
free(phi)
# Step 2: Enumerate all distinct floor(N/i) values for i >= 1
# Large values: N/i for i=1..root_n (already decreasing)
# Small values: 1..small_max (increasing, reverse to get decreasing)
# Merge them (large are all >= root_n, small are all <= root_n, possible overlap at root_n)
let root_n: i64 = isqrt(N)
let max_vals: i64 = 2 * root_n + 10
let keys: ptr<i64> = calloc(max_vals, 8)
let vals: ptr<i64> = calloc(max_vals, 8)
let mut nkeys: i64 = 0
# Large values: N/1, N/2, ..., N/root_n (decreasing)
for i in 1..(root_n + 1) {
keys[nkeys] = N / i
nkeys = nkeys + 1
}
# Small values: small_max down to 1 (decreasing)
let small_max: i64 = N / (root_n + 1)
for v in 0..small_max {
keys[nkeys] = small_max - v
nkeys = nkeys + 1
}
# Remove duplicates (consecutive equal values)
let mut unique: i64 = 0
for i in 0..nkeys {
if i == 0 || keys[i] != keys[i - 1] {
keys[unique] = keys[i]
unique = unique + 1
}
}
nkeys = unique
# Step 3: Compute S(v) for each v in INCREASING order
# (S(v) depends on S(w) for w < v, so process smallest first)
for idx in 0..nkeys {
let rev_idx: i64 = nkeys - 1 - idx
let v: i64 = keys[rev_idx]
if v <= LIMIT {
vals[rev_idx] = pref[v]
continue
}
if v == 0 {
vals[rev_idx] = 0
continue
}
# S(v) = v*(v+1)/2 - sum_{m=2..v} S(floor(v/m))
# Work mod MOD throughout
let mut result: i64 = 0
let v_mod: i64 = v % MOD
let vp1_mod: i64 = (v + 1) % MOD
# v*(v+1)/2 mod MOD = v_mod * vp1_mod / 2 mod MOD
# Since MOD is odd, /2 = * inverse(2) = * (MOD+1)/2
let mut tri_mod: i128 = (v_mod as i128) * (vp1_mod as i128) % (MOD as i128)
tri_mod = tri_mod * (INV2 as i128) % (MOD as i128)
result = tri_mod as i64
# Group by floor(v/m) values
let mut m: i64 = 2
while m <= v {
let w: i64 = v / m
if w == 0 { break }
let m_end: i64 = v / w # largest m with floor(v/m) = w
# Look up S(w)
let sw: i64
if w <= LIMIT {
sw = pref[w]
} else {
let wi: i64 = bsearch_desc(keys, nkeys, w)
if wi < 0 { return 1 }
sw = vals[wi]
}
let count: i64 = m_end - m + 1
let count_mod: i64 = count % MOD
let term: i128 = (count_mod as i128) * (sw as i128) % (MOD as i128)
let mut r128: i128 = (result as i128) - term
r128 = r128 % (MOD as i128)
if r128 < 0 { r128 = r128 + (MOD as i128) }
result = r128 as i64
m = m_end + 1
}
vals[rev_idx] = result
}
# Step 4: Compute f(N) = sum_{t>=1} (S(floor(N/2^t)) - 1)
let mut ans: i64 = 0
let mut t: i64 = 1
while (N >> t) > 0 {
let m_val: i64 = N >> t
let sm: i64
if m_val <= LIMIT {
sm = pref[m_val]
} else {
let mi: i64 = bsearch_desc(keys, nkeys, m_val)
if mi < 0 { return 1 }
sm = vals[mi]
}
ans = (ans + sm - 1 + MOD) % MOD
t = t + 1
}
printf("%lld\n", ans)
free(vals)
free(keys)
free(pref)
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);
int64_t lcm_i64_i64(int64_t a, int64_t b);
int64_t isqrt_i64(int64_t n);
int64_t mulmod_i64_i64_i64(int64_t a0, int64_t b0, int64_t mod);
int64_t mod_pow_i64_i64_i64(int64_t base, int64_t exp, int64_t mod);
bool is_prime_i64(int64_t n);
int64_t bsearch_desc_ptr_i64_i64_i64(int64_t* keys, int64_t n, int64_t key);
int32_t main(void);
static const int64_t MOD = 1000000007;
static const int64_t LIMIT = 10000000;
static const int64_t INV2 = 500000004;
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;
}
int64_t lcm_i64_i64(int64_t a, int64_t b) {
if ((a == 0 || b == 0)) {
return 0;
}
return (FLOW_CHECKED_DIV((a), (gcd_i64_i64(a, b))) * b);
}
int64_t isqrt_i64(int64_t n) {
if (n < 2) {
return n;
}
int64_t x = n;
int64_t y = FLOW_CHECKED_DIV(((x + 1)), (2));
while (y < x) {
x = y;
y = FLOW_CHECKED_DIV(((x + FLOW_CHECKED_DIV((n), (x)))), (2));
}
return x;
}
int64_t mulmod_i64_i64_i64(int64_t a0, int64_t b0, int64_t mod) {
int64_t a = FLOW_CHECKED_MOD((a0), (mod));
int64_t b = FLOW_CHECKED_MOD((b0), (mod));
int64_t result = 0;
while (b > 0) {
if (FLOW_CHECKED_MOD((b), (2)) == 1) {
result = FLOW_CHECKED_MOD(((result + a)), (mod));
}
a = FLOW_CHECKED_MOD(((a * 2)), (mod));
b = FLOW_CHECKED_DIV((b), (2));
}
return result;
}
int64_t mod_pow_i64_i64_i64(int64_t base, int64_t exp, int64_t mod) {
if (mod == 1) {
return 0;
}
int64_t result = 1;
int64_t b = FLOW_CHECKED_MOD((base), (mod));
int64_t e = exp;
while (e > 0) {
if (FLOW_CHECKED_MOD((e), (2)) == 1) {
result = mulmod_i64_i64_i64(result, b, mod);
}
b = mulmod_i64_i64_i64(b, b, mod);
e = FLOW_CHECKED_DIV((e), (2));
}
return result;
}
bool is_prime_i64(int64_t n) {
if (n < 2) {
return 0;
}
if (n < 4) {
return 1;
}
if ((FLOW_CHECKED_MOD((n), (2)) == 0 || FLOW_CHECKED_MOD((n), (3)) == 0)) {
return 0;
}
int64_t i = 5;
while ((i * i) <= n) {
if ((FLOW_CHECKED_MOD((n), (i)) == 0 || FLOW_CHECKED_MOD((n), ((i + 2))) == 0)) {
return 0;
}
i = (i + 6);
}
return 1;
}
int64_t bsearch_desc_ptr_i64_i64_i64(int64_t* keys, int64_t n, int64_t key) {
int64_t lo = 0;
int64_t hi = (n - 1);
while (lo <= hi) {
int64_t mid = FLOW_CHECKED_DIV(((lo + hi)), (2));
if (keys[mid] == key) {
return mid;
}
if (keys[mid] > key) {
lo = (mid + 1);
} else {
hi = (mid - 1);
}
}
return (-1);
}
int32_t main(void) {
int64_t N = 100000000000;
int64_t* phi = (int64_t*)(calloc((LIMIT + 1), 8));
if (phi == NULL) {
return 1;
}
phi[1] = 1;
int8_t* is_comp = (int8_t*)(calloc((LIMIT + 1), 1));
int64_t* primes_arr = (int64_t*)(calloc((FLOW_CHECKED_DIV((LIMIT), (10)) + 100), 8));
int64_t npc = 0;
int32_t __flow_step_1 = 1;
for (int32_t i = 2; (2 <= (LIMIT + 1)) ? i < (LIMIT + 1) : i > (LIMIT + 1); i += (2 <= (LIMIT + 1)) ? 1 : -1) {
if (is_comp[i] == 0) {
primes_arr[npc] = i;
npc = (npc + 1);
phi[i] = (i - 1);
}
int32_t __flow_step_2 = 1;
for (int32_t j = 0; (0 <= npc) ? j < npc : j > npc; j += (0 <= npc) ? 1 : -1) {
int64_t p = primes_arr[j];
int64_t ip = (i * p);
if (ip > LIMIT) {
break;
}
is_comp[ip] = 1;
if (FLOW_CHECKED_MOD((i), (p)) == 0) {
phi[ip] = (phi[i] * p);
break;
} else {
phi[ip] = (phi[i] * (p - 1));
}
}
}
int64_t* pref = (int64_t*)(calloc((LIMIT + 1), 8));
int64_t s = 0;
int32_t __flow_step_3 = 1;
for (int32_t i = 1; (1 <= (LIMIT + 1)) ? i < (LIMIT + 1) : i > (LIMIT + 1); i += (1 <= (LIMIT + 1)) ? 1 : -1) {
s = FLOW_CHECKED_MOD(((s + phi[i])), (MOD));
pref[i] = s;
}
free(is_comp);
free(primes_arr);
free(phi);
int64_t root_n = isqrt_i64(N);
int64_t max_vals = ((2 * root_n) + 10);
int64_t* keys = (int64_t*)(calloc(max_vals, 8));
int64_t* vals = (int64_t*)(calloc(max_vals, 8));
int64_t nkeys = 0;
int32_t __flow_step_4 = 1;
for (int32_t i = 1; (1 <= (root_n + 1)) ? i < (root_n + 1) : i > (root_n + 1); i += (1 <= (root_n + 1)) ? 1 : -1) {
keys[nkeys] = FLOW_CHECKED_DIV((N), (i));
nkeys = (nkeys + 1);
}
int64_t small_max = FLOW_CHECKED_DIV((N), ((root_n + 1)));
int32_t __flow_step_5 = 1;
for (int32_t v = 0; (0 <= small_max) ? v < small_max : v > small_max; v += (0 <= small_max) ? 1 : -1) {
keys[nkeys] = (small_max - v);
nkeys = (nkeys + 1);
}
int64_t unique = 0;
int32_t __flow_step_6 = 1;
for (int32_t i = 0; (0 <= nkeys) ? i < nkeys : i > nkeys; i += (0 <= nkeys) ? 1 : -1) {
if ((i == 0 || keys[i] != keys[(i - 1)])) {
keys[unique] = keys[i];
unique = (unique + 1);
}
}
nkeys = unique;
int32_t __flow_step_7 = 1;
for (int32_t idx = 0; (0 <= nkeys) ? idx < nkeys : idx > nkeys; idx += (0 <= nkeys) ? 1 : -1) {
int64_t rev_idx = ((nkeys - 1) - idx);
int64_t v = keys[rev_idx];
if (v <= LIMIT) {
vals[rev_idx] = pref[v];
continue;
}
if (v == 0) {
vals[rev_idx] = 0;
continue;
}
int64_t result = 0;
int64_t v_mod = FLOW_CHECKED_MOD((v), (MOD));
int64_t vp1_mod = FLOW_CHECKED_MOD(((v + 1)), (MOD));
__int128 tri_mod = FLOW_CHECKED_MOD(((((__int128)(v_mod)) * ((__int128)(vp1_mod)))), (((__int128)(MOD))));
tri_mod = FLOW_CHECKED_MOD(((tri_mod * ((__int128)(INV2)))), (((__int128)(MOD))));
result = ((int64_t)(tri_mod));
int64_t m = 2;
while (m <= v) {
int64_t w = FLOW_CHECKED_DIV((v), (m));
if (w == 0) {
break;
}
int64_t m_end = FLOW_CHECKED_DIV((v), (w));
int64_t sw;
if (w <= LIMIT) {
sw = pref[w];
} else {
int64_t wi = bsearch_desc_ptr_i64_i64_i64(keys, nkeys, w);
if (wi < 0) {
return 1;
}
sw = vals[wi];
}
int64_t count = ((m_end - m) + 1);
int64_t count_mod = FLOW_CHECKED_MOD((count), (MOD));
__int128 term = FLOW_CHECKED_MOD(((((__int128)(count_mod)) * ((__int128)(sw)))), (((__int128)(MOD))));
__int128 r128 = (((__int128)(result)) - term);
r128 = FLOW_CHECKED_MOD((r128), (((__int128)(MOD))));
if (r128 < 0) {
r128 = (r128 + ((__int128)(MOD)));
}
result = ((int64_t)(r128));
m = (m_end + 1);
}
vals[rev_idx] = result;
}
int64_t ans = 0;
int64_t t = 1;
while (FLOW_CHECKED_SHR((N), (t)) > 0) {
int64_t m_val = FLOW_CHECKED_SHR((N), (t));
int64_t sm;
if (m_val <= LIMIT) {
sm = pref[m_val];
} else {
int64_t mi = bsearch_desc_ptr_i64_i64_i64(keys, nkeys, m_val);
if (mi < 0) {
return 1;
}
sm = vals[mi];
}
ans = FLOW_CHECKED_MOD(((((ans + sm) - 1) + MOD)), (MOD));
t = (t + 1);
}
printf("%lld\n", ans);
free(vals);
free(keys);
free(pref);
return 0;
}