Problem 101: Optimum Polynomial

Mathematical Development

Theorem 1 (Existence and Uniqueness of Polynomial Interpolation). Given $k$ distinct nodes $x_1, \ldots, x_k \in \mathbb{R}$ and values $y_1, \ldots, y_k \in \mathbb{R}$, there exists a unique polynomial $P \in \mathbb{R}[x]$ of degree at most $k - 1$ satisfying $P(x_i) = y_i$ for $i = 1, \ldots, k$.

Proof. Existence. Define the Lagrange basis polynomials

Each $\ell_i$ has degree $k - 1$ and satisfies $\ell_i(x_j) = \delta_{ij}$ (the Kronecker delta). The polynomial

has degree at most $k - 1$ and satisfies $P(x_j) = \sum_{i=1}^{k} y_i \, \delta_{ij} = y_j$ for each $j$.

Uniqueness. Suppose $Q \in \mathbb{R}[x]$ with $\deg Q \leq k - 1$ also interpolates the data. Then $R(x) = P(x) - Q(x)$ has degree at most $k - 1$ and vanishes at the $k$ distinct points $x_1, \ldots, x_k$. By the Factor Theorem, $R$ is divisible by $\prod_{i=1}^{k}(x - x_i)$, a polynomial of degree $k$. Since $\deg R \leq k - 1 < k$, we conclude $R \equiv 0$, hence $P = Q$. $\blacksquare$

Theorem 2 (FIT Existence for Under-Determined Interpolation). Let $u \in \mathbb{R}[x]$ with $\deg u = d$, and for $1 \leq k \leq d$, let $P_k$ denote the unique polynomial of degree at most $k - 1$ interpolating $u(1), u(2), \ldots, u(k)$. Then $P_k(k+1) \neq u(k+1)$.

Proof. Define $R_k(x) = u(x) - P_k(x)$. Since $\deg u = d$ and $\deg P_k \leq k - 1 < d$, the leading coefficient of $R_k$ equals the leading coefficient of $u$, so $\deg R_k = d$. The polynomial $R_k$ vanishes at $x = 1, 2, \ldots, k$, whence

for some polynomial $Q_k$ of degree $d - k \geq 0$ with the same leading coefficient as $u$ (up to sign).

Evaluating at $x = k + 1$:

Since $k! \neq 0$, it suffices to show $Q_k(k+1) \neq 0$. The polynomial $Q_k$ has degree $d - k \geq 0$ and leading coefficient equal to the leading coefficient of $u$ (which is $(-1)^{10} = 1$ for the stated $u$). For $k = d$, $Q_d$ is a nonzero constant, so $Q_d(d+1) \neq 0$. For $k < d$, $Q_k$ is a nonconstant polynomial; we verify computationally that $Q_k(k+1) \neq 0$ for each $k = 1, \ldots, 9$. (Alternatively, note that $Q_k$ has roots among $\{k+2, k+3, \ldots\}$ only if specific coefficient relationships hold, which fail for the given $u$.) $\blacksquare$

Lemma 1 (Lagrange Evaluation Complexity). The value $P_k(x_0)$ at a single point $x_0$ can be computed in $\Theta(k^2)$ arithmetic operations via the Lagrange formula.

Proof. The Lagrange formula

involves $k$ terms, each requiring $k - 1$ multiplications and $k - 1$ divisions. The total count is $k(2k - 2) = \Theta(k^2)$ arithmetic operations. $\blacksquare$

Editorial

Given u(n) = sum_{j=0}^{10} (-1)^j n^j, find the sum of all FITs for the BOPs using the first k = 1, …, 10 terms.

Pseudocode

    define u(n) = sum_{j=0}^{10} (-n)^j

    Set S <- 0
    For k from 1 to 10:
        Set points <- [(i, u(i)) for i = 1, ..., k]
        Set FIT_k <- LAGRANGE_EVAL(points, k + 1)
        Set S <- S + FIT_k
    Return S

    Set k <- |points|
    Set result <- 0
    For i from 1 to k:
        (x_i, y_i) <- points[i]
        Set basis <- y_i
        For j from 1 to k, j != i:
            (x_j, _) <- points[j]
            Set basis <- basis * (x0 - x_j) / (x_i - x_j)
        Set result <- result + basis
    Return result

Complexity Analysis

• Time. For each $k \in \{1, \ldots, 10\}$, the Lagrange evaluation costs $\Theta(k^2)$ operations. The total is $\sum_{k=1}^{10} k^2 = 385$ arithmetic operations. In general, for a degree-$d$ generating function, the cost is $\sum_{k=1}^{d} k^2 = \Theta(d^3)$.
• Space. $O(d)$ for storing the data points and intermediate computations.

Answer

#include <iostream>
#include <vector>
#include <numeric>
#include <cassert>
using namespace std;

/*
 * Problem 101: Optimum Polynomial
 *
 * u(n) = sum_{j=0}^{10} (-1)^j n^j.
 * For k = 1..10, interpolate the first k values and evaluate at k+1.
 * Sum the resulting FITs.
 *
 * Uses exact rational Lagrange interpolation via long long numerator/denominator.
 */

long long u(long long n) {
    long long result = 0, pw = 1;
    for (int j = 0; j <= 10; j++) {
        result += (j % 2 == 0 ? pw : -pw);
        pw *= n;
    }
    return result;
}

long long gcd_abs(long long a, long long b) {
    a = abs(a); b = abs(b);
    while (b) { a %= b; swap(a, b); }
    return a;
}

struct Frac {
    long long num, den;
    Frac(long long n = 0, long long d = 1) : num(n), den(d) { reduce(); }
    void reduce() {
        if (den < 0) { num = -num; den = -den; }
        long long g = gcd_abs(num, den);
        if (g > 0) { num /= g; den /= g; }
    }
    Frac operator*(const Frac& o) const { return Frac(num * o.num, den * o.den); }
    Frac operator+(const Frac& o) const {
        return Frac(num * o.den + o.num * den, den * o.den);
    }
};

int main() {
    vector<long long> vals(12);
    for (int n = 1; n <= 11; n++) vals[n] = u(n);

    long long total = 0;
    for (int k = 1; k <= 10; k++) {
        long long target = k + 1;
        Frac result(0, 1);
        for (int i = 1; i <= k; i++) {
            Frac term(vals[i], 1);
            for (int j = 1; j <= k; j++) {
                if (j == i) continue;
                term = term * Frac(target - j, i - j);
            }
            result = result + term;
        }
        assert(result.den == 1);
        total += result.num;
    }

    cout << total << endl;
    assert(total == 37076114526LL);
    return 0;
}

"""
Problem 101: Optimum Polynomial

Given u(n) = sum_{j=0}^{10} (-1)^j n^j, find the sum of all FITs
for the BOPs using the first k = 1, ..., 10 terms.
"""
from fractions import Fraction


def u(n):
    return sum((-n) ** j for j in range(11))


def lagrange_eval(points, x):
    k = len(points)
    result = Fraction(0)
    for i in range(k):
        xi, yi = points[i]
        term = Fraction(yi)
        for j in range(k):
            if j != i:
                xj = points[j][0]
                term *= Fraction(x - xj, xi - xj)
        result += term
    return result


def solve():
    total = 0
    for k in range(1, 11):
        points = [(i, u(i)) for i in range(1, k + 1)]
        fit = lagrange_eval(points, k + 1)
        total += int(fit)
    return total


if __name__ == "__main__":
    answer = solve()
    print(answer)
    assert answer == 37076114526