Problem 330: Euler's Number

Detailed mathematical approach

Problem Summary

A sequence of real numbers \(a(n)\) is defined for every integer \(n\), the negative ones included, by

$$a(n)=\begin{cases}1 & n\lt 0\\ \sum_{i=1}^{\infty}\frac{a(n-i)}{i!} & n\ge 0\end{cases}$$

The infinite sum is harmless: as soon as \(i\gt n\) the index \(n-i\) is negative, so the term is simply \(1/i!\), and the series is a finite sum of earlier values plus a tail of the series for \(e\). The statement lists \(a(0)=e-1\), \(a(1)=2e-3\) and \(a(2)=\frac{7}{2}e-6\), notes that every value with \(n\ge 0\) can be written as

$$a(n)=\frac{A(n)e+B(n)}{n!}$$

with integers \(A(n)\) and \(B(n)\), and gives the checkpoint

$$a(10)=\frac{328161643\,e-652694486}{10!}.$$

We must compute

$$\bigl(A(10^9)+B(10^9)\bigr)\bmod 77\,777\,777.$$

The whole solution in one sentence: the requested quantity is just \(n!-A(n)\), and modulo each prime factor \(p\) of \(77\,777\,777\) the integers \(A(n)\) repeat with period \(p(p-1)\) once \(n\ge p\), which replaces the index \(10^9\) by one below \(p^2\). ^[1]

Mathematical Approach

Throughout we abbreviate \(A_n=A(n)\) and \(B_n=B(n)\).

1. Where the boundary bites

Fix \(n\ge 0\) and split the defining series at \(i=n\). The terms with \(1\le i\le n\) use values \(a(0),\dots,a(n-1)\) that are already known; the terms with \(i\gt n\) all have \(a(n-i)=1\). Since \(e=\sum_{i\ge 0}1/i!\), that second group is a tail of the series for \(e\):

$$\sum_{i\gt n}\frac{1}{i!}=e-\sum_{k=0}^{n}\frac{1}{k!}.$$

So the definition is really a finite recurrence with one copy of \(e\) attached,

$$\begin{aligned}a(n)&=\sum_{i=1}^{n}\frac{a(n-i)}{i!}\\ &\quad+e-\sum_{k=0}^{n}\frac{1}{k!},\end{aligned}$$

so each \(a(n)\) is well defined by induction on \(n\). At \(n=0\) the finite part is empty and the tail is \(e-1\), the value the statement gives.

2. Separating the two integer sequences

Substitute \(a(k)=(A_ke+B_k)/k!\) into the finite part and put \(k=n-i\). Multiplying the whole equation by \(n!\) and using \(n!/\bigl(k!(n-k)!\bigr)=\binom{n}{k}\) gives

$$\begin{aligned}A_ne+B_n&=\sum_{k=0}^{n-1}\binom{n}{k}(A_ke+B_k)\\ &\quad+n!\,e-\sum_{k=0}^{n}\frac{n!}{k!}.\end{aligned}$$

Collecting the coefficient of \(e\) and the rational part separately suggests the two integer recurrences

$$A_n=n!+\sum_{k=0}^{n-1}\binom{n}{k}A_k,$$

$$\begin{aligned}B_n&=\sum_{k=0}^{n-1}\binom{n}{k}B_k\\ &\quad-\sum_{k=0}^{n}\frac{n!}{k!}.\end{aligned}$$

Read in this direction nothing need be assumed about \(e\): define the integers \(A_n\), \(B_n\) by these recurrences, and the identity above shows by induction that \(a(n)=(A_ne+B_n)/n!\) really holds. Irrationality is needed only to know these are the statement's \(A(n)\) and \(B(n)\) and no others: if \(Ae+B=A'e+B'\) with integers and \(A\ne A'\), then \(e=(B'-B)/(A-A')\) would be rational, which Euler ruled out in 1737. ^[2]

3. The first few values

Running the two recurrences by hand is quick. The inner factorial sums \(\sum_{k\le n}n!/k!\) are \(1,2,5,16,65,326\) for \(n=0,\dots,5\), and the results are

