An Open Source Framework for Integrated Additive Manufacturing and Level-Set-Based Topology Optimization

Topology optimization can be considered as seeking the optimal distribution of one or more materials in a prescribed design domain according to a performance function. Since the pioneering work of Bendsøe and Kikuchi in 1988 [1], several methods of topology optimization have been developed, including homogenization [1–3], evolutionary method [4], solid isotropic material with penalization method [5,6], and level-set methods [7–10].

Even though the eye catching advantage of the level-set methods is that no postprocessing work is needed, however such method is not 100% accurate. Since dealing with computer-driven analysis, which includes optimization algorithms combined with finite element analysis, computational time and resources have to be accounted. That leads to common techniques such as downgrading three-dimensional (3D) problems to two-dimensional (2D) and applying symmetries where possible. Metamaterial examples are also worth mentioning, since they are a combination of periodically assembled unit cells. In these cases, some postprocessing work is practically needed: adding thickness, mirroring, repeating the design periodically, or even a combination of these. Scaling and rotating the design are also useful tools for achieving the desired product.

Additive manufacturing (AM) such as stereolithography, fused deposition modeling, and selective laser sintering can fabricate parts directly from computer-aided design models without specific tooling and fixtures. The main distinguishing feature of AM technologies is their layerwise fabrication approach. The ability of selectively place (multi) materials in 3D spaces provides unique design opportunities and capabilities to create any complex geometric shape and compositional material distribution. Conventional design practice relies heavily on designers' intuition and restricted by the manufacturing technologies. Topology optimization is a powerful computational tool which recasts the design problem as an optimal material distribution problem so that the system will find an efficient geometry which fulfills the requirements quantified by the objective and constraint functions in the optimization scheme. However, topology optimization results exist in complicated freeform geometries, which pose great challenges for conventional manufacturing technologies. The direct rationale implies that the topology optimization would be one of the best technologies for design for additive manufacturing.

Existing work, including Liu and Tovar [11] and Zegard and Paulino [12], have already developed tools for generating STereoLithography (STL) files based on topology optimization results. On the other hand, there are plenty of software packages providing options for postprocessing the design before printing.

However, this paper provides a code that can be integrated with the optimization algorithms, offering the designer a variety of postprocessing tools and generating the STL file. This integrated way reduces or eliminates the need for use of external tools. The code is mainly developed for the needs of a level-set-based topology optimization method, but it is also applicable to other methods, such as solid isotropic material with penalization. In the rest of the paper, an overview of the STL format is provided, the main parts of the code are explained, and several demonstration examples are presented.

STereoLithography is a standard file widely used by commercial 3D printers and computer-aided design software. It contains information only about the boundary surface of the 3D object in a specific way which will be explained in this section.

Since the surface of the design can describe only the boundaries of the object, the normal vectors play an important role on defining on which side of the facet there will be material and which not. The normal vector should always point outward to the void. Hence, on the other direction, the existence of material is implied.

A normal vector along with three vertices contains all the information for each facet with 12 numbers. These data can be used in two different ways in order to generate an STL file: ASCII and binary mode [14]. The main advantage of the ASCII mode generated STL files is that the data can be readable by the designer, and a typical example is shown as follows:

• facet normal 0.257496774 −0.25384599 0.932339907

•  outer loop

•  vertex 100.810745 51.1502495 16.9487667

•  vertex 100.845886 50.9629555 16.8880672

•  vertex 99.8509369 52.9121056 17.6935463

•  endloop

• endfacet

But the size of these files is large, making their use sometimes impractical. On the other hand, if the binary mode is selected, the size of the file is reduced, but the data can only be read by software—a drawback that is not considered important. Moreover, since the designs obtained by topology optimization can be complicated, a high resolution may be needed. Following the binary mode, the size of the generated files is kept relatively small.

At the beginning of the STL file, the user can insert an 80-character title (80 bytes), which should be followed by the total number of the facets (n). The file ends with the prepared data in a 25 × n structure. Even though the dimensions in the STL file are dimensionless, the metric system (mm) will be used here for reference.

In this section, the main parts of the matlab code (Appendix) are explained. The code, apart from creating the STL file, also provides some useful tools for preparing the design before printing. There are two main categories of inputs: the ones that are given while calling the function (direct inputs) and those that can be changed inside the code (built-in parameters). The numbers of the lines that each tool is referring to are given in brackets.

The direct inputs are split into two categories.

