Visualization¶
Triangle Mesh¶
[1]:
import plotly.io as pio
from lapy import Solver, TetMesh, TriaMesh, io, plot
pio.renderers.default = "sphinx_gallery"
This tutorial will show you some of our visualization functionality. For that we load a larger mesh of the cube and compute the first three eigenvalues and eigenvectors. We also show how to save the eigenfunctions to disk.
[2]:
tria = TriaMesh.read_vtk("../data/cubeTria.vtk")
fem = Solver(tria)
evals, evecs = fem.eigs(k=3)
evDict = dict()
evDict['Refine'] = 0
evDict['Degree'] = 1
evDict['Dimension'] = 2
evDict['Elements'] = len(tria.t)
evDict['DoF'] = len(tria.v)
evDict['NumEW'] = 3
evDict['Eigenvalues'] = evals
evDict['Eigenvectors'] = evecs
io.write_ev("../data/cubeTria.ev", evDict)
TriaMesh with regular Laplace-Beltrami
Solver: spsolve (LU decomposition) ...
Let’s look at the result by visualizing the first non-constant eigenfunction on top of the cube mesh. You can see that the extrema localize in two diametrically opposed corners.
[3]:
plot.plot_tria_mesh(tria, vfunc=evecs[:,1], xrange=None, yrange=None, zrange=None, showcaxis=False, caxis=None)
We can also adjust the axes and add a color scale.
[4]:
plot.plot_tria_mesh(tria, vfunc=evecs[:,1], xrange=[-2, 2], yrange=[-2, 2], zrange=[-2, 2], showcaxis=True, caxis=[-0.3, 0.5])
Tetrahedral Mesh¶
Next we load a tetrahedral mesh and again compute the first 3 eigenvectors.
[5]:
tetra = TetMesh.read_vtk("../data/cubeTetra.vtk")
fem = Solver(tetra)
evals, evecs = fem.eigs(k=3)
evDict = dict()
evDict['Refine'] = 0
evDict['Degree'] = 1
evDict['Dimension'] = 2
evDict['Elements'] = len(tetra.t)
evDict['DoF'] = len(tetra.v)
evDict['NumEW'] = 3
evDict['Eigenvalues'] = evals
evDict['Eigenvectors'] = evecs
io.write_ev("../data/cubeTetra.ev", evDict)
--> VTK format ...
--> DONE ( V: 9261 , T: 48000 )
TetMesh with regular Laplace
Solver: spsolve (LU decomposition) ...
The eigenvector defines a function on all vertices, also inside the cube. Here we can see it as a color overlay on the boundary.
[6]:
plot.plot_tet_mesh(tetra, vfunc=evecs[:,1], xrange=None, yrange=None, zrange=None, showcaxis=False, caxis=None)
Flipped 24000 tetrahedra
Found 4800 triangles on boundary.