• \(A_n\) for \(n=0,\dots,6\): \(1, 2, 7, 34, 211, 1596, 14259\);
• \(B_n\) for \(n=0,\dots,6\): \(-1, -3, -12, -62\), \(-398, -3072, -27798\).

For instance \(A_2=2!+\binom20 A_0+\binom21 A_1\)\({}=2+1+4=7\), while \(B_2\) is \(\binom20 B_0+\binom21 B_1-(2+2+1)\)\({}=-1-6-5=-12\), so \(a(2)=(7e-12)/2!=\frac72 e-6\), as the statement says. Continuing to \(n=10\) reproduces the statement's \(A_{10}=328161643\) and \(B_{10}=-652694486\), and a high-precision evaluation of the original series agrees with \((A_ne+B_n)/n!\) to forty decimal places for every \(n\le 12\).

4. The sum \(A_n+B_n\) collapses

Only the sum of the two sequences is wanted, so set \(C_n=A_n+B_n\). Adding the two recurrences of Section 2 gives

$$\begin{aligned}C_n&=n!-\sum_{k=0}^{n}\frac{n!}{k!}\\ &\quad+\sum_{k=0}^{n-1}\binom{n}{k}C_k.\end{aligned}$$

Claim. \(C_n=n!-A_n\) for every \(n\ge 0\).

Proof. For \(n=0\) both sides are \(0\). Assume \(C_k=k!-A_k\) for all \(k\lt n\). Because \(\binom{n}{k}k!=n!/(n-k)!\), substituting the hypothesis gives

$$\begin{aligned}&\sum_{k=0}^{n-1}\binom{n}{k}C_k\\ &=\sum_{j=1}^{n}\frac{n!}{j!}\\ &\quad-\sum_{k=0}^{n-1}\binom{n}{k}A_k,\end{aligned}$$

where \(j=n-k\) was substituted in the first sum. Since \(\sum_{k=0}^{n}n!/k!=n!+\sum_{j=1}^{n}n!/j!\), the two factorial sums cancel and only

$$C_n=n!-n!-\sum_{k=0}^{n-1}\binom{n}{k}A_k$$

survives. The recurrence for \(A_n\) says that the last sum equals \(A_n-n!\), so \(C_n=n!-A_n\), which closes the induction. \(\square\)

The values of \(C_n\) for \(n=0,\dots,6\) are \(0, -1, -5, -28\), \(-187, -1476, -13539\), which indeed equal \(n!-A_n\). From here on \(B_n\) can be forgotten: only \(A_n\) and \(n!\) matter, and \(n!\) is a multiple of every prime below \(n\).

5. From a convolution to a one-step recurrence

The recurrence for \(A_n\) is a convolution over all earlier indices, which is useless at \(n=10^9\). Generating functions turn it into a single step. Let \(f(x)=\sum_{n\ge 0}A_nx^n/n!\) be the exponential generating function. Including the \(k=n\) term in the binomial sum and rearranging gives

$$\sum_{k=0}^{n}\binom{n}{k}A_k=2A_n-n!\qquad(n\ge0).$$

The factorial denominators are what makes multiplication produce binomial coefficients. In \(e^xf(x)\), a term of degree \(n\) pairs degree \(k\) from \(f\) with degree \(n-k\) from \(e^x\), so the coefficient of \(x^n\) is

$$\sum_{k=0}^{n}\frac{A_k}{k!(n-k)!},$$

which multiplied by \(n!\) becomes the binomial sum on the left. ^[3] For \(n=2\) the coefficient is \(\frac12+2+\frac72=6\), and \(2!\cdot 6=12=2A_2-2!\). These are formal power series, so every coefficient involves only finitely many terms. On the right, \(\sum_{n}n!\,x^n/n!=1/(1-x)\), so \(e^xf(x)=2f(x)-1/(1-x)\), that is

$$f(x)=\frac{1}{(1-x)(2-e^{x})}.$$

