metaball meshing

For another project I have been working with metaball scalar fields, and I thought it would be fun to have a go at meshing them using Paul Bourke's marching tetrahedra methods. I know that Daniel Piker has done a little of this before more generally, and that there is an excellent implementation of the marching cubes algorithm as a part of Millipede...this is really just more specifically set up for generating metaball meshes.

On its own, the marching tetrahedra algorithm as outlined by Paul Bourke produces some odd results for metaballs...basically, the points that drive the mesh are assigned along the edges of the tetrahedra according to a linear interpretation of the iso value...what this means is that if you have an "ideal" iso value say of 0.29, and the sample points at each end of your tetrahedron edge return values of 0.2 and 0.3, the point will be located along that edge at 10% of the length away from the point with the 0.3 sample. In general, this works fine, except (from what it appears to me at least...could be wrong here) that the metaball field itself isn't really linear, and the tetrahedra used for sampling are carved from cubes such that the lengths of each edge differ significantly...so that there are a few spikes and dips here and there. Sampling resolution has a lot to do with these results as well. I address this issue here by embedding laplacian smoothing as a part of the process, and also by an additional field value sampling run for each mesh vertex after the initial mesh has been constructed.

So the inputs are your metaball centers (as a list), their radii (also as a list, which is the location at which their isosurface value contribution is 0), the iso level threshold, the grid resolution for sampling (here a lower number produces a higher resolution...basically is the rough size of the grid in the x, y and z directions), and the number of smoothing passes run after the mesh is generated. The outputs include the metaball meshes as well as the bounding box that's used, as well as the centers of the sample cells (both of these last are superfluous but perhaps interesting). I'd suggest starting with a low resolution (higher res value) for sketching, as it works pretty quickly...then if you want, you can increase the the sampling density (and incur the resulting exponential slowdown).

I've used a couple of tricks for reducing the number of sampled points (nothing so sophisticated as an octree subdivision), but there are also definitely inefficiencies in points shared between sampling cells having their fields calculated multiple times. There's really plenty of room for improvement...this is just a rough pass really, but I figured some people might have some fun with it and it seems that every now and again the issue of making metaball meshes comes up, and I figured that it'd be useful to have a quick way of building (relatively) clean meshes laying about.

Also, probably lots of bugs everywhere...feel free to point them out...and the code's pretty ugly, but poke around, take it, do what you will with it.

Comment

You need to be a member of Grasshopper to add comments!

Comment by Taimaa Barakat on October 17, 2020 at 1:57pm: Hello, I just wanted to thank you, this helped a couple of classmates and I to figure out how to model connected circles.

