Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Question: plot_trisurf (matplotlib) directly from qhull #343

Closed
orbisvicis opened this issue May 17, 2022 · 3 comments
Closed

Question: plot_trisurf (matplotlib) directly from qhull #343

orbisvicis opened this issue May 17, 2022 · 3 comments

Comments

@orbisvicis
Copy link

I want to plot the Delaunay triangulation directly with plot_trisurf from mplot3d because the grid interpolation consumes too much memory, but the resulting surface is a jumble of triangles even though learner.npoints is 71016.

interp = learner.interpolator()
qhull = interp.tri
mptri = mpl.tri.Triangulation\
    ( x = qhull.points[:, 0]
    , y = qhull.points[:, 1]
    , triangles = qhull.vertices
    )
fig = plt.figure()
ax = plt.axes(projection="3d")
ax.plot_trisurf\
    (mptri, interp.values[:, 0], cmap="viridis")

The issue appears related to the way adaptive triangulates the mesh, because I can create a triangulation from a regular grid without any issues:

x = np.linspace(-100, 100, 1000)
y = np.linspace(-100, 100, 1000)
X, Y = np.meshgrid(x, y)
Z = func(X, Y, ...)
xr = np.ravel(X)
yr = np.ravel(Y)
zr = np.ravel(Z)
interp = scipy.interpolate.LinearNDInterpolator\
    (np.vstack([xr, yr]).T, zr)
qhull = interp.tri
mptri = mpl.tri.Triangulation\
    ( x = qhull.points[:, 0]
    , y = qhull.points[:, 1]
    , triangles = qhull.vertices
    )
fig = plt.figure()
ax = plt.axes(projection="3d")
ax.plot_trisurf\
    (mptri, interp.values[:, 0], cmap="viridis")

It's interesting that the plot and interpolated_on_grid results return the expected surface, except that the output is transposed.

def func_wrapper(xy):
    return func(xy[0], xy[1], ...)

learner = adaptive.Learner2D\
    (test_func, bounds=[(-100,100), (-100, 100)])
runner = adaptive.BlockingRunner\
    (learner, goal = lambda l: l.loss() < 0.00001)
x, y, Z = learner.interpolated_on_grid(4000)
# Un-transpose the output
Z = Z.T
X, Y = np.meshgrid(x, y)
fig = plt.figure()
ax = plt.axes(projection="3d")
ax.plot_surface(X, Y, Z, cmap="viridis")
@orbisvicis
Copy link
Author

orbisvicis commented May 18, 2022

I am able to plot the adaptive triangulated mesh of other functions. I'm not sure why the triangulation of this particular adaptive surface is so incorrect while the resulting interpolated grid (interpolated_on_grid) presents a smooth surface. You're looking at the sum of several shifted step approximation functions.

@basnijholt
Copy link
Member

Hi @orbisvicis, thanks for your interest in Adaptive.

I am having trouble understanding the exact problem. Would it be possible to post a full running example that reproduces the problem?

@orbisvicis
Copy link
Author

orbisvicis commented May 18, 2022

The problem was that the domains differ in scale by about 1e4 while interpolated_on_grid creates a scaled interpolator but unscales the result. So I extract my data from adaptive using _data_in_bounds() - perhaps this should be a public method - and scale, triangulate, then unscale the triangulation before plotting. The result is very good; the triangulation before scaling consisted of nearly-parallel lines.

Question. Since the learner builds an incremental triangulation and the default loss function depends on the triangle area, would I get better results if I scale my bounds before running the learner? For example, scaling [0,10] to [0,1] but reverting the scale within the function to learn. Or is this done automatically?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants