restart;
P:=z->a*(z-p2)^2+b*(z-p2)+fp2;
eq1:=P(p1)=fp1;
eq2:=P(p0)=fp0;
S:=solve({eq1,eq2},{a,b});
simplify(S[1]);
simplify(S[2]);
with(CodeGeneration);
C(S[1]);
C(S[2]);
plot(z^3-2*z^2-3*z+10,z=-3..4);