PAVEL: Decorative Patterns with Packed Volumetric Elements

Packing 3D shapes in an arbitrary volume is a well-studied problem, and it is known to be NP-complete. It is relevant in graphics, engineering, operational research, and many other fields. Recent works such as (Egeblad et al., 2009; Liu et al., 2015; Romanova et al., 2018) achieve an excellent volume coverage. Still, they are heavily limited in the orientation of objects, are significantly time-consuming, and can only deal with a limited number of items with few facets. Overcoming these limitations, Ma and colleagues (2018) propose an algorithmic approach well suited to complex models typically used in graphics. Their algorithm initially scales down all the objects to avoid intersections, then iteratively increases their size, optimizes rotations, and tries some object swapping to arrange local solutions. The algorithm can obtain good results but can take up to an hour to finish the computation to place less than 100 objects. Our approach differs from both since we pack elements on surfaces and since we decide element sizes upfront, which makes the problem generally harder.

The packing of 2D shapes is a relevant field of research for visualization and illustration purposes. Reinert and colleagues (2013) proposes an algorithm based on Centroidal Voronoi Diagram relaxation to fill an arbitrary-shaped planar domain with multiple shapes. Saputra and colleagues (2019) also describe an algorithm to pack shapes in a planar domain, with a greater focus on leaving as little space as possible. To achieve their results, they allow deformations of the decorative shapes and adds small ones, in the end, to obtain better coverage and reduce space. For similar purposes, Kwan and colleagues (2016) propose a new shape descriptor and use it to measure the touching shapes’ quality. Another relevant method, by Chen and colleagues (2017), packs tiles on top of a curved surface for fabrication and decorative purposes. They describe an attract-and-repulse algorithm in which tiles are allowed to move and rotate, but not to deform, which can take up to fifteen minutes to place more or less 150 tiles. Two works focusing on decorations are (Lu et al., 2015), and (Chen et al., 2016), which investigate the possibility of transforming a surface carving holes following a pattern. In (Lu et al., 2015) a texture guides the building of a pattern of connected stripes taking into account the structural problems that can arise when fabricating the final object. In (Chen et al., 2016), the pattern generates filigrees using the skeletal representations of simple starting elements. The described pipeline can then place them on a surface, greedily refine their position, and slightly modify them to ensure contact between close elements without overlaps. However, all these methods are inherently different from ours since they treat 2D shapes rather than 3D elements. Furthermore, most of these works employ specific representations that cannot be generalized to 3D.

Placing decorations with circular symmetry over the surface of an object can be interpreted as a packing problem as well as a point sampling one. Blue noise point distributions and their variations can undoubtedly be used to roughly pack a surface with circular shapes. Bridson (2007) proposes a simple and efficient algorithm to generate a blue noise distribution. Due to its greedy nature, it tends to leave significant gaps, particularly undesirable for artistic purposes. A more uniform distribution is given by a Centroidal Voronoi Diagram, which can be computed over arbitrary surfaces as described in (Xin et al., 2016). Knöppel and colleagues (2015) introduce a special kind of entity distribution over a surface putting regular striped patterns that can adapt to the local characteristics of the shape, curvature included. This method is useful guidance for scattering points on the surface of the form only onto the stripes. We make use of the latter two methods to place seeds for decorative elements in our pipeline.

Algorithms to control deformable surfaces based on level set formulation are popular in Visual Computing after the seminal works by Osher, Sethian, and Malladi (1988; 1995). Using implicit formulations, it is possible to model complex interface behaviors, including elastic and plastic deformations, but with huge computational complexity and stability issues (Barton et al., 2011). A simplified approach is represented by the fast marching algorithm (Sethian, 1996), which allows a fast estimation of the evolution of closed surfaces as a function of an underlying velocity field assuming unidirectional front propagation along the surface normals. By controlling the underlying speed, it is possible to drive the surface evolution towards the desired position. Simplified implicit deformable models have been used in Computer Graphics to simulate soft bodies’ contact, and the resulting deformation (Gascuel, 1993) and to animate deformable shapes (Desbrun and Cani, 1998). A similar approach has also been recently exploited in a generalized framework to interactively control the blending of two shapes represented by signed distance fields (SDF) on regular grids (Angles et al., 2017). We adapt the fast marching approach to creating the effect of plastic-like volume-preserving deformation of multiple objects efficiently in our work.

