>	restart:with(plots):with(linalg):

Warning, the name changecoords has been redefined

Warning, the protected names norm and trace have been redefined and unprotected

>	with(VectorCalculus):SetCoordinates('cartesian'[x,y,z]):

>	C:=NULL:B:=NULL:F:=NULL:for n from 1 to 9 do xval:=sqrt(10)*(n-1)/16:yval:=sqrt(10)*(n-1)/16:zval:=xval^2+yval^2+3:tan1:=TangentPlane((x,y)->x^2+y^2+3, xval,yval):pt1:=subs(x=0,y=0,tan1(x,y)):pt2:=subs(x=1,y=0,tan1(x,y)):pt3:=subs(x=0,y=1,tan1(x,y)):

>	v1:=convert(pt1,vector):v2:=convert(pt2,vector):v3:=convert(pt3,vector):

>	c1:=crossprod(evalm(v3-v1),evalm(v2-v1)):C:=C,plot3d({x^2+y^2+3,c1[1]*x+c1[2]*y-(c1[1]*xval+c1[2]*yval+c1[3]*zval)},x=-3..3,y=-3..3,view=-16..16,axes=normal,orientation=[20,75]):

xval2:=sqrt(5)*cos(Pi/4-Pi*(n-1)/32):yval2:=sqrt(5)*sin(Pi/4-Pi*(n-1)/32):zval2:=xval2^2+yval2^2+3:tan2:=TangentPlane((x,y)->x^2+y^2+3, xval2,yval2):pt4:=subs(x=0,y=0,tan2(x,y)):pt5:=subs(x=1,y=0,tan2(x,y)):pt6:=subs(x=0,y=1,tan2(x,y)):

>	v4:=convert(pt4,vector):v5:=convert(pt5,vector):v6:=convert(pt6,vector):

>

c2:=crossprod(evalm(v6-v4),evalm(v5-v4)):B:=B,plot3d({x^2+y^2+3,c2[1]*x+c2[2]*y-(c2[1]*xval2+c2[2]*yval2+c2[3]*zval2)},x=-3..3,y=-3..3,view=-16..16,axes=normal,orientation=[20,75]):xval3:=sqrt(5)*(8-n)/8:yval3:=0:zval3:=xval3^2+yval3^2+3:tan3:=TangentPlane((x,y)->x^2+y^2+3, xval3,yval3):pt7:=subs(x=0,y=0,tan3(x,y)):pt8:=subs(x=1,y=0,tan3(x,y)):pt9:=subs(x=0,y=1,tan3(x,y)):

>	v7:=convert(pt7,vector):v8:=convert(pt8,vector):v9:=convert(pt9,vector):

>	c3:=crossprod(evalm(v9-v7),evalm(v8-v7)):F:=F,plot3d({x^2+y^2+3,c3[1]*x+c3[2]*y-(c3[1]*xval3+c3[2]*yval3+c3[3]*zval3)},x=-3..3,y=-3..3,view=-16..16,axes=normal,orientation=[20,75]):od:display(C,B,F,insequence=true);

>

>

>