A 3D Φ consists of a discretized 2D domain and a value on each node, and hence, the input phi will be a 2D matrix (or 3D matrix for a 3D design). For a 3D matrix, the initial grid will be i × j × k. In a general case, if one dimension has been discretized by n, the grid will consist of n + 1 nodes. Thus, the original resolution will be calculated from the largest dimension as 1/(max([i, j, k]) − 1). Moreover, the code, by default, will scale the largest dimension to 1 mm—a dimension that can be changed later.

In contrast with phi, the following inputs are not vital. Even if they are not defined, the default settings will be enabled. In case of an input phi ready to be printed without any postprocessing, the designer still has to change the scaling factor r in order to get the desired size. That will drive the largest dimension of the initial phi to be r mm, and the rest of the dimensions will be changed accordingly.

Parameter res can change the accuracy of the design [320–327]. Larger res means higher accuracy and by result larger-in size-STL file. For speed improvement of the code, the resolution change happens before [115–122] or after [171–176] the modifications of Φ, based on the value of res (less or greater than 1).

Throughout the code, several figures are generated to visualize the design before and after the changes. It can help the designer decide the changes that need to be done, but it increases significantly the computational time. When a high resolution is preferred, it is recommended that the visualization is turned off (vis = 0). A resolution should be selected wisely: low resolution may reduce the quality of the design (Table 2) and high resolution will severely slow down the code and may lead to huge STL files (Table 3), especially for designs with a lot of complicated features. An ideal resolution has to be found for every design separately, in order to balance the size of the STL file, the speed of the code, and the quality of the design.

Unlike 3D designs that may be ready-to-print, 2D designs need, at least, the necessary step of giving thickness [189–192]. Since the code changes the scale of the desired twice: automatically to 1 mm and by designer's input to desired r mm, the thickness should be provided in terms of a reference dimension. In this implementation, the initial x dimension has been selected. If this dimension was 0.8 mm and t = 20%, the result will be a thickness of 0.16 mm. If r = 100, then x = 80 mm and the thickness will be resized to 16 mm. On the other hand, periodicity or other processes, which have been applied to x direction, will not affect the thickness. In the latter example, if per_x = 2, then x = 160 mm but the thickness will remain 16 mm. The thickness is used to convert a 2D design to 3D, which actually can be described as 2.5D since it only extrudes the initial 2D design. This conversion is done by the convert_2d_to_3d function [297–318] which is included at the end of the code.

In summary, Table 4 presents the five different modes, in which the input can be given, with the corresponding default settings.

Apart from the direct inputs, there are several parameters, each of which is responsible for a specific action that is needed for some designs (Table 5). Since they are not as necessary as the basic inputs, they can only be changed inside the code. By default, these parameters are set to values, keeping these tools disabled. The designer can always modify the code to alter these parameters and can further include them as direct inputs. That will be useful in case of optimizations that need symmetry, periodicity, or change on their placement orientation. By including them as direct inputs, the communication between the code and the optimization algorithm becomes automated. The reader is also referred to Figs. 2–4 for the technical implementation, and Figs. 5–14 of the examples for a demonstration of these tools.

Flipping the design can be useful and is usually selected in order to correctly mirror the design (symmetry). It can also change the desired orientation along with the orientation parameters (which will be discussed later). fl_x, fl_y, and fl_z can be either 0, keeping the design unaffected, or 1 (flip in specific direction) (Figs. 9(a) and 9(b)).

If the optimization result was obtained by having considered symmetries, symx, symy, and symz can mirror the design in the corresponding directions. In 2D designs, by default, the initial two dimensions are set to be x and y; hence, symmetry in z direction has no physical meaning and is being ignored. Value 1 will keep the initial design, and 2 will apply the symmetry in the specific direction (Figs. 9(b) and 9(c)).

Some designs, including metamaterials, consist of periodically assembled unit cells. If a unit cell is generated by the optimization procedure, the result should get repeated accordingly, in order to build the desired structure. For computational time economy, providing thickness to a 2D design comes after applying the periodicity. So, as in the symmetry, per_z has no effect on a 2D design. Instead, increasing the thickness parameter t, from the direct inputs, is an option. If all per_x, per_y, and per_z are equal to 1, no periodicity is applied in any direction (Figs. 9(a) and 9(b)).

All the changes, up to this point, will be saved as a new data file containing the updated Φ (phi.mat) [269–272].

