ACM ICPC Dhaka Regional 2008: Stopping Doom's Day

ACM ICPC Dhaka Regional 2008: Stopping Doom's Day

2017, May 05    

Tags: Polynomial Interpolation, Gaussian Elimination, Modular Inverse

Problem Statement:

Solution: One can probably find a closed form for this. But this is not necessary and also it can cause a lot of time-waste during a live contest. Instead, one can notice that, the term inside the summands is of degree 5. So, after all the summations are done, the final value is a 10 degree polynomial of $n$.
So, we can assume,

$T = f(n) = a_{0} + a_{1}n^{1} + a_{2}n^{2} + a_{3}n^{3} + a_{4}n^{4} + a_{5}n^{5} + a_{6}n^{6} + a_{7}n^{7} + a_{8}n^{8} + a_{9}n^{9} + a_{10}n^{10}$

We can manually calculate the value of $f(n)$ for $n = 1, 2, \cdots, 11$. Then, it comes down to solving a system of 11 equations with 11 variables. We solve for the coefficients and use them to determine $f(n)$ for future values of $n$.
Implementation:

#include <bits/stdc++.h>
#define MOD 10007

using namespace std;



#define MAX         3
/// n rows, m columns
/// MAX is the maximum number of equations or variables
/// X holds the solution
/// returns 1 if unique solution, -1 if inconsistent, inf if infinite solutions
struct matrix{
    long long arr[MAX+10][MAX+10];
    int n, m;
    void print()
    {
        for(int i = 0; i<n; i++)
        {
            for(int j = 0; j<=m; j++)
                cout << arr[i][j] << " ";
            cout << endl;
        }
    }
}A;

int where[MAX+10];
long long X[MAX+10];
const int inf = INT_MAX;


long long modinverse(long long a, long long n){
    if(n==0) return 1LL;
    long long ret = modinverse(a,n/2);
    ret = (ret*ret) % MOD;
    if(n%2==1){
        ret = (ret*a) % MOD;
    }
    return ret;
}


int gauss()
{
    memset(where,-1,sizeof(where));
    int row, col;
    for(row = col = 0; row<A.n && col<A.m; col++)
    {
        int pivot = row;
        for(int i = row+1; i<A.n; i++)
        {
            if(abs(A.arr[pivot][col]) < abs(A.arr[i][col]))
                 pivot = i;
        }
        if(pivot != row)
        {
            for(int i = 0; i<=A.m; i++)
                swap(A.arr[row][i], A.arr[pivot][i]);
        }
        if(A.arr[row][col]==0)
            continue;
        where[col] = row;
        for(int i = row+1; i<A.n; i++)
        {
            if(A.arr[i][col])
            {
                long long c = (A.arr[i][col]*modinverse(A.arr[row][col],MOD-2))%MOD;
                for(int j = col; j<=A.m; j++){
                    A.arr[i][j] -= (c*A.arr[row][j])%MOD;
                    A.arr[i][j] = (A.arr[i][j]+MOD)%MOD;
                }
            }
        }
        row++;
    }


    for(int i = 0; i<A.n; i++)
    {
        long long total = 0;
        for(int j = 0; j<A.m;j++){
            total += abs(A.arr[i][j]);
            total %= MOD;
        }
        if(abs(total)==0 && abs(A.arr[i][A.m])==0)
            return -1;
    }
    for(int i = A.n-1; i>=0; i--)
    {
        if(where[i] == -1)  return inf;
        int row = where[i];
        long long sltn = A.arr[row][A.m];
        for(int j = i+1; j<A.m; j++){
            sltn -= (A.arr[row][j]*X[j])%MOD;
            sltn %= MOD;
            if(sltn < 0) sltn+= MOD;
        }
        X[i] = (sltn *modinverse(A.arr[row][i],MOD-2))%MOD;
    }
    return 1;
}


void buildmatrix(){
    A.n = 11, A.m = 11;
    for(int i = 0; i < A.n; i++){
        A.arr[i][0] = 1LL;
        for(int j = 1; j < A.m; j++){
            A.arr[i][j] = (A.arr[i][j-1] * (i+1)) %MOD;
        }
    }
    for(int n = 1; n <= 11; n++){
        for(int i = 1; i <= n; i++){
            for(int j = 1; j <= n; j++){
                for(int k = 1; k <= n; k++){
                    for(int l = 1; l <= n; l++){
                        for(int m = 1; m <= n; m++){
                            A.arr[n-1][11] += ((long long)abs(i-j)*abs(j-k)*abs(k-l)*abs(l-m)*abs(m-i))%MOD;
                            A.arr[n-1][11] %= MOD;
                        }
                    }
                }
            }
        }
    }
}


long long Fun(long long n){
    long long ret = X[10];
    for(int i = 9; i >= 0; i--){
        ret = (ret * n )%MOD+ X[i];
        ret %= MOD;
    }
    return ret;
}


int main(){
    buildmatrix();
    int x = gauss();
    long long n;
    while(1){
        scanf("%lld",&n);
        if(n==0) break;
        printf("%lld\n",Fun(n));
    }
    return 0;
}