k1=5;k2=3;k3=0;
P1=(1,3);P2=(3,2);P3=(5,k1);P4=(-1,k2);P5=(7, k3)
R1=[P1[0]^3,P1[0]^2,P1[0],-1]
R2=[P2[0]^3,P2[0]^2,P2[0],-1]
R3=[P3[0]^3,P3[0]^2,P3[0],-1]
R4=[P4[0]^3,P4[0]^2,P4[0],-1]
R5=[P5[0]^3,P5[0]^2,P5[0],-1]
M=matrix(QQ,5,4,[R1,R2,R3,R4,R5])
B=matrix(QQ,5,1,[P1[1],P2[1],P3[1],P4[1],P5[1]])
tu=transpose(M)*M
mi=tu^(-1)
yo=mi*transpose(M)*B
a0=float(yo[0][0]);a1=float(yo[1][0]);a2=float(yo[2][0]);a3=float(yo[3][0])
print B
print a0,a1,a2,a3
f(x)=a0*x^3+a1*x^2+a2*x-a3
tetu=plot(f(x),x,-4,8)
trit=points([P1,P2,P3,P4,P5],color='red')
tutuu=tetu+trit
tutuu.show()