This last tool can be used if the orientation of the printed product has to be altered. Given the variety of 3D printers in the market, the capabilities vary. More economic 3D printers or high-tech large-scale production printer may require the product to be placed in a specific orientation for better result or manufacturing convenience. Parameters x_new, y_new, and z_new are, by default, set to 1, 2, and 3 accordingly. Setting x_new = 2 will move the initial y direction to x direction and so on. As a result, the values 1, 2, and 3 have to be used once, in order to get a correctly rotated 3D design. If more rotation is required, these parameters can be combined with the flip parameters (Figs. 9(c) and 9(d)).

Writing the STL file has to follow a specific way. At the beginning, the user can insert an 80-character title, which should be followed by the total number n of the facets. Then, the prepared data can be written in a 25×n structure.

The code offers a framework of exporting an STL file from a specified Φ. In topology optimization, multimaterial designs are of great interest. They contain more than one level-set functions, depending on the method followed: color level-set method, piecewise constant level-set method, or reconciled level-set method. If the reconciled method [15–17] is followed, the optimization will include one level-set function for each material. In case of a design with m materials, the current code should be utilized m times for Φ₁ to Φ_m. If other methods are followed, the code can still fit in the needs of each framework. Each STL file will represent a single material and be a part of the full multimaterial assembly. Having more than one level-set functions Φ implies that each STL file should be saved under a separate name [258].

While the designer is trying to tune up the parameters for getting the desired output, the visualization should be on and the resolution relatively low, which will accelerate the code. Designs with high resolution may lead to calculations difficult to be handled by the computer. If high resolution is required, a wise selection of res is recommended in accordance to the computer's capabilities. Unlike res, scale factor r does not affect the computational time.

The main errors that can occur are the following:

• Phi structure input is not the expected or a parameter has an illegal value.

• The computer is unable to handle complicated postprocessing work on phi and simultaneously high resolution res. This can make the STL generation time-consuming and may use the entire computer's power.

The code completes the level-set-based topology optimization and creates a design framework integrated with additive manufacturing. In this section, the code is demonstrated in action with several optimization examples. The first one is a benchmark optimization known as Michell beam. The design is obtained through an in-house level-set-based topology optimization algorithm. The second one is a microgripper; a design which is the result of a robust topology optimization process [18]. In these examples, the functionality of the code is presented in terms of postprocessing operations and generation of STL. Both examples need data that can be downloaded online in order to function as input. For convenient testing, the code contains a third example with a simple 3D Φ representing a random design as well. At the end of this section, a metamaterial with negative Poisson's ratio is also presented to demonstrate the code's capability of handling more than one material.

For the Michell beam, the design domain consists of an 11 × 8 plate with a circle filled with material. It is loaded at the center of the right edge with F = 30 N in the vertical direction [19]. The mean compliance problem is solved using level-set-based topology optimization method, and the initial and final design can be seen in Fig. 5.

Since an integrated framework is used, at the end of the optimization the resulted phi is send directly to phi2stl(phi), and the STL file is generated. By using the basic input only (phi), the design will follow the default settings, transforming the 11 mm dimension to 1 mm and using 0.22 mm as thickness (20%) (Fig. 6). No extra tools will be enabled.

A more complete use of the code is to include more direct inputs phi2stl(phi, 110, 2, 0, 10). That will produce a final design with the x dimension scaled to 110 mm, the resolution doubled, the visualization disabled, and a thickness of 11 mm (10%) (Fig. 7).

The STL generated by the code is, then, ready to be printed. Assuming that the part has to be printed just for realizing the design-result, an economic and fast way can be selected. In this example, the object was built with polylactic acid material using a Flashforge Dreamer 3D printer, with low accuracy settings for fast printing. The design procedure, from the initial Φ to the printed model, is shown in Fig. 8.

In this example, some parameters from inside the code will be changed, as well, for a better demonstration of the capability of the code to adapt to the problem, designer's intent, and manufacturing needs. The example contains a half microgripper, which is the result of an optimization [18]. For computational reasons, symmetry conditions were considered during the problem setting.

Since phi contains information only for half of the desired design, the data have to be modified before printing. By calling phi2stl(phi), where phi is the only input, the result would be a 1 mm long half microgripper (Fig. 9(a)). In order to get an STL file with the desired result, direct inputs and built-in parameters have to be set. In this case, the code is called as phi2stl(phi, 50, 1, 1) in order to request a half design that is 50 mm long. At the same time, symmetry has to be enabled in z direction, after the initial design has been flipped. To do so, built-in parameters fl_z = 1 (Fig. 9(b)) and symz = 2 (Fig. 9(c)) have to be set.