This identifies the sequence as A337000. ^[8] Multiplying \(f\) by \(1-x\) leaves \(1/(2-e^x)\), whose combinatorial meaning drives all that follows. Let \(F_n\) count the ways to arrange \(n\) labelled items in a sequence of nonempty groups: the groups are ordered, the items inside a group are not. Two items can sit together, or in two singleton groups in either order, so \(F_2=3\); the empty arrangement gives \(F_0=1\). For \(n\ge 1\), choosing the \(k\) items of the first group and arranging the rest gives

$$F_n=\sum_{k=1}^{n}\binom{n}{k}F_{n-k}.$$

If \(G\) is the exponential generating function of these counts, the multiplication rule just proved turns this into \(G=1+(e^x-1)G\), that is \(G=1/(2-e^x)\). So \(F_n\) is the sequence \(1,1,3,13,75,541,4683\) of Fubini numbers, also called ordered Bell numbers. ^{[4] [5]} Hence \((1-x)f(x)=G\), and since the coefficient of \(x^n/n!\) in \(xf(x)\) is \(nA_{n-1}\), comparing coefficients yields

$$A_n=nA_{n-1}+F_n\qquad(n\ge 1),$$

with \(A_0=F_0=1\). A check: \(A_4=4\cdot 34+75=211\) and \(A_5=5\cdot 211+541=1596\). One index of \(A\) now costs one multiplication and one Fubini number.

6. The Fubini numbers modulo an odd prime

Fix a prime factor \(p\) of the modulus. All five are odd, so \(2\) is invertible modulo each: dividing by \(2\) means multiplying by that inverse, which is \(4\) modulo \(7\) because \(2\cdot 4\equiv 1\), and \(2^{-j}\) below means the inverse of \(2^{j}\).

Group the arrangements counted by \(F_n\) by their number \(k\) of groups. Forgetting the order of the groups leaves a partition into \(k\) nonempty unordered blocks, counted by the Stirling number of the second kind \(S(n,k)\), so there are \(k!\,S(n,k)\) arrangements with \(k\) groups. Equivalently, label the groups \(1,\dots,k\) and send each item to its group: that is a map hitting all \(k\) labels, a surjection. Inclusion–exclusion counts those by discarding the maps that miss labels — three items, two labels: \(2^3\) maps minus the \(2\) constant ones leaves \(6\) — and allowing only \(j\) labels gives \(j^n\) maps with \(\binom{k}{j}\) choices: ^[6]

$$k!\,S(n,k)=\sum_{j=0}^{k}(-1)^{k-j}\binom{k}{j}j^{n}.$$

Hence \(F_n=\sum_{k\ge 0}k!\,S(n,k)\), where the terms with \(k\gt n\) vanish. Modulo \(p\) the terms with \(k\ge p\) vanish as well, because \(k!\,S(n,k)\) is a multiple of \(k!\) and \(k!\) is a multiple of \(p\). Only \(k\le p-1\) is left, so exchanging the two sums gives

$$F_n\equiv\sum_{j=0}^{p-1}w_j\,j^{n}\pmod p,$$

$$w_j=\sum_{k=j}^{p-1}(-1)^{k-j}\binom{k}{j},$$

with the convention \(0^0=1\). The weights \(w_j\) do not depend on \(n\), and they have a closed form. Collected in a polynomial, interchanging the sums once more,

$$\sum_{j=0}^{p-1}w_jx^{j}=\sum_{k=0}^{p-1}(x-1)^{k}.$$

Multiplying the geometric sum on the right by \((x-1)-1=x-2\) gives \((x-1)^p-1\), which modulo \(p\) is \(x^p-2\) because \((u+v)^p\equiv u^p+v^p\). The same product appears in \(x^p-2^p=(x-2)\sum_{j=0}^{p-1}2^{\,p-1-j}x^{j}\), and \(2^p\equiv 2\) by Fermat's little theorem. ^[7] So over \(\mathbb{F}_p\), the residues modulo \(p\), both polynomials have the same product with \(x-2\). A product of non-zero polynomials over \(\mathbb{F}_p\) is non-zero, since leading coefficients multiply to something non-zero modulo a prime, so \(x-2\) simply cancels; no value is substituted for \(x\). Comparing coefficients gives \(w_j\equiv 2^{\,p-1-j}\equiv 2^{-j}\), so for every \(n\ge 0\)

