// Create Vandermonde Matrix
#include <stdio.h>
#include <stdlib.h>
#include <time.h>


#define LONG long long int
#define ULONG unsigned long long int


void VandermondeMatrix(LONG ** A, int n, LONG * x, LONG p){

   for(i=0;i < n; i++){
      t=1;
         for(j=0;j<n;j++){
             A[i*n+j] = t;
             t = t*x[j] mod p;
         }
   }
}

LONG det64s( LONG * A, int n, LONG p ) {
// M = [[a,b,c],[d,e,f],[1,2,3],[4,5,6]]  M[i][j]
// A = [a,b,c,d,e,f,1,2,3,4,5,6]          A[m*i+j]   m = #columns
    int i,j,k;
    LONG d,t;
    recint P;
    P = recip1(p);
    d = 1;
    for( k=0; k<n; k++ ) {
        //if( k>0 && k%100 == 0 ) printf("elimination at row %d\n",k);
        for( i=k; i<n && A[n*i+k]==0; i++ ); // look for non-zero pivot
        if( i>=n ) { d = 0; break; }
        if( i!=k ) { // interchange row k with row i
             for( j=k; j<n; j++ ) { t = A[k*n+j]; A[k*n+j] = A[i*n+j]; A[i*n+j] = t; }
             d = neg64s(d,p);
        };
        d = mulrec64(d,A[k*n+k],P);
        if( A[k*n+k]==0 ) { printf("division by 0\n"); return -1; }
        t = inv64s(A[k*n+k],p);
        A[n*k+k] = 1;
        for( j=k+1; j<n; j++ ) { A[n*k+j] = mulrec64(A[n*k+j],t,P); }
        for( i=k+1; i<n; i++ ) { // row i
            if( A[i*n+k]!=0 )
                for( j=k+1; j<n; j++ )
                     // A[i,j] = A[i,j] - A[i,k] A[k,j] mod p
                     A[i*n+j] = sub64s(A[i*n+j],mulrec64(A[i*n+k],A[k*n+j],P), p );
            A[i*n+k] = 0;
        }
    }
    return d;
}