For manufacturing purposes, the design might be expected to have a different orientation, such as having the original xz plane attached to the ground. Again, built-in parameters have to be set to: x_new = 1, y_new = 3, and z_new = 2 (Fig. 9(d)). The modifications of the design along with the creation of the STL file are shown in Fig. 9.

The design was further printed using a high-quality 3D printer (Objet260 Connex, Stratasys) with a soft material (TangoBlackPlus) (Fig. 10).

• Boundary [329–340]

  The generation of the STL can be done only on 3D designs that do not have any openings on the external surface. It is not uncommon that some optimization results have open faces at the locations where the original design meets the boundaries. That is why a boundary function is used at the end of all the modifications [194–196]. That ensures that no open faces exist and the creation of the STL file will not fail (Fig. 11).

• Limit [206–210]

  By setting a small offset (±0.0001) on Φ instead of the zero level, it ensures that the boundaries will be included in the final design. In the current implementation, positive values of Φ correspond to material and negative values to void. If the optimization follows the opposite way (positive values of Φ correspond to void), the signs have to be reversed for keeping the normal vectors in the correct directions, pointing outward (phi_stl = ±phi) (Fig. 11).

• Smoothing Φ [212–213]

  Smoothing Φ is a technique of getting good quality surface with a relatively small resolution. That may lead to differences between the original design and the printed object, and therefore, it has to be used with caution (Fig. 12).

Since the previous examples require downloading data, one more example is presented using a simple Φ created inside matlab. Assuming that the current design is a unit cell, it has to be repeated periodically. Here, it will be repeated three times in the x direction and twice in y direction resulting to a 2 × 3 structure (Fig. 13). As in the previous example, apart from having phi as an input, built-in parameters per_x and per_z have to be set to 3 and 2, accordingly. Since, having the design repeated is not common, these parameters are not set as direct inputs. It is up to designer's option to modify the code, in order to make them direct inputs, if that is necessary.

In contrast to the single material designs, a multimaterial design, consisting of m materials, needs m different STL files to be printed. Each STL file will be assigned with a material during the printing procedure. In this last subsection, a multimaterial design is used to demonstrate the use of the code with more than one material. The current example is a negative Poisson's ratio metastructure. A 2D reconciled level-set-based topology optimization method is followed with two functions Φ₁ and Φ₂ associated with one hard and one soft material. The result of the optimization is a 2D unit cell which has to be periodically assembled to get the desired structure. In this example, the code has to be called twice: once for each Φ. The structure is periodic (3 × 3), and the dimension of each unit cell is chosen to be 10 mm × 10 mm and the thickness 100% (10 mm).

The design was printed using a high-quality 3D printer (Objet260 Connex, Stratasys) with a hard (VeroWhite) and a soft material (TangoPlus) (Fig. 14).

This paper presented a method of connecting the topology optimization with the additive manufacturing. By using the existing code, an STL file is getting created based on the binary mode. That makes the STL file smaller-in size-than following the alternative way of ASCII mode. The gap between these two fields was bridged by a matlab code, which apart from generating an STL file, provides also some useful postprocessing tools: scaling, changing resolution, adding thickness to 2D designs, flipping, symmetry, periodicity, and changing the orientation. It can work with the existing optimization algorithms based on level-set method in a straightforward way, but it can also be used with other topology optimization methods. Inputting direct inputs for the basic use of the code and changing built-in parameters were explained in detail. The only necessary input for generating the STL was the level-set function Φ that corresponds to a 2D or 3D design, and the rest of the inputs are used to create the desired design. The built-in parameters are generally less frequently utilized, but the designer still has the freedom to change them or even include them as direct inputs. The speed aspect of the code has been analyzed in terms of the design's resolution and visualization of the results. The paper concluded with the demonstration of three examples-designs that have been obtained using a level-set-based topology optimization method.

The purpose of this work was to create a tool that can be integrated with the existing topology optimization algorithms. However, a graphical user interface is also available (Fig. 15) for postprocessing already generated designs.

The authors acknowledge the support from the Region 2 University Transportation Research Center (UTRC).

• Ford Motor Company (Integrated Topology).

• National Science Foundation (Grant No. CMMI1462270).

• State University of New York at Stony Brook (Start-Up Funds).

The supplementary data code used in this work to covert the level-set-based topology optimization results into STL files can be found in the website link.²