'''
This code was written in 15-113 lecture on Tue 31-Jan and Thu 2-Feb.
It is only for demonstrational purposes, and may contain
dubious style and even an occasional bug.
'''

import numpy as np
from numpy.linalg import inv
from fractions import Fraction

def getA(p):
    maxPower = p+1
    rows = cols = maxPower+1
    A = [([0]*cols) for _ in range(rows)]
    for row in range(rows):
        for col in range(cols):
            A[row][col] = (row+1)**(cols-1-col)
    return A

def getB(p):
    maxPower = p+1
    rows = maxPower+1
    cols = 1
    B = [([0]*cols) for _ in range(rows)]
    colSum = 0
    for row in range(rows):
        colSum += (row+1)**p
        B[row][0] = colSum
    return B

def powerSumPolynomial(p):
    A = getA(p)
    b = getB(p)
    Ainv = inv(A)
    result = np.matmul(Ainv, b)
    return [convertToFraction(v[0]) for v in result]

def convertToFraction(f):
    return str(Fraction(f).limit_denominator(10**3))


print(powerSumPolynomial(4))