Question

Having a volume implicitly defined by

x*y*z <= 1 

for

-5 <= x <= 5 
-5 <= y <= 5 
-5 <= z <= 5 

how would I go about plotting its outer surface using available Python modules, preferably mayavi?

I am aware of the function mlab.mesh, but I don't understand its input. It requires three 2D arrays, that I don't understand how to create having the above information.

EDIT:

Maybe my problem lies with an unsufficient understanding of the meshgrid()-function or the mgrid-class of numpy. I see that I have to use them in some way, but I do not completely grasp their purpose or what such a grid represents.

EDIT:

I arrived at this:

import numpy as np

from mayavi import mlab
x, y, z = np.ogrid[-5:5:200j, -5:5:200j, -5:5:200j]
s = x*y*z

src = mlab.pipeline.scalar_field(s)

mlab.pipeline.iso_surface(src, contours=[1., ],)
mlab.show()

This results in an isosurface (for x*y*z=1) of a volume though, which is not quite what I was looking for. What I am looking for is basically a method to draw an arbitrary surface, like a "polygon in 3d" if there is such a thing.

I created the following code, which plots a surface (works with mayavi, too). I would need to modify this code to my particular problem, but to do that I need to understand why and how a 3d surface is defined by three 2d-arrays? What do these arrays (x, y and z) represent?

import numpy as np
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import axes3d, Axes3D

phi, theta = np.mgrid[0:np.pi:11j, 0:2*np.pi:11j]
x = np.sin(phi) * np.cos(theta)
y = np.sin(phi) * np.sin(theta)
z = np.cos(phi)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_wireframe(x,y,z)
fig.show()
Was it helpful?

Solution

The outer surface, implicitly defined by

x*y*z = 1,

cannot be defined explicitly globally. To see this, consider x and y given, then:

z = 1/(x*y),

which is not defined for x = 0 or y = 0. Therefore, you can only define your surface locally for domains that do not include the singularity, e.g. for the domain

0 < x <= 5
0 < y <= 5

z is indeed defined (a hyperbolic surface). Similarly, you need to plot the surfaces for the other domains, until you have patched together

-5 <= x <= 5
-5 <= y <= 5

Note that your surface is not defined for x = 0 and y = 0, i.e. the axis of your coordinate system, so you cannot patch your surfaces together to get a globally defined surface.

Using numpy and matplotlib, you can plot one of these surfaces as follows (adopted from http://matplotlib.org/mpl_toolkits/mplot3d/tutorial.html#surface-plots):

from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import matplotlib.pyplot as plt
import numpy as np

fig = plt.figure()
ax = fig.gca(projection='3d')

X = np.arange(0.25, 5, 0.25)
Y = np.arange(0.25, 5, 0.25)
X, Y = np.meshgrid(X, Y)
Z = 1/(X*Y)

surf = ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.coolwarm,
        linewidth=0, antialiased=False)
ax.set_zlim(0, 10)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')    
fig.colorbar(surf, shrink=0.5, aspect=5)
plt.show()

I'm not familiar with mayavi, but I would assume that creating the meshes with numpy would work the same.

OTHER TIPS

The test case in the Mayavi docs where the function test_mesh() is defined is capable of producing a sphere. This is done by replacing

r = sin(m0*phi)**m1 + cos(m2*phi)**m3 + sin(m4*theta)**m5 + cos(m6*theta)**m7

with r = 1.0 say.

However, your problem is you need to understand that the equations you are writing define a volume when you want to draw a sphere. You need to reformulate them to give a parametric equation of a sphere. This is essentially what is done in the above example, but it may be worth your while to try it yourself. As a hint consider the equation of a circle and extend it.

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top