$$F_n\equiv\sum_{j=0}^{p-1}\frac{j^{n}}{2^{\,j}}\pmod p.$$

A sum of \(p\) terms therefore replaces the quadratic recurrence. For \(p=7\) it gives \(F_n\equiv 1,1,3,6,5,2,0\) for \(n=0,\dots,6\), which is \(1,1,3,13,75,541,4683\) modulo \(7\).

7. The period of \(F_n\) is \(p-1\)

For \(n\ge 1\) the term \(j=0\) of that sum vanishes, and \(j^{p-1}\equiv 1\) for \(1\le j\le p-1\) by Fermat's little theorem. ^[7] Replacing \(n\) by \(n+p-1\) multiplies every surviving term by \(j^{p-1}\equiv 1\), so

$$\begin{aligned}F_{n+p-1}\equiv F_n\pmod p\\ (n\ge 1).\end{aligned}$$

The restriction \(n\ge 1\) cannot be dropped. Putting \(y=2^{-1}\), the identity \((y-1)\sum_{j\lt p}y^{j}=y^{p}-1\) together with \(y^p\equiv y\) gives \(\sum_{j=0}^{p-1}2^{-j}\equiv 1\), so \(F_{p-1}\equiv 1-1=0\), whereas \(F_0=1\). Index \(0\) is the single exception, and it is what shapes the final algorithm. For \(p=7\) the residues of \(F_n\) are \(1\), then \(1,3,6,5,2,0\) repeating from \(n=1\) on; a direct computation confirms that \(p-1\) is the minimal period for each of the five primes.

8. Multiples of \(p\) erase the past

Now look at \(A_n=nA_{n-1}+F_n\) modulo \(p\). When \(p\) divides \(n\) the factor \(n\) vanishes and the whole history is wiped out:

$$\begin{aligned}A_n\equiv F_n\pmod p\\ (p\mid n,\ n\ge 1).\end{aligned}$$

Write \(n=qp+r\) with \(0\le r\lt p\) and \(q\ge 1\). Unrolling the recurrence from the index \(qp\) upwards, the multiplier that carries the contribution of \(F_{qp+s}\) up to \(n\) is

$$\begin{aligned}&(qp+s+1)\cdots(qp+r)\\ &\equiv(s+1)\cdots r=\frac{r!}{s!}\pmod p,\end{aligned}$$

the empty product for \(s=r\) being \(1\). Therefore

$$A_n\equiv\sum_{s=0}^{r}\frac{r!}{s!}\,F_{qp+s}\pmod p.$$

Only \(r+1\le p\) Fubini numbers, at consecutive indices starting at a multiple of \(p\), decide the value; nothing before \(qp\) is relevant.

9. The period of \(A_n\) modulo \(p\)

Replace \(n\) by \(n+p(p-1)\). The remainder \(r\) is unchanged, since \(p(p-1)\) is a multiple of \(p\), and \(q\) grows by \(p-1\). Every index \(qp+s\) in the formula of Section 8 therefore grows by \(p(p-1)\), a multiple of \(p-1\), and every one of these indices is at least \(qp\ge p\ge 1\), so Section 7 applies to each of them and no Fubini value changes. The coefficients \(r!/s!\) do not change either, hence

$$\begin{aligned}A_{n+p(p-1)}\equiv A_n\pmod p\\ (n\ge p).\end{aligned}$$

The bound \(n\ge p\) is not cosmetic: with \(q=0\) the formula would reach index \(0\), the one place where Section 7 fails, and the first \(p\) residues really are different. For \(p=7\), with period \(42\), the values \(A_0,\dots,A_6\) are \(1,2,0,6,1,0,0\) modulo \(7\) while \(A_{42},\dots,A_{48}\) are \(0,1,5,0,5,6,1\) — different in every position. From \(n=7\) on they agree: \(A_7,\dots,A_{13}\) and \(A_{49},\dots,A_{55}\) are both \(1,4,0,5,1,5,3\). A direct search shows \(p(p-1)\) is in fact the minimal period for all five primes: \(42, 110, 5256, 10100, 18632\).

