ring r = (2, A), x, lp;

minpoly = A^3 + A + 1;

// This is a comment
// Lagrange's interpolation formula
// Ni = numerator of the i'th sum term, Di = denominator
// See formula from my slides/book chapter

poly N1 = (x-1)*(x-A)*(x-A+1)*(x-A^2)*(x-A^2+1)*(x-A^2+A)*(x-A^2+A+1);
poly D1 = (0-1)*(0-A)*(0-A+1)*(0-A^2)*(0-A^2+1)*(0-A^2+A)*(0-A^2+A+1);
poly f1 = (N1/D1)*0;

poly N2 = (x-0)*(x-A)*(x-A+1)*(x-A^2)*(x-A^2+1)*(x-A^2+A)*(x-A^2+A+1);
poly D2 = (1-0)*(1-A)*(1-A+1)*(1-A^2)*(1-A^2+1)*(1-A^2+A)*(1-A^2+A+1);
poly f2 = (N2/D2)*0;

poly N3 = (x-0)*(x-1)*(x-A+1)*(x-A^2)*(x-A^2+1)*(x-A^2+A)*(x-A^2+A+1);
poly D3 = (A-0)*(A-1)*(A-A+1)*(A-A^2)*(A-A^2+1)*(A-A^2+A)*(A-A^2+A+1);
poly f3 = (N3/D3)*1;

poly N4 = (x-0)*(x-1)*(x-A)*(x-A^2)*(x-A^2+1)*(x-A^2+A)*(x-A^2+A+1);
poly D4 = (A+1-0)*(A+1-1)*(A+1-A)*(A+1-A^2)*(A+1-A^2+1)*(A+1-A^2+A)*(A+1-A^2+A+1);
poly f4 = (N4/D4)*1;

poly N5 = (x-0)*(x-1)*(x-A)*(x-A+1)*(x-A^2+1)*(x-A^2+A)*(x-A^2+A+1);
poly D5 = (A^2-0)*(A^2-1)*(A^2-A)*(A^2-A+1)*(A^2-A^2+1)*(A^2-A^2+A)*(A^2-A^2+A+1);
poly f5 = (N5/D5)*A;

poly N6 = (x-0)*(x-1)*(x-A)*(x-A+1)*(x-A^2)*(x-A^2+A)*(x-A^2+A+1);
poly D6 = (A^2+1-0)*(A^2+1-1)*(A^2+1-A)*(A^2+1-A+1)*(A^2+1-A^2)*(A^2+1-A^2+A)*(A^2+1-A^2+A+1);
poly f6 = (N6/D6)*A;

poly N7 = (x-0)*(x-1)*(x-A)*(x-A+1)*(x-A^2)*(x-A^2+1)*(x-A^2+A+1);
poly D7 = (A^2+A-0)*(A^2+A-1)*(A^2+A-A)*(A^2+A-A+1)*(A^2+A-A^2)*(A^2+A-A^2+1)*(A^2+A-A^2+A+1);
poly f7 = (N7/D7)*(A+1);

poly N8 = (x-0)*(x-1)*(x-A)*(x-A+1)*(x-A^2)*(x-A^2+1)*(x-A^2+A);
poly D8 = (A^2+A+1-0)*(A^2+A+1-1)*(A^2+A+1-A)*(A^2+A+1-A+1)*(A^2+A+1-A^2)*(A^2+A+1-A^2+1)*(A^2+A+1-A^2+A);
poly f8 = (N8/D8)*(A+1);

poly f = f1 + f2 + f3 + f4 + f5 + f6 + f7 + f8;

// C-like printf
printf("The interpolated polynomial is: %s", f);

quit;