The most common approach to extract a mesh representation from a discrete scalar field is the Marching Cubes algorithm (Lorensen and Cline, 1987). When the scalar field is Boolean, an additional smoothing step is often necessary. Coeurjolly and colleagues (2018) propose an optimization-based approach to extract a smooth manifold from a voxel representation to avoid this. This method works directly on the voxel field’s boundary quads, moving the vertices to optimize fidelity, fairness, and smoothness. We do not use this method in our reconstruction since it is unacceptably slow for resolving our voxel grids. It produces overly smoothed results, which removes details in small decorative elements. We instead use marching cubes and a smoothing remesher to convert the deformed decorations back to a surface.

To obtain the deformed elements, we cannot use existing physics engines for several reasons. First, most engines cannot accept an initial condition with large overlaps between objects, which in our case is as high as 30% in volume. Second, complex custom forces should be defined to ensure a reasonable equilibrium, e.g., the forces that sculptors apply to get the precise deformation they desire while making sure that decorations adhere to each other. We tested multiple physics solvers, including Houdini’s FEM and XPBD solvers (SideFX, 2020), Blender’s solver (Blender, 2020), and a projective dynamics solvers (Fratarcangeli et al., 2016). In all cases, we could not successfully remove the overlap for the reasons listed above. Figure 4 shows a failure case. While physics engines’ behavior is reasonable, as we start the simulation in a non-physically-plausible configuration, we cannot rely on these tools for our goals.

Since we cannot use physical based simulators, we propose a novel method based on an implicit formulation, inspired by works on interactive visual simulation of volumetric object contact (Gascuel, 1993). We rely on fast marching (Sethian, 1996) using a specifically designed parametric function for the propagation field that determines the deformed shape appearance and enforce volume preservation.

We voxelize each decoration $d_{i}$ in its own dense grid $g_{i}$ that is $50\%$ larger than the decoration’s dimension to ensure that the computed deformation will fit in the grid. Compared to the whole surface, this is a sparse discretization that keeps memory consumption relatively low since we are effectively discretizing a slab around the surface. While each decoration has its own representation, we ensure that the voxel dimension and the grid rotation is consistent across all voxelizations. As a consequence, we can perform Boolean operations just by keeping into account their mutual offsets. Additionally, separate voxelizations, instead of one global labeled grid, allow for simple handling of the overlaps and compatibility with the fast marching algorithm.

Given as input a set of voxelized decorations $D=\{d_{0},d_{1},\dotsc,d_{n}\}$, our first step is to resolve the overlaps. We do this by detecting the voxels shared by multiple decorations and removing them from all but the one with the shortest distance from the decorative element center. Furthermore, we remove all the voxels of the decorations shared with the base objects as well. These removals generate a new set of decorations $D^{\prime}$, where each decoration has associated a measured volume loss $v_{i}$, equal to the sum of the voxels removed in the previous steps. We then ensure that each element recovers its original volume by expanding it in a way that simulates a plastic deformation.

To compute the expansion, we use a fast marching approach (Sethian, 1996) to grow each element in the empty regions until they recover the original volume. To simulate a plastic deformation, which typically generates bulges near the soft bodies’ contact points, we use a non-constant velocity field to drive the fast marching algorithm. For each shape $d^{\prime}_{i}$, we find the contact surfaces with the other shapes, and from them, we build the respective distance fields. The distance fields are then converted in velocity fields $F_{i}$ through a parametric formula, the volume recovery function, that can be adapted to obtain the desired behavior. The recovery function is a piece-wise polynomial, obtained by connecting three second-degree polynomials and enforcing $C^{1}$ continuity, and defined by:

In the inset, we show an example plot of the function. We keep the velocity null in contact points to simulate the constraints related to the attachment of surfaces. We then have two control points $a,b$ used to define the distance for maximal bulging and the maximal distance reached by the plastic wave. This allows effective control of “material behavior” driven by the amount and locality of the

deformation.

The velocity fields are then weighted proportionally to the volume losses to simulate a faster growth of the elements more “squeezed” by the neighboring ones and used to compute the new arrival times for each decoration. To recover the volume lost by each decoration, we evolve the front of each decoration with the fast marching method. For each elements’ sub-grid $g_{i}$, a set of active front voxels is initialized, taking those adjacent to the initial voxelization border. In contrast, we set all voxels belonging to either the initial shape, the base shape, or other elements to inactive. Each decoration can recover its volume only by occupying empty voxels.

