var('x,y,x1,y1,z,t')
f(x,y)=(24-x^2-y^2)/8
a=implicit_plot3d(f(x,y)==z,(x,-10,10),(y,-10,10),(z,-10,10),color='gold')
x1=2;y1=4
po=point([x1,y1,f(x1,y1)],size=20,color='red')
fx(x,y)=diff(f(x,y),x)
fy(x,y)=diff(f(x,y),y)
print("fx=",fx(x,y))
print("fy=",fy(x,y))
print("fx=",fx(x1,y1))
print("fy=",fy(x1,y1))
plane(x,y)=fx(x1,y1)*(x-x1)+fy(x1,y1)*(y-y1)+f(x1,y1)
plot1=parametric_plot3d([x1+t,y1,f(x1,y1)+t*fx(x1,y1)],(t,-10,10),thickness=5,color='maroon')
plot2=parametric_plot3d([x1,y1+t,f(x1,y1)+t*fy(x1,y1)],(t,-10,10),thickness=5,color='maroon')
pl=implicit_plot3d(plane(x,y)==z,(x,-10,10),(y,-10,10),(z,-10,10),color='gray')
a+po+pl+plot1+plot2
