var('p')
w=vector([0,0,p])
print w
var('x')
var('y')
var('z')
r=vector([x , y , z])
print r
w.cross_product(r)
p=10
f=(x**2+y**2)**(1/2)
print f
plot3d(p*f, (x, -100, 100), (y,-100,100))