The algorithm then iterates over all the elements. At each step, it evaluates the arrival time for all the voxels of the active front. It selects the voxel with a minimal arrival time according to the Eikonal equation’s solution estimated on the base of the velocity field defined in the corresponding grid.

The selected voxel is then assigned to the current element, and it is set as inactive for the other elements’ front propagation. The active front is finally updated, removing the selected voxel from it and adding to it the voxel’s free neighbors. This procedure stops when all the decorations have recovered their original volume. It is theoretically possible that the algorithm could not find enough voxels to recover the volume when the evolution of neighboring elements occludes the front, but this is not likely to happen for a layer of small decorations over a base surface and never happened in our experiments.

The method ensures by construction that all volume is correctly maintained while deforming the decorations. The algorithm is fast since it greedily iterates over elements instead of estimating all the fronts’ arrival times at once. Fast marching is not easily parallelizable to speed up the computation. Still, as we simplify the decorations’ distance fields as independent of each other, we can concurrently compute them for a significant performance gain. We use the fast marching of (Scikit FMM, 2020) and handle IO and voxelization with (Trimesh, 2020).

The method described above can simulate the elements’ plastic deformation if the decorative elements are initialized with a reasonable overlap. To understand what is a “reasonable” overlap, we performed several tests to define control parameters for the point distribution. To control the decorations’ overlap, we choose to use the surface coverage, that is, the ratio between the sum of the surface areas of the decorations’ footprints and the base object area to be decorated. Fig. 5 show results of packing small spheres on a base sphere with different surface coverage values. As we can see, a ratio of 1.2 is barely enough to cause some deformation, while a ratio of 1.6 produces excessive and unrealistic distortions of the elements’ shape.

We also performed extensive tests to choose the parameters of the volume recovery function (Eqn. 4). Fig. 6 shows the results obtained with the same overlapping spheres but different choices for the values of the $a$ and $b$ parameters. By shifting the control parameters, it is possible to simulate different plastic behaviors. All the results shown in the paper, but Fig. 6, are obtained with the same settings $(a=0.1,b=0.3)$.

As it is possible to see, for example, in Fig. 7, this choice provides results that are visually quite close to real plastic deformation.

While manual placement provides complete artistic control, the decoration of complex surfaces with small decorative elements would be very time consuming, just like in real-world crafts. To simplify creating complex decorations, we also support the automatic placement of decorative elements using two procedural methods.

Our placement methods are based on point sampling algorithms. We support the placement of objects with approximate 2D radial symmetry, with the symmetry axis aligned to the base surface’s local normal. We also assume that the decorations are small compared to the base surface’s curvature to ignore the variations in the local volume extruded from the base surface. Both of these are reasonable assumptions found in the handcrafted models that inspired our work.

Instead of a packing algorithm, our use of point sampling is motivated by a trade-off between speed and packing quality. It’s worth reminding that optimal packing is NP-hard in the Euclidean domain and just as hard in the intrinsic manifold metric. This means that all published packing algorithms are approximations that trade-off element shapes, packing quality, and computation time. In general, it is well known in both the literature and the practice that most packing algorithms fail or are not optimal for all but the simplest shapes. Note also that we do not allow to change the decoration size, which is the most used approximation for soft packing in graphics. For all these reasons, our point sampling approaches are just another heuristic within their stated constraints. If better placement heuristics are available, we could integrate them into our pipeline with ease.

Our first automatic placement method produces isotropic distributions that work well for decorations with cylindrical symmetry. The user controls the density of decorations by setting the desired surface coverage, and the algorithm instantiates the correct number of seeds to guarantee the target density.

To ensure a high-quality distribution, we rely on optimization using a method based on the Centroidal Voronoi Tessellation (CVT) and its efficient implementation described in (Yan et al., 2009) and available in the Geogram library (Levy, 2020). We generate the base surface’s tessellation with the desired number of vertexes and use them as the centers of the decorations. When converged, this iterative approach produces a point distribution with good statistical properties, resulting in a high-quality decorated surface. To obtain a more “hand-made” look, we allow users to reduce the number of iterations. Fig. 8 shows an example results obtained with this method that closely resembles the hand-made objects shown in Fig. 2.

