Problem: Professor Abacus has just built a new computing engine for making numerical tables. It was designed to calculate the values of a polynomial function in one variable at several points at a time. With the polynomial function f(x) = x2 + 2x + 1, for instance, a possible expected calculation result is 1 (= f(0)), 4 (= f(1)), 9 (= f(2)), 16 (= f(3)), and 25 (= f(4)).

It is a pity, however, the engine seemingly has faulty components and exactly one value among those calculated simultaneously is always wrong. With the same polynomial function as above, it can, for instance, output 1, 4, 12, 16, and 25 instead of 1, 4, 9, 16, and 25.

You are requested to help the professor identify the faulty components. As the first step, you should write a program that scans calculation results of the engine and finds the wrong values.



Input

The input is a sequence of datasets, each representing a calculation result in the following format.

d
v0
v1
...
vd+2


Here, d in the first line is a positive integer that represents the degree of the polynomial, namely, the highest exponent of the variable. For instance, the degree of 4x5 + 3x + 0.5 is five and that of 2.4x + 3.8 is one. d is at most five.

The following d + 3 lines contain the calculation result of f(0), f(1), ... , and f(d + 2) in this order, where f is the polynomial function. Each of the lines contains a decimal fraction between -100.0 and 100.0, exclusive.

You can assume that the wrong value, which is exactly one of f(0), f(1), ... , and f(d+2), has an error greater than 1.0. Since rounding errors are inevitable, the other values may also have errors but they are small and never exceed 10-6.

The end of the input is indicated by a line containing a zero.

Output

For each dataset, output i in a line when vi is wrong.

Example

Input

2
1.0
4.0
12.0
16.0
25.0
1
-30.5893962764
5.76397083962
39.3853798058
74.3727663177
4
42.4715310246
79.5420238202
28.0282396675
-30.3627807522
-49.8363481393
-25.5101480106
7.58575761381
5
-21.9161699038
-48.469304271
-24.3188578417
-2.35085940324
-9.70239202086
-47.2709510623
-93.5066246072
-82.5073836498
0


Output

2
1
1
6

Solution (in CPP):

Solution:
#include <cassert>// c
#include <ctime>
#include <iostream>// io
#include <iomanip>
#include <fstream>
#include <sstream>
#include <vector>// container
#include <map>
#include <set>
#include <queue>
#include <bitset>
#include <stack>
#include <algorithm>// other
#include <complex>
#include <numeric>
#include <functional>
#include <random>
#include <regex>

using namespace std;
typedef long long ll;

#define ALL(c) (c).begin(),(c).end()
#define FOR(i,l,r) for(int i=(int)l;i<(int)r;++i)
#define REP(i,n) FOR(i,0,n)
#define FORr(i,l,r) for(int i=(int)r-1;i>=(int)l;--i)
#define REPr(i,n) FORr(i,0,n)
#define EACH(i,c) for(__typeof((c).begin()) i=(c).begin(); i!=(c).end(); ++i)
#define IN(l,v,r) ((l)<=(v) && (v)<(r))
#define UNIQUE(v) v.erase(unique(ALL(v)),v.end())
//debug
#define DUMP(x)  cerr << #x << " = " << (x)
#define LINE()    cerr<< " (L" << __LINE__ << ")"

template<typename T,typename U> T pmod(T x,U M){return (x%M+M)%M;}

typedef long double D;
const D INF = 1e12,EPS = 1e-8;
int sig(D a,D b=0){return a<b-EPS?-1:a>b+EPS?1:0;}

#define H(mat) mat.size()
#define W(mat) mat[0].size()


class Ambiguous{};//many
class NoSolution{};//none