10. Reducing the index, and the modulus

The modulus factorises as

$$77\,777\,777=7\cdot 11\cdot 73\cdot 101\cdot 137,$$

five distinct odd primes, so the Chinese remainder theorem reconstructs a residue modulo \(77\,777\,777\) uniquely from the five residues modulo the factors. ^[9] For each prime and each \(n\ge p\), put

$$n'=p+\bigl((n-p)\bmod p(p-1)\bigr),$$

so that \(p\le n'\le p^2-1\) and \(A_n\equiv A_{n'}\). For \(n\ge p\) the factorial \(n!\) contains the factor \(p\) and vanishes modulo \(p\), so Section 4 reduces to

$$A_n+B_n\equiv -A_{n'}\pmod p,$$

while the tiny case \(n\lt p\) keeps \(n!-A_n\) at the unreduced index. The largest index ever reached is \(137^2-1=18768\), not \(10^9\).

11. A worked example: \(n=300\)

Take \(n=300\), large enough for two of the five primes to need the reduction. Listing \(p\), \(p(p-1)\), the reduced index \(n'\), the residue \(A_{n'}\bmod p\) and \(-A_{n'}\bmod p\):

• \(7\): \(42\), \(n'=48\), \(A_{48}\equiv 1\), residue \(6\);
• \(11\): \(110\), \(n'=80\), \(A_{80}\equiv 8\), residue \(3\);
• \(73\): \(5256\), \(n'=300\), \(A_{300}\equiv 39\), residue \(34\);
• \(101\): \(10100\), \(n'=300\), \(A_{300}\equiv 35\), residue \(66\);
• \(137\): \(18632\), \(n'=300\), \(A_{300}\equiv 37\), residue \(100\).

The five residues \((6,3,34,66,100)\) are merged one at a time, exactly as the code does it. Start with \(x=6\) modulo \(7\). For the second residue write \(x=6+7t\) and ask for \(x\equiv 3\pmod{11}\); as \(7^{-1}\equiv 8\pmod{11}\), this gives \(t\equiv(3-6)\cdot 8\equiv 9\), so \(x=69\) modulo \(77\). In general, with the primes used so far multiplying to \(M\), the next residue \(a\) modulo \(p\) needs \(t\equiv(a-x)M^{-1}\pmod p\) and replaces \(x\) by \(x+Mt\). The last three steps take \(t=46,58,65\) and pass through \(3611\) and \(329629\) to \(37231494\), which is also what \(300!-A_{300}\) gives in exact big integers. The same pipeline returns \(9410000\) at \(n=100\) and \(64356042\) at \(n=10\), the latter being \(10!-328161643\) modulo \(77\,777\,777\), so it confirms the statement's \(A_{10}+B_{10}\) as well.

How the Code Works

All three programs implement Sections 4, 9 and 10 directly, computing \(A_n\bmod p\) with the quadratic convolution of Section 2; the one-step recurrence of Section 5 only justifies the index reduction they rely on.

The entry points differ. Python's solve sets the index, modulus and prime list locally, calls solve_mod_prime once per prime and passes the five residues to crt_combine; Java's solve does the same from the constants DEFAULT_N, MODULUS and PRIME_FACTORS, with camelCase names. In C++, main parses its Options, runs the checkpoints unless switched off, then calls solve_mod, which gathers the residues and combines them.

reduced_index_for_prime is the three-line heart of the method: it returns \(n\) unchanged when \(n\lt p\), and otherwise \(p+((n-p)\bmod p(p-1))\), exactly the formula of Section 10.

compute_A_mod_prime(p, max_n) builds the table \(A_0,\dots,A_{n'}\) modulo \(p\). It keeps one row of Pascal's triangle in the array comb and advances it in place from right to left, so comb[k] holds \(\binom{n}{k}\bmod p\) while row \(n\) is processed; going backwards is what lets the update run without a second array. A variable fact tracks \(n!\bmod p\) beside it, multiplied by \(n\) while \(n\lt p\) and set to \(0\) from \(n=p\) on — Section 10 in one line. The inner loop adds \(\binom{n}{k}A_k\) for \(k\lt n\) to it and stores the result.

solve_mod_prime reduces the index, builds the table up to it, takes \(n!\bmod p\) from factorial_mod_small when \(n\lt p\) and \(0\) otherwise, and returns \(n!-A_{n'}\) modulo \(p\). crt_combine then merges the residues incrementally, as in Section 11: ^[9] it carries the partial solution x and the product current_modulus of the primes used so far, and inverts that product with inverse_mod_prime, which is Fermat's little theorem as \(a^{p-2}\equiv a^{-1}\). ^[7] The inverse exists because the earlier primes all differ from the current one.

The C++ program adds the self-checks. parse_arguments accepts --n= to run a different index and --skip-checkpoints to suppress the tests. compute_exact_ab evaluates both recurrences of Section 2 in exact 64-bit integers up to \(n=12\), with a full triangle of binomial coefficients, and run_checkpoints compares \(A_{10}\) and \(B_{10}\) with the statement, compares the modular pipeline with the exact one for every \(n\le 12\), and spot-checks Section 9 by building the whole table from \(A_0\) to \(A_{p^2+80}\), since \(p+p(p-1)=p^2\), and testing \(A_{p+t}\equiv A_{p+t+p(p-1)}\) for \(80\) offsets \(t\).

The implementations also differ in arithmetic and storage types. Python relies on arbitrary-precision integers and on the three-argument pow for the modular inverse; its normalize_mod and pow_mod helpers are never called. C++ and Java reduce after every product, which is safe because \(\binom{n}{k}A_k\) stays below \(p^2\le 137^2\) and the accumulator is 64-bit. Java holds the table in an int[] and casts the index bound to int, sound because Section 10 bounds it by \(p^2-1\le 18768\). C++ counts the Pascal row down with an unsigned type and stops at \(k\ge 1\), so the counter never wraps.

Complexity and Verification

The index reduction is one division per prime. Writing \(m=n'\), compute_A_mod_prime then performs \(m(m-1)/2\) Pascal-row updates and \(m(m+1)/2\) convolution steps, so exactly \(m^2\) inner iterations — an addition and a conditional subtraction in the first kind, a multiplication and a reduction in the second: time \(O(m^2)\), space \(O(m)\). With \(m\le p^2-1\) the worst case per prime is \(O(p^4)\); for \(p=137\) the extreme index \(18768\) gives \(352\,237\,824\) inner iterations, and the arrays values and comb hold \(18769\) residues each, stored as 32-bit integers in C++ and Java and as references to integer objects in Python. The reduction is the whole reason \(10^9\) is irrelevant to the cost.

As wall-clock time on an Apple M4 Max, with C++ built by Homebrew clang 22.1.6 at -O2 -std=c++17, OpenJDK 17 and CPython 3.11, each the median of repeated runs including start-up: C++ takes \(0.13\) s for the five residues and \(0.89\) s with the checkpoints, whose tables up to \(p^2+80\) dominate it; Java takes \(0.24\) s, Python \(3.16\) s. Sections 5 and 6 would allow \(O(m\,p)\) instead, maintaining the \(p\) powers \(j^n\) with one multiplication each per index; computing one period of Fubini residues in \(O(p^2)\) and then sweeping \(A_n=nA_{n-1}+F_n\) brings it to \(O(p^2+m)\). The convolution is already instantaneous here, so the programs keep it.

The packaged checks are those listed above: the statement's \(A_{10}\) and \(B_{10}\), agreement of the modular pipeline with the exact recurrences for \(0\le n\le 12\), and \(80\) periodicity offsets per prime. The last is a spot check, not a proof; Sections 6 to 9 supply the proof.

Every number printed above was recomputed for this article. The recurrences of Section 2 reproduce the statement's \(a(0)\), \(a(1)\), \(a(2)\) and \(a(10)\), and a 120-digit evaluation of the defining series agrees with \((A_ne+B_n)/n!\) for all \(n\le 12\). The identities \(A_n+B_n=n!-A_n\), \(A_n=nA_{n-1}+F_n\), \(A_n=\sum_{j\le n}(n!/j!)F_j\) and the first-group recurrence for \(F_n\) hold for every \(n\le 400\), as does the exponential sum of Section 6 modulo all five primes; the unrolled formula of Section 8 matches \(A_n\bmod p\) for all \(p\le n\le 4p^2\). The periods \(p-1\) from \(n=1\) and \(p(p-1)\) from \(n=p\) were confirmed over three full periods per prime and found minimal, and for every \(n\lt p\) the congruence \(A_n\equiv A_{n+p(p-1)}\) indeed fails. Finally a pipeline resting on Sections 5 and 6 alone, neither of which the published code uses, agrees with the exact big-integer \(A_n+B_n\) modulo \(77\,777\,777\) for every \(n\le 320\); the C++ program reproduces both worked examples of Section 11, and all three programs agree at the requested index.

Footnotes and References

The references below supply the background facts used above. The two integer recurrences, the identity \(A_n+B_n=n!-A_n\), the closed form of \(w_j\) and both periodicity statements are derived explicitly here.

1. Project Euler 330 — Euler's Number. The official statement sets a(n) = 1 for negative n and the infinite sum of a(n − i)/i! otherwise, lists a(0) = e − 1, a(1) = 2e − 3, a(2) = (7/2)e − 6 and a(10) = (328161643e − 652694486)/10!, writes a(n) = (A(n)e + B(n))/n! with integers A(n), B(n), and asks for A(10⁹) + B(10⁹) modulo 77777777.
2. Wikipedia — Proof that e is irrational. The page gives Euler's continued-fraction argument of 1737, published in 1744, and Fourier's proof, in which assuming e = a/b makes a scaled remainder of the series an integer strictly between 0 and 1. Section 2 needs only the conclusion.
3. Wikipedia — Binomial transform. The page defines the transform as a weighted sum with binomial coefficients and gives its effect on exponential generating functions: for its signed convention the transformed function is eˣ·f(−x), hence eˣ·f(x) for the unsigned sum of Section 5.
4. OEIS A000670 — Fubini numbers. The entry lists 1, 1, 3, 13, 75, 541, 4683, 47293, gives the exponential generating function 1/(2 − exp(x)), the recurrence a(n) = Σ_{k=1..n} binomial(n, k)·a(n − k) for n ≥ 1 with a(0) = 1, and the Stirling-number form a(n) = Σ k!·S(n, k) used in Section 6.
5. Wikipedia — Ordered Bell number. The page explains that these numbers count the weak orderings of an n-element set, derives the generating function 1/(2 − eˣ), and states both the Stirling sum Σ k!·S(n, k) and the double sum Σ_k Σ_j (−1)^(k−j)·C(k, j)·jⁿ that Section 6 starts from.
6. Wikipedia — Stirling numbers of the second kind. The page gives S(n, k) = (1/k!)·Σ (−1)^(k−j)·C(k, j)·jⁿ, derives it by inclusion–exclusion on the surjections from an n-set onto a k-set, and states that those surjections number k!·S(n, k) — the formula and the divisibility used in Section 6.
7. Wikipedia — Fermat's little theorem. For a prime p and any integer a, p divides aᵖ − a, and when a is coprime to p this says a^(p−1) ≡ 1 (mod p). Section 6 uses the first form for 2ᵖ ≡ 2, Section 7 the second for the period p − 1, the code the variant a^(p−2) ≡ a⁻¹.
8. OEIS A337000. The entry is the sequence A(n) of this problem: generating function 1/((1 − x)(2 − exp(x))), described as the exponential convolution of the Fubini numbers with the factorials, with a(n) = n!·Σ A000670(k)/k!; of its 21 terms, beginning 1, 2, 7, 34, 211, 1596, the one at index 10 is 328161643, the A(10) of the statement.
9. Wikipedia — Chinese remainder theorem. For pairwise coprime moduli the congruences always have a solution, and any two solutions agree modulo the product of the moduli. The page also gives the incremental construction, one congruence at a time with Bézout coefficients, which is what Section 11 and crt_combine do with a modular inverse.