How to smooth pyvista.StructuredGrid?

Question

Here is a Hopf torus made in Python with PyVista:

import numpy as np
import pyvista as pv

A = 0.44
n = 3
def Gamma(t):
    alpha = np.pi/2 - (np.pi/2-A)*np.cos(n*t)
    beta = t + A*np.sin(2*n*t)
    return np.array([
      np.sin(alpha) * np.cos(beta),
      np.sin(alpha) * np.sin(beta),
      np.cos(alpha)
    ])

def HopfInverse(p, phi):
    return np.array([
      (1+p[2])*np.cos(phi),
      p[0]*np.sin(phi) - p[1]*np.cos(phi), 
      p[0]*np.cos(phi) + p[1]*np.sin(phi),
      (1+p[2])*np.sin(phi)
    ]) / np.sqrt(2*(1+p[2]))

def Stereo(q):
    return 2*q[0:3] / (1-q[3])

def F(t, phi):
    return Stereo(HopfInverse(Gamma(t), phi))

angle = np.linspace(0, 2 * np.pi, 300)
theta, phi = np.meshgrid(angle, angle)
x, y, z = F(theta, phi)

# Display the mesh
grid = pv.StructuredGrid(x, y, z)
grid.plot(smooth_shading=True)

The color is not entirely smooth: on the lobe at the bottom right, you can see a line which separates pale gray and dark gray. How to get rid of this line?

Answer 1

I think what's going on here is that there's no connectivity information where the two ends of your structured grid meet. One way to fix this is to turn your grid into a PolyData using the extract_geometry() method, and then using clean with a larger tolerance. This will force pyvista to realise that there's a seam in the mesh where points are doubled, causing the points to be merged and the seam closed:

import numpy as np
import pyvista as pv

A = 0.44
n = 3
def Gamma(t):
    alpha = np.pi/2 - (np.pi/2-A)*np.cos(n*t)
    beta = t + A*np.sin(2*n*t)
    return np.array([
      np.sin(alpha) * np.cos(beta),
      np.sin(alpha) * np.sin(beta),
      np.cos(alpha)
    ])

def HopfInverse(p, phi):
    return np.array([
      (1+p[2])*np.cos(phi),
      p[0]*np.sin(phi) - p[1]*np.cos(phi), 
      p[0]*np.cos(phi) + p[1]*np.sin(phi),
      (1+p[2])*np.sin(phi)
    ]) / np.sqrt(2*(1+p[2]))

def Stereo(q):
    return 2*q[0:3] / (1-q[3])

def F(t, phi):
    return Stereo(HopfInverse(Gamma(t), phi))

angle = np.linspace(0, 2 * np.pi, 300)
theta, phi = np.meshgrid(angle, angle)
x, y, z = F(theta, phi)

# Display the mesh, show seam
grid = pv.StructuredGrid(x, y, z)
grid.plot(smooth_shading=True)

# convert to PolyData and clean to remove the seam
cleaned_poly = grid.extract_geometry().clean(tolerance=1e-6)
cleaned_poly.plot(smooth_shading=True)

Your mileage for the tolerance parameter may vary.

Just as a piece of trivia, we can visualize the original seam by extracting the feature edges of your original grid:

grid.extract_feature_edges().plot()

These curves correspond to the open edges in your original grid:

>>> grid.extract_surface().n_open_edges
1196

Since your surface is closed and watertight, it should have 0 open edges:

>>> cleaned_poly.n_open_edges
0