cint

cint university / 201 / ced, the cint editor

201-09: functions

This is lesson 201-09 of CU-201.2, checked against cint 0.2.0.

You will pass numbers and arrays to functions, write a recursive function, read the call stack in a fault report and meet the limit on how deep calls can go.

the program

The harmonic number H(n) is 1 + 1/2 + 1/3 + ... + 1/n. This program adds it up as an exact fraction and keeps the fraction in lowest terms with Euclid's algorithm:

// Exact sums of fractions, reduced with Euclid's algorithm.

// Euclid: gcd(a, b) is gcd(b, a % b), and gcd(a, 0) is a.
I64 gcd(I64 a, I64 b) {
    if (b == 0) {
        return a;
    }
    return gcd(b, a % b);
}

// Add 1/k to the fraction f[0]/f[1] and reduce it to lowest terms.
void add_unit(inout I64[2] f, I64 k) {
    I64 num = f[0] * k + f[1];
    I64 den = f[1] * k;
    I64 g = gcd(num, den);
    f[0] = num / g;
    f[1] = den / g;
}

void main() {
    "gcd(1071, 462) = {gcd(1071, 462)}\n";
    I64[2] h;
    h[0] = 0;
    h[1] = 1;
    for k in 1..=50 {
        add_unit(h, k);
        if (k % 10 == 0) {
            "H({k}) = {h[0]}/{h[1]}\n";
        }
    }
}

Run this program in ced

run it

$ cint run harmonic.ci
gcd(1071, 462) = 21
H(10) = 7381/2520
H(20) = 55835135/15519504
H(30) = 9304682830147/2329089562800
H(40) = 2078178381193813/485721041551200
fault[E_OVERFLOW]: the exact result is outside the type's range
  --> harmonic.ci:13:20
   |
13 |     I64 num = f[0] * k + f[1];
   |                    ^ mul.checked.i64
   |
   = operation: mul.checked.i64
   = operand-count: 2
   = operand: I64 532145396070491417
   = operand: I64 44
   = exact:   23414397427101622348
   = limit:   I64 9223372036854775807
   = revision: 32960ad2bbac044b12ed41dcf78b45746985bc46b3eb99160afef5baf46c7225
   = source-map: none
   = address: none
   = stack-depth: 1
   = called from: harmonic.ci:26:9
   = fuel-consumed: 920 (not part of the fault record)
   = state: the program stopped; stdout holds what it wrote before the fault

The four fractions agree with fractions.Fraction in Python 3.13 on Ubuntu 24.04. Then the program stopped while adding 1/44.

reading the call stack

The fault happened inside add_unit, at line 13. stack-depth: 1 says one call was in progress below main, and called from: harmonic.ci:26:9 says where it was made: the line add_unit(h, k); in main. When calls are nested deeper, the report has one called from line for each call, starting from the call made in main.

The operands tell the rest. The first is the numerator of H(43), 532145396070491417, and the second is 44. Their product, 23414397427101622348, is more than twice the largest I64, so f[0] * k could not be computed.

parameters

A number passed to a function is copied. When gcd calls itself with gcd(b, a % b), its new a and b are new copies, and the caller's variables are unchanged.

An array is passed as a view of the caller's array, with a mode that says what the function may do with it:

  • in I64[2] f lets the function read the array but not write it; a write is error C2060.
  • inout I64[2] f lets it read and write the caller's array. add_unit writes the new fraction into h this way.

A number cannot be passed inout: void exchange(inout I64 a, inout I64 b) is error C2064. To give back one number, return it. To give back several, pass an array inout or return a struct, as lab 401-01 does.

A function's name cannot be the name of a built-in function. swap is a built-in function, so a function called swap is error C3001, like a variable called sum.

a smaller way to add

The program multiplies both denominators together and only then divides out the common part, so its numbers grow faster than the answer does. Adding over the least common denominator keeps them smaller. Replace add_unit with this version:

// Add 1/k to the fraction f[0]/f[1], over the least common denominator.
void add_unit(inout I64[2] f, I64 k) {
    I64 den = f[1] / gcd(f[1], k) * k;
    I64 num = f[0] * (den / f[1]) + den / k;
    I64 g = gcd(num, den);
    f[0] = num / g;
    f[1] = den / g;
}