The previous method works well if the decorative elements are small compared to the surface curvature. If that assumption is not valid, the decorated surface has visible artifacts, shown in Fig. 9. The main problem is that our method generates equally spaced seeds that represent the bottom of the decorations but does not take into account decorations’ volume and surface curvature. Decorations overlap more in concave regions than convex ones. The volume above the surface is smaller in the former case and can result in heavily compressed areas in concave zones and non-contacting in convex regions. When decorations are small compared to the surface curvature, this effect is negligible. But it becomes more noticeable as the decorations’ size increases.

For this reason, we also provide a method that considers the volume of the elements when point sampling their locations. The basic idea is to extract an isosurface at a signed distance $d$ from the base mesh, perform the sampling on it, and then map back the seed points on the original base surface, projecting them along the local normal direction.

An appropriate choice of $d$ is the height of the element’s z-section with the maximal area. In this way, we can handle large decoration sizes even on surfaces with varying curvature. We show a clear example of this in Fig. 9. Here the bunny features gulfs and bulges. They are handled without problems by sampling the offset surface. Conversely, sampling on the surface reveals significant artifacts. Since decorations are small in general, and the modified sampling has a non-negligible computational cost (see Section 7), we leave the choice of whether to use the offset surface to the user.

To support the automatic generation of anisotropic element patterns, we also provide a stripe-based sampling approach. This method can create gridded patterns of decorations with approximate cylindrical symmetry and even with elongated shapes.

Our method is based on (Knöppel et al., 2015) that allows us to generate patterns of evenly spaced stripes over a generic base surface. The method is based on a simple optimization problem and approximately maintains a user-specified orientation and spacing of the lines by automatically inserting branch points where necessary. Furthermore, the generated patterns are globally continuous, which avoids visible seams after the applications of our decorative elements, and with frequencies not depending on the mesh resolution. We exploit the stripe-based method by generating two crossing stripe patterns with a user-defined crossing angle and sampling seed points at the line crossings. In this manner, we obtain lines with points sampled at constant distances and with a controlled offset among the decorative elements of two adjacent stripes. Fig. 10 shows an example with crossing stripes creating an angle of 60 degrees and elements placed with a constant orientation concerning the main stripe.

Since we can independently set the sample spacing on the main stripes and the distance between adjacent stripes, we are no longer forced to pack objects with approximate radial symmetry, but we can insert elements with arbitrary elongation. Furthermore, since the stripes encode a local directional field, we can use non-isotropic decorative elements oriented in a controlled way with respect to the stripes’ directions. Fig. 11 shows an example with perpendicular stripes-based sampling and non-isotropic elements with main directions aligned with the stripes fields.

In our stripe-based algorithm, we create the sampling points and estimate and store the point connectivity of the sample along with the two stripe directions. In this way, we can develop decorations alternating different elements and elements’ orientations, as shown in Fig. 12.

Aside from placement, we can easily integrate controls for decorations’ rotations and filling different areas with different decorations. The rotation for each element can be set randomly or aligned to a field defined over the object surface, giving the possibility to create different artistic effects that work particularly well when placing elongated objects over a gridded seed distribution. Figure 12 shows an example of an artistic pattern obtained with alternate rotations of 180 degrees of the elements along the principal stripe field direction.

The average elements’ size is chosen to guarantee the desired decorative effect, and the number of elements can be tuned to guarantee the desired amount of initial overlap. However, we leave the possibility of applying a further random perturbation to the size of each decoration to increase the number of design options. We allow variations up to $\pm 10\%$, as larger ones would lead to inconsistent overlaps and too many empty spaces.

Lastly, it is possible to load more than one decoration and randomly select one of those for each seed. We can exploit this to increase the final results’ variety loading a multitude of slightly different decorations (e.g., stones of various types, flowers, etc.).

Fig. 13 shows an example of decoration mimicking the one in Fig. 2 obtained by using three different elements shapes and random rotations with respect to the main field direction.

To generate a 3D mesh of the whole object, we first independently create meshes for each deformed decorative element using the Marching Cubes algorithm (Lorensen and Cline, 1987). Then, we merge the elements and the base object into a single model and remesh the result using the Central Voronoi Tessellation approach (Yan et al., 2009), which provides feature-preserving results and is robust and efficient for large meshes.

It is possible to use our pipeline also to create decorative shells for real-world objects. To decorate an object, we first scan it to obtain an accurate digital model and then create a PAVEL decoration for it. This produces a mesh representing all the decorative elements fused together. We then decompose this final mesh into printable parts and wrap it around the original object. Fig. 14 shows the whole process of decorating an object.