template <typename T> class Mat :public vector<vector<T>> {
	public:

	// tropical
	// static const T Ae=INF;
	// T A(T a,T b){return min(a,b);}
	// static const T Me=0;
	// T M(T a,T b){return a+b;}
	
	// double
	static constexpr T Ae=0.0;
	T A(T a,T b){return a+b;}T Ainv(T a){return -a;}
	static constexpr T Me=1.0;
	T M(T a,T b){return a*b;}T Minv(T a){return 1/a;}
	// mod
	// const int MD=2;
	// static constexpr T Ae=0;
	// T A(T a,T b){return pmod(a+b,MD);}T Ainv(T a){return pmod(-a,MD);}
	// static constexpr T Me=1;
	// T M(T a,T b){return pmod(a*b,MD);}T Minv(T a){return pow(a);/*mod inv */}
	T S(T a,T b){return a+Ainv(b);}

	Mat(int H,int W) : vector<vector<T>>(H){
		REP(y,H)this->at(y)=vector<T>(W);
	}
	Mat(int H,int W,int v) : vector<vector<T>>(H){
		REP(y,H)this->at(y)=vector<T>(W,v);
	}
	const Mat I(int H,int W){
		Mat A(H,W);
		REP(y,H)REP(x,W)A[y][x]=Ae;
		REP(i,min(H,W))A[i][i]=Me;	
		return A;	
	}
	//O(n^2)
	Mat operator+(const Mat& r){
		const Mat& l=*this;
		assert(H(r)==H(l) && W(l)==W(r));
		Mat res(H(l),W(r));
		REP(y,H(l))REP(x,W(r)) res[y][x]=A(l[y][x],r[y][x]);
		return res;
	}
	//O(n^2)
	Mat operator-(const Mat& r){
		const Mat& l=*this;
		assert(H(r)==H(l) && W(l)==W(r));
		Mat res(H(l),W(r));
		REP(y,H(l))REP(x,W(r)) res[y][x]=S(l[y][x],r[y][x]);
		return res;
	}
    //O(LH*N*RW) O  (n^3)
	Mat operator*(const Mat& r){
		const Mat& l=*this;
		assert(H(r)==W(l));
		Mat res(H(l),W(r));
		REP(y,H(l))REP(x,W(r))
			REP(k,H(r)) res[y][x]=A(res[y][x],M(l[y][k],r[k][x]));
		return res;
	}
	
	//O(n^2)
	vector<T> operator*(const vector<T>& r){
		const Mat& l=*this;
		vector<T> res(r.size());
		REP(y,H(l))REP(x,W(l))res[y]=A(res[y],M(l[y][x],r[x]));
		return res;
	}

	// O(n^3*log(v))
	// Mat x must be n*n size
	Mat pow(T v){
		Mat x=*this,res=I(H(x),W(x));

		while(v!=0){
			if(v%2)res=res*x;
			x=x*x;
			v/=2;
		}
		return res;
	}

	// //verified by SRM306 LightSwitches
	// //次数を求める
	// //O(n^3)
	// template<typename T>
	// int rank(Mat<T> A){
	// 	 int H = H(A),W=W(A);
	// 	int py =0,px = 0;
	// 	while(py < H && px < W){
	// 		for(int y=py+1;y<H;y++){ // pivot
	// 			if(abs(A[y][px]) > abs(A[py][px])){
	// 				swap(A[y],A[py]);swap(b[y],b[py]);
	// 			}
	// 		}
	// 		if(abs(A[py][px]) > 0){
	// 			T d = inv(A[py][px]);
	// 			for(int x=0;x<W;x++) A[py][x] = mod(A[py][x] * d);
	// 			for(int y=py+1;y<H;y++){
	// 				T k = A[y][px];
	// 				for(int x=0;x<W;x++) A[y][x] = mod(A[y][x] - k * A[py][x]);
	// 			}
	// 			py++;
	// 		}
	// 		px++;
	// 	}
	// 	return py;
	// }
		
	// //連立方程式
	//O(n^3)
	//Ax = b
	void gauss(vector<T>& b){
		Mat A=*this;
		const int H = H(A),W=W(A);
		int py =0,px = 0;
		while(py < H && px < W){
			FOR(y,py+1,H){ // pivot
				if(sig(abs(A[y][px]),abs(A[py][px]))>0){
					swap(A[y],A[py]);swap(b[y],b[py]);
				}
			}
			if(sig(A[py][px]) != 0){
				T d = Minv(A[py][px]);
				FOR(x,0,W) A[py][x] =M(A[py][x],d);
				b[py] = M(b[py],d);
				FOR(y,py+1,H){
					T k = A[y][px];
					FOR(x,0,W) A[y][x] = S(A[y][x],M(k,A[py][x]));
					b[y] = S(b[y],M(k,b[py]));
				}
				py++;
			}
			px++;
		}
		FOR(y,py,H)if(sig(b[y]) != 0)throw NoSolution();
		if(py < W || px < W)throw Ambiguous();

		FORr(x,0,W)FOR(y,0,x){
			b[y] = S(b[y],M(b[x],A[y][x]));
		}
	}


	// 行列A の逆行列を求める．
	// 制約:N*N 行列
	// O(n^3)
	Mat inv(){
		Mat A=*this;int n=H(A);Mat B=I(n,n);
		REP(i,n){
			T d = A[i][i];
			REP(j,n){
				A[i][j]=M(A[i][j],Minv(d)); B[i][j]=M(B[i][j],Minv(d)); 
			}
			REP(r,n)if(r!=i){
				T dif = A[r][i];
				REP(j,n){
					A[r][j]=S(A[r][j],M(dif,A[i][j]));
					B[r][j]=S(B[r][j],M(dif,B[i][j]));
				}
			}
		}
		return B;
	}
};
int main() {
	cerr << fixed <<setprecision(20);
	cout << fixed <<setprecision(20);
  while (true) {
	int d;cin >> d;if(d==0)break;

	vector<D> ps(d+3);
	REP(i,d+3)cin >> ps[i];

	vector<int> par(d+3);
	bool ex=false;
	REP(i,d+3)REP(j,d+3)if(!ex)if(i!=j){
		Mat<D> mat(d+1,d+1);
		{int id=0;REP(k,d+3)if(k!=i && k!=j)par[id++]=k;}

		REP(y,d+1)REP(x,d+1)mat[y][x]=pow(par[y],x);
		vector<D> bs(d+1);REP(y,d+1)bs[y]=ps[par[y]];

		try{
			mat.gauss(bs);
			bool ok=false;int ind=-1;
			REP(x,d+3){
				D v=0;
				REP(r,d+1)v+=bs[r]*pow((D)x,r);
				if(abs(v-ps[x])>1e-3){
					if(!ok && ind==-1){
						ok=true;
						ind=x;
					}else{
						ok=false;
					}
				}
			}
			if(ok){
				cout << ind <<endl;
				ex=true;
			}
		}catch(NoSolution e){}
		catch(Ambiguous e){}
	}
	if(!ex)assert(false);
  }
}