Run this program in ced

It gets three numbers further and faults while adding 1/47:

fault[E_OVERFLOW]: the exact result is outside the type's range
  --> harmonic.ci:13:35
   |
13 |     I64 den = f[1] / gcd(f[1], k) * k;
   |                                   ^ mul.checked.i64
   |
   = operation: mul.checked.i64
   = operand-count: 2
   = operand: I64 1345655451257488800
   = operand: I64 47
   = exact:   63245806209101973600
   = limit:   I64 9223372036854775807
   ...

This time the exact line is the answer itself. 63245806209101973600 is the denominator of H(47), and it does not fit in an I64. No way of arranging the arithmetic can hold H(47) in this type. A wider type can hold it, and 301-04 carries the sum to H(100) in I256.

recursion and the depth limit

gcd calls itself until b is 0. Euclid's algorithm takes at most five steps for each decimal digit of the smaller number, so for any two I64 values it calls itself fewer than 100 times.

Every program has a limit on how deep calls can go: 256 calls in progress at once, counting main. A recursion that goes deeper faults with E_DEPTH:

I64 sum_to(I64 n) {
    if (n == 0) {
        return 0;
    }
    return n + sum_to(n - 1);
}
void main() {
    "{sum_to(100)}\n";
    "{sum_to(100_000)}\n";
}

Run this program in ced

sum_to(100) returns 5050, and sum_to(100_000) faults. sum_to(254) is the largest call that runs without a fault.

fault[E_DEPTH]: call depth limit exceeded
  --> deep.ci:5:16
   |
 5 |     return n + sum_to(n - 1);
   |                ^ call.enter
   |
   = operation: call.enter
   = operand-count: 0
   = exact:   none
   = limit:   I64 256
   ...
   = stack-depth: 255
   = called from: deep.ci:9:7
   = called from: deep.ci:5:16
   ...

main counts as one of the 256, so stack-depth stops at 255, and the report lists those 255 calls, with one called from line for each. A loop computes the same sum without adding to the call depth. Recursion suits work that shrinks by at least half every step or two, such as Euclid's algorithm or a binary search, because then the depth grows with the number of digits and not with the number itself.

elsewhere

The first version of the program in C, with int64_t, built with GCC 13.3 and -O2 -Wall -Wextra on Ubuntu 24.04 for x86-64, gave no warning and exited with status 0. Its last line was H(50) = -396549522021549939/211166340108011552: a negative sum of positive fractions. Signed overflow in C is undefined behavior, so the C standard places no requirement on what the program does next.

fractions.Fraction in Python 3.13 on Ubuntu 24.04 has integers of any size, so it computes H(50) exactly. Python limits recursion instead: sum_to(100_000) written in Python stopped with RecursionError: maximum recursion depth exceeded. C has no depth limit of its own. A recursion deep enough to use up the stack can crash the program, unless the compiler turns the recursion into a loop. With GCC 13.3 or Clang 18.1 on the same machine, sum_to(100000000) in C crashed with a segmentation fault at -O0 and printed 5000000050000000 at -O2.

walk-through

if (b == 0) { return a; } is the case that stops the recursion. Every recursive function needs one, and every chain of calls must reach it.

I64[2] h; declares an array of two numbers, both 0. h[0] = 0; changes nothing and is there so the starting fraction reads plainly, and h[1] = 1; makes it 0/1.

f[1] / gcd(f[1], k) * k divides before it multiplies. Written as f[1] * k / gcd(f[1], k), the product would overflow sooner.

num and den are local to each call of add_unit, so each call starts with its own num and den.

try this

  1. Print every H(k) instead of every tenth one. What is the last harmonic number each version of add_unit prints?
  2. Write lcm(a, b) as a function that calls gcd, and use it in add_unit.
  3. Rewrite sum_to with a loop and check that it gives 5050 for 100 and 5000050000 for 100,000.

reference

functions, parameters and views, recursion and depth, fault record, checked arithmetic

Next: lab 401-02, hunting an overflow.