f := 3*z^2*cos(x*z)+x^2/(y^2+sin(y)^2*z); Dfx := Diff(f,x) = diff(f,x); fx := unapply(op(Dfx)[2],x,y,z); `fx(1,1,0)` = fx(1,1,0); `fx(1,1,0)` = subs({x=1,y=1,z=0},op(Dfx)[2]); Diff('f',x,y$2) = diff(f,x,y$2); f1 := (x,y) -> 16-4*x^2-y^2; g := unapply(f1(x,2),x); fx := D(g); `fx(1)` = fx(1); gx := subs({x=2,y=1},diff(f1(x,y),x)); S := plot3d(f1(x,y),x=-5..5,y=-5..5,axes=framed): with(plots): L1 := spacecurve(evalm([1,2,8]+t*[1,0,-8]),t=-5..1,color=yellow): L2 := spacecurve(evalm([2,1,-1]+t*[1,0,-16]),t=-10..1/16,color=yellow): display3d({S,L1,L2},axes=framed); u := (x,t) -> exp(-a^2*k^2*t)*sin(k*x); # ut := diff(u(x,t),t); Rhs := a^2*diff(u(x,t),x$2); ut - Rhs;