R2=matrix([[4,-1,0],[-1,4,-1],[0,-1,4]])
print R2
p2=vector([21,-45,33])
x2=R2.solve_right(p2)
print x2