Comment by David Stasiuk on September 18, 2015 at 7:31am: You can also check out Cytoskeleton (it's within the Exoskeleton component group: http://www.grasshopper3d.com/group/exoskeleton) It has a good way to create a thickened mesh from the wireframe of another mesh, and you can have it work on the Dual of the other mesh as well...

Comment by David Stasiuk on September 18, 2015 at 6:46am: You want to get a "Mesh Dual" Weaverbird has it...

Comment by Victor Lu on September 18, 2015 at 6:11am: Thank you David! One more question, how could you get your mesh in hexagon. Mine is always triangular.

Comment by David Stasiuk on September 15, 2015 at 9:41am: Hi Victor. You may want to check out the Cocoon plug-in I've developed, which advances this metaball meshing strategy quite a bit further.

But as to the patterning, this project by TheVeryMany appears to use some sort of structural analysis on the shell to erode portions of the mesh from the final structure (it's actually very difficult to know with their projects, as they don't disseminate their work beyond images/built objects...no papers on their methods are out there). You could use a plug-in like Millipede to get an understanding of the shell's structural properties.

But for sure, it's not a trivial thing to make this type of structure...

Comment by Victor Lu on September 13, 2015 at 10:23pm: Hi, David, I tried to use your fantastic definition to reverse engineering this project.

But it seems that it is impossible to apply pattern on this complicated mesh, especially on the edge two spheres join together. Could you please do me a favor to achieve it? Thanks!

Comment by David Stasiuk on February 21, 2015 at 6:34am: Hi Rodion- sorry I missed your post on Wednesday! I would really stay away from MeshFromPoints...it's super heavy and totally unnecessary, as the marching tetrahedrons algorithm organizes all of your points perfectly for creating your own mesh faces to begin with. You really don't need to build surfaces at all, or get points...

I don't use Python, so I can't describe the syntax there, but the C# syntax for creating a mesh from those points would be something like this:

First, create your final mesh object:

Mesh MyFullMesh = new Mesh();

Then, run your marching tetrahedrons algorithm, and for each Triangle output (which should give you three points) do this:

Mesh NewTriangle = new Mesh();

NewTriangle.Vertices.Add(TrianglePoint0);

NewTriangle.Vertices.Add(TrianglePoint1);

NewTriangle.Vertices.Add(TrianglePoint2);

NewTriangle.Faces.AddFace(0,1,2);

MyFullMesh.Append(NewTriangle);

MyFullMesh.Vertices.CombineIdentical(true, true);

Sorry that it's late...hope it helps!

Comment by David Stasiuk on February 6, 2015 at 1:37am: Hi Rodion-

I think you're confusing marching cubes with marching tetrahedrons...you're right that marching cubes has 256 cases, but the lookup you're using for marching tetrahedrons is correct. There are just the 8 cases, so you're actually doing it right here (scroll down on this page, and you'll see a separate subset all about marching tetrahedrons http://paulbourke.net/geometry/polygonise/). The benefit to using marching tetrahedrons is exactly this: that the number of possible "cuts" through the tetrahedron are dramatically smaller in number than those through a cube.

However, I have found that also what you're seeing that the linear interpolation creates some odd distortions (which is why I went ahead and later did the marching cubes implementation). Some of this comes from the density of the sampling grid: the more dense, the fewer distortions.

What I would suggest, if you want a (relatively) quick way to improve this outcome:

1) build up a full mesh rather that bunch of surfaces, and use Rhinocommon to combine identical vertices, and rebuild the vertex normals

2) run a couple rounds of laplacian smoothing on the mesh to better distribute your vertices (for each vertex, make it equal in location to the average of its neighbours)

3) create a line normal to each vertex roughly the length of your sampling grid and test the endpoints of it against your scalar field formula, and then do one final linear interpolation between those two points for your vertex.

This should give you a smoother mesh for sure.

But good work getting this far!

Comment by David Stasiuk on February 4, 2015 at 4:43am

Sounds like you're on the right track, Rodion. I really haven't dug into this particular code very closely, but from what it appears to me, yes, the PolygoniseTri should give you the vertex locations for any new mesh faces derived from that particular sampling tetrahedron. And it looks to me like it's correctly cycling through each tetrahedron in each cube:

triangles.extend(PolygoniseTri(grid,isolevel,0,2,3,7))
triangles.extend(PolygoniseTri(grid,isolevel,0,2,6,7))
triangles.extend(PolygoniseTri(grid,isolevel,0,4,6,7))
triangles.extend(PolygoniseTri(grid,isolevel,0,6,1,2))
triangles.extend(PolygoniseTri(grid,isolevel,0,6,1,4))
triangles.extend(PolygoniseTri(grid,isolevel,5,6,1,4))

So in grasshopper, this is where you would define a new mesh face, using the points indicated from this function (looks like they return either one or two triangles depending on the intersection state of the tetrahedron, which is right). By iterating through your entire sample grid, you'll create your continuous mesh.

I have found that the easiest way to do it is to create a new Mesh for each sample tetrahedron, and then to append them to a complete mesh as you go along. Hope this helps...it'll be cool to make a Python implementation for GH!

Comment by David Stasiuk on February 2, 2015 at 1:44am: This looks like it is a pretty reasonable Python implementation, based on the same source that I used for mine:

http://michelanders.blogspot.dk/2012/02/marching-tetrahedrons-in-py...

You'd have to adapt it to Rhinocommon, of course, but all of the key ingredients are there.