For objects with a simple topology, the decomposition works as follows. First, we segment the base object. Any mesh partitioning can be used for this purpose. In the example shown, we cut the object with two mutually perpendicular planes passing through the model’s symmetry axis. Given the base segmentation, we cluster the deformed elements placed on each patch and remesh them into single shell elements to be 3D separately printed. One can use a more appropriate decomposition for more complex base objects while the rest of the pipeline remains the same.

In the inset, we show the final result obtained by gluing the parts of the tangible shell created from the design in Fig. 14 over the base object. In this case, the object was a simple paper cup covered with modeling clay. We scanned it with a low-cost device (Intel Realsense F200) and automatically decorated with an isotropic distribution of elements with different shapes and randomized size and orientation. We finally printed the shell elements with a Stratasys J826 3D printer.

Using PAVEL, we can create a variety of decorations that mimic real-world artworks. Fig. 15 shows side-by-side views of original artworks, courtesy of Heather Knight, and objects generated with our pipeline and rendered with clay-like material properties.

PAVEL decorations are not limited to the replicas of clay pieces but can be used to generate various artworks of different kinds. Figs. 16 and 17 shows a few examples, demonstrating potential uses in jewelry and pottery, exploiting manual elements positioning and automatic pattern creation.

To understand if PAVEL could be a useful design tool, it is necessary to evaluate the time required by the different steps and the global effort required to create the decorated models. We obtain all our results on a 12-core Ryzen 3900x with 32GB of RAM. Tab. 1 shows the time required by the automatic seeding procedure. All the times for the non-curvature adaptive methods ranges between 0.6 and 1.5s, with minimal differences between CVT and striped seeding. The most computationally expensive sampling strategy is, as expected, the curvature-dependent one. Our test took 54s to extract the isosurface but only 2.5s to sample it. This speed allows a reasonably fast interactive design procedure. The digital artist can quickly test different object choices and seeding options visualizing the seeded configurations superimposed on the base model.

Tab. 2 shows the times required by the volume recovery and mesh generation. While the whole procedure is clearly not interactive, the time required is reasonable for an offline rendering procedure. It is not critically variable varying the surface coverage and the elements’ volumes overlap within the parameters ranges discussed in Section 4. The time required by the volume recovery step is directly dependent on both the number of decorations and their resolution. The times range from 20 to 360 seconds, with the faster times achieved by the models with around 50 decorative elements. Still, the method scales relatively well, handling nearly 2,000 decorations in just under 6 minutes. Another critical factor in performances is the resolution, which is also crucial to maintaining the decoration’s details. Our results show how the method can handle millions of voxels to represent the decorations, even in intersections up to 33% of the volume. Lastly, another time-consuming step is CVT smoothing. While it ranges from 10 to 200 seconds to remesh our test models, it’s worth noticing that the marching cubes meshes’ resolution reaches up to 40 million triangles, and we output smooth meshes with up to 4 million triangles.

Our method has two significant limitations. First, decorations’ automatic placement is limited to approximately isotropic elements or elongated shapes placed on stripes. As stated in Section 1, this is a willing design choice to simplify the packing procedure leaving to the deformation step to fill the gaps between elements. However, we plan to improve the initial packing step with novel options, adding a greedy refinement of elements’ position and orientation after the initial placement to enhance the quality of the packing of non-radially symmetric elements.

The second limitation is our use of fast marching that limits the complexity of the deformed objects’ details since the number of voxels quickly grows w.r.t. surface details. In the future, we plan to explore sparse representations based on spatial hashing and wide-branching trees that are of common use in 3D animation.

In addition to the improvements listed above, we plan to improve the pipeline by adding more options to control the design, namely the possibility of performing automatic seeding for specific decorations on different parts of the mesh.

Another planned work is the optimization of the shell decomposition for printing. We plan to develop methods for the automatic decomposition of the PAVEL decoration in patches of controlled size and printer-friendly. We will also investigate the possibility of adding printability constraints in the fast marching evolution.

We also plan to better exploit the elements’ deformation’s parallelizability to obtain a faster implementation. We believe that it should be possible to parallelize both the computation of the elements deformations and the merging of the final result and obtain significant speed up in the case of large numbers of small decorations.