y^2+z^2 = 9; x^2+z^2 = 9; with(linalg): i := vector([1,0,0]): j:=vector([0,1,0]): k:=vector([0,0,1]): r = x*i + y*j + sqrt(9-y^2)*k; z := sqrt(9-y^2): dS := sqrt(1+diff(z,x)^2+diff(z,y)^2); Answer := 2*Int(Int(dS,x=-abs(y)..abs(y)),y=-3..3) = 2*int(int(dS,x=-abs(y)..abs(y)),y=-3..3); evalf(rhs("));