Interactive infill topology optimisation guided by user drawn patterns

ABSTRACT

Widespread use of topology optimisation as a design tool for additive manufacturing faces major inhibiting obstacles, such as high computational costs and complexity, concern for other failure modes, and manufacturability. Interactive infill topology optimisation presents an alternative approach to circumvent some of these barriers. The novel contribution of the present work prompts the user to draw a tailored infill pattern, specify regions of interest to locate the infill, and control how strictly the pattern is replicated in the material layout of the design using appearance constraints. This approach improves engineering metrics not directly included in the optimisation formulation by incorporating the user’s engineering experience, thereby avoiding increased computational costs, parameter tuning, and numerical artifacts associated with complex objective functions and constraints. Two 2D benchmark examples increase the linear buckling resistance and energy absorption, respectively, and a 2.5D example minimises compliance while reducing the quantity of overhang supports for additive manufacturing.

KEYWORDS:

• Topology optimisation
• iterative user interaction
• infill
• appearance constraint

1. Introduction

Topology optimisation is a computational design method to generate material efficient, high performing structures according to a user defined objective function and constraints. The free form design characteristics of topology optimisation yield unique, unconventional structures that are frequently difficult to manufacture using traditional means; therefore, additive manufacturing is often essential to production. The layer-by-layer approach of most additive manufacturing technologies, such as Fused Deposition Method (FDM) or Selective Laser Melting (SLM), are feasible manufacturing methods for the topology optimised designs. Especially as plastic extrusion and resin photolithography 3D printers are more affordable and widely available, it is less difficult to produce topology optimised structures. The field of topology optimisation early on recognised its reliance upon 3D printing technology and therefore, the research expanded in the last decade to develop algorithms tailored to additive manufacturing [1–8]. Examples include topology optimisation codes with an overhang angle constraint [4–6,9], modified stiffness for filament bonding strength [7,10] and infill optimisation [11–16].

Density-based topology optimisation discretizes the design domain into finite elements and imposes loads and boundary conditions. The elemental densities are the design variables and traditionally, the optimisation minimises the compliance, or maximises the stiffness, subject to a volume constraint. Additionally, a filter is applied for the regularisation of the design variables, imposing a length scale that suppresses checkerboards and ensures the optimisation formulation is well-posed [17]. Initially, the algorithm distributes the material uniformly across the design domain. A finite element analysis determines the primal solution of the density field for a given optimisation iteration, calculating the objective function and constraint values, as well as determining the derivatives of the objective function and constraint with respect to each elemental design variable of the density field. The elemental derivatives of the objective function and constraint, called sensitivities, are applied in a mathematical programming, such as the Method of Moving Asymptotes (MMA) [18], redistributing the density field in each optimisation iteration to minimise the objective function and satisfy the constraints. The next optimisation iteration resumes with the new density field associated with a lower objective function value (higher performance) and feasible constraint values (satisfying the volume constraint). Topology optimisation algorithms iteratively approach a local minimum for the material layout of the design space.

The flexibility when applying gradient-based optimisation allows for a variety of objective functions or constraints. The field of topology optimisation has grown to reduce stress concentrations [19–25], maximise linear buckling load [26–29], and increase energy absorption [30–33] as objective function and constraints. However, these enhanced algorithms are generally associated with higher computational costs as additional optimisation iterations, evolved modelling, prerequisite knowledge of feasible and well-posed topology optimisation formulations, and more design parameters. Interactive topology optimisation for guided geometrical patterns does not attempt to match the performance of specialised optimisation codes, but instead provide a computationally cheap, user-guided alternative approach to improve alternative engineering metrics while minimising compliance or similar simple objective functions. For example, the 250-line topology optimisation code maximising the buckling load [26] takes 420 minutes to increase the buckling load by a factor of 5 [34], while the present work can increase the buckling load by 1.47 in 20 minutes.

Interactive topology optimisation schemes fill a practical gap by allowing users to leverage their engineering experience and expertise to increase other performance metrics, while keeping the optimisation formulation complexity low. An example of guided interactive infill topology optimisation is shown in Figure 1, where an initial optimisation is done in Figure 1(a) and afterwards, the user draws the infill pattern in Figure 1(b) and specifies the region of interest (ROI) for the infill in Figure 1(c) on the previous optimised design. Subsequently, the user is actively involved in the interactive infill optimisation, guiding the optimised design toward desired properties or design targets, as shown in Figure 1(d).

Figure 1. Interactive infill topology optimisation initially applies (a) 50 optimisation iterations using traditional compliance topology optimisation, after which the user (b) draws an infill pattern and (c) defines the region of interest (ROI) to locate the infill pattern. The algorithm resumes using an activated appearance constraint (d) to geometrically guide the optimised material layout for the defined region.

2. Background and motivation

Interactive infill topology optimisation builds upon the emerging research in design optimisation using a computer graphics-based image analysis method called Patch Match [35,36]. Patch Match is a tool for image editing determining the nearest neighbour matches between patches of source and target images to reshuffle, complete, or manipulate the target image [35,36]. The foundation of Patch Match lies in calculating the Euclidean distance between the pixels of images [37] and is applied to image and video editing, computer graphics, and art to create exciting, unique imagery [38,39]. A complete review of patch-based synthesis for image editing is outside the scope of this paper but can be found in [40] and references therein. Due to the discretization of the design domain for density-based topology optimisation, a patch-based approach is a natural extension for image editing using pixel patches to design optimisation using finite element patches. Recent work explores the use of Patch Match for aesthetic design [41–46], the generation of infill for additive manufacturing [47], and the use of machine learning to assist neural style transfers between a source image and topology optimised material layout [48]. However, interactive infill topology optimisation introduces a unique approach allowing one to generate the source image and decide where the appearance constraint should be applied within the design domain.

The introduction of infill in topology optimisation is a proven method for improving manufacturability [11,13–15] or alternative engineering performance metrics such as buckling [12]. Topology optimisation of shell/infill or coated structures creates porous, lightweight, interconnected structures that do not require support material and that have lower residual stresses during additive manufacturing [49]. Designs that can be printed without supports avoid the potential defects caused by support connection points, as well as the extra manufacturing time, sacrificial material, tools, and labour required for the removal of support parts. Recent work in multiscale topology optimisation explores the simultaneous design of the macro and microscale for lattice structures [50], additive manufacturing [51], and robust load cases [52]. Unlike traditional topology optimisation, these multiscale methods explicitly control feature sizes by mapping infill patterns within a coated shell or macrostructure [51]. Engineering the microscale of structures through the design of mechanical metamaterials for unique functions, such as maximising specific energy absorption, can be done using triply periodic minimal surfaces (TPMS) [53–55], or through multiscale topology optimisation, which is discussed in detail in [56].

Another established method for generating these infill structures and improving manufacturability in topology optimisation is through length scale control. Increasing the minimum length scale of void and solid phases allows for the prevention of small void regions or thin members that are difficult to additively manufacture [57–62]. Control of the maximum length scale generates porous infill structures reducing the need for support material [58,62–65]. Length scale modifications to improve manufacturability and generate infill structures can also be generated for topology optimisation using the level set method [13,66,67]. Few algorithms explore the spatial variation of length scale controls, which imposes a maximum or minimum feature size in the design for localised regions to target certain geometrical material layouts, such as reducing the maximum stress at a sharp corner [68]. This paper builds upon the spatial variation of length scale control by allowing control over where the infill will be located and locally, what the material layout should resemble.

There exist multiple methods to incorporate geometrical user defined design targets in topology optimisation. Algorithms such as iBESO (Bidirectional Evolutionary Structural Optimisation) and HiTop (Human Informed Topology Optimisation) applies human engineering experience and expertise to the design process using two different optimisation methods. iBESO uses subjective scoring of designs to guide the next set of designs generated by the BESO algorithm until the design satisfies the targets [69]. There is an additional feature within iBESO that allows the users to draw a feature, which the BESO algorithm will attempt to replicate in the optimised design [69]. HiTop is a density-based topology optimisation approach that pauses pre-convergence to allow the user to modify the length scale requirements for certain regions of interest (ROIs) across the domain [34,70]. Commercially available software, such as Tosca or Fusion, have topology optimisation plug-ins including manufacturability features that allow for similar spatial variations in length scale requirements [71]. Unlike the interactive topology optimisation approaches, one must define the prescribed regions for the modified feature sizes prior to running the optimisation and the feature sizes cannot directly be modified during the optimisation iterations.

2.1. Existing methods

The interactive infill topology optimisation presented herein is a novel approach that builds upon the principles of interactive topology optimisation, first presented as iBESO and HiTop. The iBESO approach enables users to draw patterns that the algorithm will explicitly match in the final design [69]. In contrast, HiTop [34,70] enables user selection of ROIs where local feature-size specifications are updated during the optimisation iterations. These interactive approaches have been shown to provide an alternative method of improving performance metrics not directly included in the optimisation formulation without increasing complexity. Topology optimisation implementations that are specialised, for example dynamics or large deformations applications, will attain designs having higher optimised performance, but these solutions frequently require high computational costs and heuristic approaches as artificial numerical stabilisation of void regions [72,73]. Patch Match topology optimisation algorithms that universally replicate source images across the design domain also experience increased computational costs for the void regions [41].

Established numerical implementations and commercial software accommodating manufacturing constraints often lack the flexibility of spatial local variations in the manufacturing requirements. Frequently for commercial software, one must partition the design domain as preprocessing and select partitioned regions to impose an overhang constraint, minimum and/or maximum feature size requirement prior to running the optimisation. In these scenarios, one must predefine and determine how to influence the optimised design without prior knowledge of how an optimised design might look or which regions could be most critical to other design requirements.

2.2. Key contributions

The present development of guided optimised designs having selective regions of drawn infill is novel, as unlike established academic approaches or commercial solutions, it attains a higher degree of design freedom to consider other design requirements that are not directly formulated in the objective function or constraints. Interactive infill topology optimisation increases geometrical control of the design by the user and does not only allowing users to draw a pattern to guide the design layout as in [69] but also to decide the specific regions for the infill and how strongly the drawn pattern is replicated in the final optimised design. Using the guided interactive infill pattern, one chooses where to locate the infill and will likely not choose all void regions, thereby avoiding a significant increase in computational costs, as in [41]. This additional design control allows one to target other performance metrics according to anticipated load types, failure criteria, and manufacturing considerations not included directly in the optimisation formulation without increasing the computational costs and modelling complexity for non-linear static structural transient analysis and the corresponding optimisation formulations. For example, if one anticipates an impact load near the surface of a structure, then one may enforce a porous infill pattern only at the region close to the surface to increase the crashworthiness energy absorption without considerably affecting the overall stiffness and strength near the supports for the static load cases.

The guided interactive infill approach differs from previous interactive topology optimisation works by placing the drawn infill patterns using an appearance constraint during the optimisation iterations. Using a simple binary field description and appearance constraint, without the complexity of applying nonlinear weighting functions, the interactive infill approach avoids the undesirable artifacts and intermediate densities seen by many length scale projection schemes. Thereby, the present solution yields clearer solutions using less optimisation iterations. It is also more intellectually approachable for novice topology optimisation users looking to generate customised designs engineered with specific infill to improve desired performance metrics. Unlike established academic solutions or commercial software, interactive infill topology optimisation provides control to the human designer, often without a significant loss in optimality for the objective function and constraints compared to the inherent optimisation formulation. Methods and costs for post-processing the optimised topology for manufacturability can be eliminated or reduced as manufacturing considerations and decisions can be added and modified during the interactive infill topology optimisation process [70].

Initially, the paper will describe a passive approach for customised infill topology optimisation in Section 3, then introduce the workflow, mathematical formulation, and numerical implementation for the guided interactive infill topology optimisation in Section 4. Section 5 demonstrates the interactive design process on two 2D benchmark examples using the drawn customised infill for improving engineering performance metrics which are not directly included in the objective function and constraints. Section 6 discusses how the present approach can reduce overhang support material for additively manufactured test beams.

3. Passive infill topology optimisation

Passive infill topology optimisation defines elements within a selected Region of Interest (ROI) to be solid or void according to a user-drawn pattern as shown in Figure 2. These elements are not included as design variables in the optimisation; therefore, this approach is better suited for fulfilling aesthetic design requirements or accommodating strict mechanical assembly requirements, for example at interfaces. The passive infill topology optimisation approach is demonstrated on an aesthetic application using a half symmetric, MBB beam, where the aim is to incorporate the logo shown in Figure 2(a). The optimisation begins by running 50 iterations employing the 88-line topology optimisation code [74], after which the user selects an ROI, shown by the blue circle in Figure 2(b) and draws a design, shown in Figure 2(c). Afterwards, the optimisation workflow resumes having the elements within the enclosed red drawing, the ‘cg,’ being set to solid, and all other elements within the blue ROI are set to void. Thereby, both the solid and void elements in the ROI are not defined as design variables. The user is enforced to draw a pattern containing information on how the pattern can connect to the surrounding design elements outside the ROI. If the pattern does not include connective members, then the design will float within the structure, unable to influence the structural properties of the optimised design. The final output is shown in Figure 2(d), preserving the passive drawing region and still minimising the compliance as the objective function, while satisfying the volume constraint.

Figure 2. Aesthetic application of passive infill implementation: (a) design inspiration, (b) user drawn region of interest (ROI), (c) passive infill pattern within the ROI having predefined interface connections to the design elements, and (d) converged optimised design.

Although this example is for purely aesthetic design exploration, it is analogous to the use of ‘frozen regions’ in industrial topology optimisation software. The passive approach can be useful for designing optimised parts in an assembly requiring certain connecting members and holes at predefined interfaces.

4. Guided interactive infill topology optimisation

Unlike the passive infill topology optimisation implementation described in Section 3, the interactive infill topology optimisation applies an appearance constraint to replicate the pattern while still including the elements of the infill pattern as design variables. Having the infill elements as design variables allows the optimisation algorithm to balance between satisfying the appearance constraint for the infill pattern and minimising the objective function and fulfilling the constraints. A specified parameter for the infill topology optimisation interactively controls the tradeoff between the geometrical appearance of the infill pattern and structural performance defined by the objective function and constraints, allowing one to determine which metrics are prioritised during the optimisation process. This approach can improve other performance targets not directly present in the objective function and constraints, for example stress, buckling, or energy absorption, in a fraction of the time. In the present work, the interactive infill topology optimisation applies the classic topology optimisation formulation minimising compliance subject to a volume constraint. The present approach is not limited to compliance topology optimisation but can be applied for minimising mass as the objective function, subject to stiffness or strength requirements. However, as previously emphasised, it may be desirable to influence other design requirements not directly included in the optimisation formulation. Therefore, for demonstration purposes, the paper focuses on minimising compliance subject to a volume fraction and influences alternative design requirements, which in the literature would often be addressed by having these directly formulated in the objective function or constraints.

4.1. Optimisation workflow

The workflow for the guided interactive infill topology optimisation is shown in Figure 3. Initially, the design domain is defined in in Figure 3(a) for the specific loads, boundary conditions, and non-designable regions. Next, the user draws half of a symmetric infill pattern shown in Figure 3(b), being converted into a full infill pattern in Figure 3(c) and the full infill pattern is repeated across the entire design domain as shown in Figure 3(d). The user only draws half of an infill pattern to generate a symmetric infill patch, which ensures continuity across the infill and optimised design. The conventional compliance topology optimisation runs for 50 optimisation iterations and the present interim topology optimised result is overlaid using the infill pattern allowing the user to draw a ROI, shown in red at Figure 3(e). The interactive infill topology optimisation applies the appearance constraint only for the specified ROI enforcing the predefined pattern, yielding the final optimised design shown in Figure 3(f). As expected, when the user enforces the drawn infill pattern, the appearance constraint is activated and the compliance increases, shown in Figure 3(g) at the highlighted region. At the start of the optimisation and at iteration 50 which activates the appearance constraint, the compliance value is lower than the converged compliance value due to those optimisation iterations violating the volume constraint. As the volume and appearance constraint are satisfied, the compliance lowers to a value comparable to the 50 iteration output. The irregularities in the compliance value, appearance constraint, and volume constraint at iterations 100 and 200 are due to continuation of parameters controlling the clarity of the optimised design.

Figure 3. (a) Defined design space for compliance optimisation. (b) Half drawing of a symmetric infill pattern, (c) generate full infill patch, (d) reflect the infill patch across the entire design domain, and (e) at 50 optimisation iterations, interactively overlay the interim topology optimised design to draw a ROI. (f) The guided interactive infill topology optimisation resumes and converges to a design applying the pattern in the ROI. (g) The convergence plot for the compliance, appearance constraint, and volume constraint, where the highlighted portion and hand figure denotes the optimisation iteration where the user enforced the drawn infill pattern.

4.2. Modifications for compatible infill

The performance of designs is highly dependent upon the applied drawings of the compatible infill patterns. Therefore, the present optimisation algorithm prompts the user to draw half of an infill pattern defined across the width of the infill patch, as shown in Figure 3(b,c). The authors found that drawing half of an infill pattern is a reliable method to prevent disconnected material and promote continuity across the infill pattern and the rest of the optimised design.

For enforcing a minimum feature size requirement, the infill pattern is passed through the classic density filter used to ensure a minimum length scale and to prevent checkerboarding of the elemental density distribution [16]. Figure 4(a) shows an example of a drawn infill pattern, while Figure 4(b) shows the filtered infill density distribution ‘smearing’ the density distribution of the pattern to satisfy the minimum length scale requirement. The infill pattern shown in Figure 4(b) is blurry due to this smearing effect upon the design layout. The filtered pattern is then passed through a Heaviside’s function to extrapolate the filtered density field toward a solid/void distribution, as shown in Figure 4(c). Afterwards, the extrapolated density distribution for the infill pattern is decisively set to 0 or 1 for each element, using a cutoff threshold of 0.5. Defining a binary 0,1 infill pattern is critical to achieve clear solutions for the optimised topologies as the appearance constraint aims to replicate the density values of the infill pattern. The infill pattern must be repeated across the entire design domain prior to filtering to ensure the connection points between the patches of the entire infill pattern also satisfy the minimum length scale requirement.

Figure 4. (a) Infill as drawn by the user, (b) passed through the density filter to enforce a minimum length scale feature size, (c) passed through a Heaviside's filter to extrapolate toward a solid/void distribution, and (d) the final binary infill pattern using a threshold cutoff of 0.5.

4.3. Topology optimisation formulation

Guided interactive infill topology optimisation follows the classic formulation shown in Equation (1), which aims to minimise the compliance subject to an equilibrium constraint, volume fraction constraint, and the lower and upper design variable bounds. However, the guided interactive infill optimisation includes the additional appearance constraint for the appearance value, $A$, which must be less than or equal to 1. Hence, the optimisation formulation is as follows [41,74]: $\underset{\mathit{x}}{\mathrm{Min}:}c\left(\mathit{\rho }\right)={\mathbf{U}}^{\mathrm{T}}\mathbf{KU}$ (1) $s.t.:\left\{\begin{array}{c}\mathbf{K}\left(\mathit{\rho }\right)\mathbf{U}=\mathbf{F}\\ {\displaystyle \frac{V\left(\mathit{\rho }\right)}{V_0}\le f}\\ \begin{array}{c}A\left(\mathit{\rho }\right)\le 1\\ 0\le x_e\le 1\mathrm{\forall }e\in \mathbf{\Omega }\end{array}\end{array}\right\}$(1) Where x is the vector of non-physical design variables, $\mathit{\rho }\left(x\right)$ is the vector of physical relative elemental densities as a function of x, simplified to ρ for notation, the objective function c is the compliance, or inverse of stiffness for a single load case, U is the vector for the nodal displacements, K is the assembled stiffness matrix, F is the vector for the nodal forces, V (ρ) is the volume of the design defined by the physical relative elemental densities, V₀ is the volume of the full design domain, f is the relative volume fraction constraint, A (ρ) is the appearance constraint value, Ω is the design domain. To encourage a solid-void solution, the Solid Isotropic Material Penalisation (SIMP) Method [17] is implemented to penalise intermediate densities using Equation (2): (2) $E\left({\rho }_e\right)=E_{min}+\left({\rho }_e{)}^{\eta }\right(E_0-E_{min})\mathrm{\forall }e\in \mathbf{\Omega }$(2) where $E\left({\rho }_e\right)$ is the Young’s modulus as a function of the physical relative elemental density ${\rho }_e$, $E_{\min }$ is a small number to ensure the stiffness matrix is positive definite, $E_0$ is the initial Young’s modulus, and $\eta$ is a penalisation factor. For this paper, $E_{\min }=1\times {10}^{-9}$, $E_0=1$, and $\eta =3$. Continuation of the penalisation factor $\eta$ increases to 5 at the 50th optimisation iteration and 7 at the 150th optimisation iteration. The design variable $x_e$, is first passed through the classic density filter to enforce a minimum feature size in the optimised design, yielding ${\tilde{x}}_e$. The filter is shown in Equations (3) and (4) [74]. (3) ${\tilde{x}}_e={\displaystyle \frac{{\sum }_{i\in {\mathit{N}}_e}{\mathit{H}}_{\mathit{e}i}\cdot x_i}{{\sum }_{i\in {\mathit{N}}_e}{\mathit{H}}_{\mathit{e}i}}\mathrm{\forall }e\in \mathrm{\Omega }}$(3) (4) ${\mathit{H}}_{\mathit{e}i}=max\left(0,r_{min}-\mathrm{\Delta }\left(e,i\right)\right)\mathrm{\forall e}\in \mathbf{\Omega },\mathrm{\forall i}\in {\mathit{N}}_e$(4) Where ${\tilde{x}}_e$ is the filtered design variable with enforced minimum length scale, ${\mathit{H}}_{\mathit{e}i}$ is a matrix of weighting factors, $r_{\min }$ is the radius defining the minimum solid feature size, and ${\mathit{N}}_{\mathit{e}}$ is the neighbourhood of elements within $r_{\min }$ of element $e$. The filter also prevents checkerboards in the optimised design caused by the SIMP method’s penalisation of intermediate densities. The filtered design variable ${\tilde{x}}_e$, is then passed through a Heaviside projection to further push its value towards 0, void, or 1, solid. For this paper, the following Heaviside projection is used [75]: (5) ${\rho }_e\left({\tilde{x}}_e\right)=1-e^{-\text{β}\cdot {\tilde{x}}_e}+{\tilde{x}}_e\cdot e^{-\text{β}}\mathrm{\forall }e\in \mathbf{\Omega }$(5) Where${\rho }_e$ is the physical relative elemental density, $\text{β}$ is the steepness factor and set initially to 5 for this work. Continuation increases the steepness factor $\text{β}$ to 15 at the 50th optimisation iteration and 25 at the 150th optimisation iteration. The physical relative elemental density ${\rho }_e$ is used for calculating the appearance constraint, volume constraint, the finite element modelling, and the corresponding compliance.

The appearance constraint is calculated using a distance formula [41] summing the squared difference between the physical elemental density and the corresponding user defined infill pattern described in Section 4.2 as the following: (6) $D_e\left(\mathit{\rho },\mathit{\alpha }\right)={\displaystyle \frac{{\sum }_{j\in {\mathit{\omega }}_{\mathit{e}}}{\left({\rho }_{e,j}-{\alpha }_{e,j}\right)}^2}{\left|{\mathit{\omega }}_{\mathit{e}}\right|}\mathrm{\forall e}\in {\mathbf{\Omega }}_{\mathbf{ROI}}}$(6) Where $D_e$ is the distance value for element e, ${\rho }_{e,j}$ is the physical relative elemental density of element j in the neighbourhood ${\mathit{\omega }}_{\mathit{e}}$ of element e, ${\alpha }_{e,j}$ is the density value for element j of the defined infill pattern, and ${\mathbf{\Omega }}_{\mathbf{ROI}}$ is the collection of elements located in the user drawn ROI. The neighbourhood containing elements $j$ for element e, ${\mathit{\omega }}_{\mathit{e}}$, is defined by one half of the width for the patch upon which the infill is drawn. For example, if the user draws the infill pattern on a patch that is 20 elements wide in the x-direction and 15 elements tall in the y-direction, then ${\mathit{\omega }}_{\mathit{e}}$ includes all elements $j$ within 10 elements of element e. A normalised sum of the squared difference between the physical relative elemental densities and the infill pattern densities allows for elements outside of the drawn ROI to influence the appearance constraint. If an element ${\rho }_{e,j}$ is not located in the drawn ROI, ${\mathbf{\Omega }}_{\mathbf{ROI}}$, but is in the neighbourhood, ${\mathit{\omega }}_{\mathit{e}}$, of element e which lies in the drawn ROI, it still impacts the elemental distance $D_e$. This facilitates continuity across the boundary of the ROI and ensures connection between the replicated drawn infill and the rest of the optimised design. The difference $D_e$ is calculated for each element e in the ROI, ${\mathbf{\Omega }}_{\mathbf{ROI}}$, generating a vector $\mathit{D}$, that is ultimately summed into $A\left(\mathit{D}\right)$. Thereby, the appearance constraint [41] yields: (7) $A\left(\mathit{D}\right)={\displaystyle \frac{1}{A^{\ast }\left|{\mathbf{\Omega }}_{\mathbf{ROI}}\right|}\sum _{e\in {\mathbf{\Omega }}_{\mathbf{ROI}}}{D}_e\left(\mathit{\rho },\mathit{\alpha }\right)}$(7) Where the appearance constraint $A\left(\mathit{D}\right)$ is a function of $\mathit{D}$, the vector of elemental distances of the drawn infill pattern, ${\mathbf{\Omega }}_{\mathbf{ROI}}$ is the ROI, and $A^{\ast }$ is a parameter controlling how strictly the appearance constraint is enforced. For low $A^{\ast }$, the distance between $\mathit{\rho }$ and $\mathit{\alpha }$. is small, generating a very exact replication of the infill pattern. For this work, $A^{\ast }$ is calculated using the scheme outlined in [41] defined as: (8) $A^{\ast }=\gamma \cdot A_{50}$(8) Where $A^{\ast }$ is the parameter controlling the rigidity of the appearance constraint, $\gamma$ is a user defined value, set to 0.3 for the numerical examples in Section 5 and varied for the numerical examples in Section 6, and $A_{50}$ is the appearance value for the defined drawn infill pattern and the interim topology optimised degn at 50 optimisation iteration using the 88 line code [74].

4.4. Sensitivity calculations

Appearance constrained topology optimisation is a gradient-based approach that requires the calculation of sensitivities for the objective and constraints with respect to the design variable, $x_e$. The sensitivity of the compliance with respect to the design variables can be found using the chain rule as shown in Equation (9). (9) ${\displaystyle \frac{\mathrm{\partial }c}{\mathrm{\partial }{x}_e}=\sum _{j\in \mathbf{\Omega }}{\displaystyle \frac{\mathrm{\partial }c}{\mathrm{\partial }{\rho }_j}{\displaystyle \frac{\mathrm{\partial }{\rho }_j}{\mathrm{\partial }{x}_e}\mathrm{\forall e}\in \mathbf{\Omega }}}}$(9) The first term is the derivative of the compliance with respect to the physical relative elemental density and yields Equation (10) for the SIMP method [74]. (10) ${\displaystyle \frac{\mathrm{\partial }c}{\mathrm{\partial }{\rho }_j}=-\eta \cdot {\rho }_j^{\left(\eta -1\right)}\cdot \left(E_0-E_{min}\right)\cdot {\mathit{u}}_{\mathit{j}}^{\mathit{T}}{\mathit{k}}_0{\mathit{u}}_{\mathit{j}}\mathrm{\forall j}\in \mathbf{\Omega }}$(10) The second term in Equation (9) is the derivative of the physical elemental density with respect to the design variable, which is similarly found using the chain rule in Equation (11). (11) ${\displaystyle \frac{\mathrm{\partial }{\rho }_j}{\mathrm{\partial }{x}_e}={\displaystyle \frac{\mathrm{\partial }{\rho }_j}{\mathrm{\partial }{\tilde{x}}_j}{\displaystyle \frac{\mathrm{\partial }{\tilde{x}}_j}{\mathrm{\partial }{x}_e}\mathrm{\forall }e\in \mathbf{\Omega },j\in \mathbf{\Omega }}}}$(11) The first term in Equation (11) applies the derivative of the Heaviside’s projection with respect to the filtered density${\tilde{x}}_e$ and is calculated using Equation (12) [74]. (12) ${\displaystyle \frac{\mathrm{\partial }{\rho }_j}{\mathrm{\partial }{\tilde{x}}_j}=\text{β}\cdot e^{-\text{β}\cdot {\tilde{x}}_j}+e^{-\text{β}}\mathrm{\forall j}\in \mathbf{\Omega }}$(12) The second term in Equation (11) takes the derivative of the classic density filter with respect to the design variable $x_e$ and is calculated using Equation (13) [76,77]. (13) ${\displaystyle \frac{\mathrm{\partial }{\tilde{x}}_j}{\mathrm{\partial }{x}_e}={\displaystyle \frac{{\mathit{H}}_{\mathit{e}j}}{{\sum }_{i\in {\mathit{N}}_{\mathit{e}}}{\mathit{H}}_{\mathit{e}i}}\mathrm{\forall e}\in \mathbf{\Omega },j\in \mathbf{\Omega }}}$(13) Similarly, the appearance constraint is calculated using the chain rule, shown in Equation (14), where the derivative of the appearance is found with respect to the physical elemental density. (14) ${\displaystyle \frac{\mathrm{\partial }A}{\mathrm{\partial }{\rho }_e}=\sum _{j\in {\mathbf{\Omega }}_{\mathbf{ROI}}}{\displaystyle \frac{\mathrm{\partial }A}{\mathrm{\partial }{D}_j}{\displaystyle \frac{\mathrm{\partial }{D}_j}{\mathrm{\partial }{\rho }_e}\mathrm{\forall }e\in {\mathbf{\Omega }}_{\mathbf{ROI}}}}}$(14) The first term of Equation (14) is found using Equation (15), applying the derivative of the appearance constraint with respect to the elemental distances [41]. The sensitivity of the appearance constraint with respect to elements outside of the ROI is set to 0 to ensure compatibility with Equations (11–13) without activating the appearance constraint outside of ${\mathbf{\Omega }}_{\mathbf{ROI}}$. (15) ${\displaystyle \frac{\mathrm{\partial }A}{\mathrm{\partial }{D}_j}={\displaystyle \frac{1}{A^{\ast }\left|{\mathbf{\Omega }}_{\mathbf{ROI}}\right|}\mathrm{\forall }j\in {\mathbf{\Omega }}_{\mathbf{ROI}}}}$(15) The second term of Equation (14) is determined using Equation (16), which calculates the derivative of the elemental distance with respect to the physical elemental density. (16) ${\displaystyle \frac{\mathrm{\partial }{D}_j}{\mathrm{\partial }{\rho }_e}={\displaystyle \frac{{\sum }_{j\in {\mathit{\omega }}_{\mathit{e}}}2\left({\rho }_{j,e}-{\alpha }_{j,e}\right)}{\left|{\mathit{\omega }}_{\mathit{e}}\right|}\mathrm{\forall e}\in {\mathbf{\Omega }}_{\mathbf{ROI}}}}$(16) Afterwards, the optimisation problem in Equation (1) is solved using the derivatives of the objective function and constraints defined in Equations (9–16) and the Method of Moving Asymptotes (MMA) [18] which is a mathematical programming method for structural optimisation.

5. Numerical examples

5.1. Increase linear buckling resistance

The guided interactive infill topology optimisation is demonstrated on a 2D benchmark short cantilever example as shown in Figure 5(a). The objective is to minimise compliance subject to a volume fraction, $f=0.15$. The number of elements in the x  – and y  – direction is 90 and 210, respectively. The number of elements in the pattern patch in the x  – and y  – direction is 20 for both directions and the appearance constraint parameter$\gamma$ is 0.3.

Figure 5. (a) Design domain for a short cantilever design domain. (b) traditional compliance optimised structure, (c) user drawn infill for the patch, (d) full infill pattern for the entire design domain, (e) user drawn ROI on overlaid infill pattern and the interim 50-iteration optimised design, and (f) converged optimised design for guided interactive infill optimisation.

The short cantilever is first designed using traditional compliance topology optimisation, as shown in Figure 5(b). However, the structure is prone to buckling due to the slenderness ratio of the lower member being in compression. The user would like to reduce this failure mode and thereby, increase the maximum buckling carrying capacity of the structure by enforcing a customised infill in the lower member of the structure. Therefore, one draws an infill pattern shown in Figure 5(c, d) and chooses the lower region ROI shown in Figure 5(e). The final optimised design is shown in Figure 5(f) for the entire workflow using the guided interactive infill optimisation approach. The small holes in the final, converged design of Figure 5(f) could pose a challenge to manufacturability but can be mitigated by imposing a two-phase minimum length scale of the solid and void material phases.

Both optimised structures in Figure 5(b, f) are post-processed as 2D shell finite element models using the finite element solver Abaqus [71], subject to the initial load, material properties, and boundary conditions. The guided interactive infill optimised cantilever structure in Figure 6(b) has a 39.3% higher compliance than the traditional compliance optimised structure in Figure 6(a). However, the traditional compliance optimised design is highly prone to buckling due to its high slenderness ratio. A linear eigenvalue buckling analysis demonstrates that the guided interactive infill topology optimised design yields a 147% increase in the lowest buckling eigenvalue for the dimensionless, normalised displacement plots for the given linear buckling eigenmodes in Figure 6. Although the coating material surrounding the infill pattern provides the main structural capacity for the linear buckling stability, rather than the infill itself, the material distribution shown in Figure 6(b) is the direct result of the changes by the user enforced on the material layout of the optimised design.

Figure 6. Lowest linear buckling eigenvalue and normalised displacements of the eigenmodes showing the failure location for the (a) traditional compliance optimised short cantilever and the (b) guided interactive infill topology optimised short cantilever, respectively.

Interactive infill topology optimisation can also be applied to more challenging design problems to increase the buckling resistance, such as the L-bracket shown in Figure 7(a). The objective is to minimise compliance subject to a volume fraction, $f=0.25$. The number of elements in the x  – and y  – direction is 200. The number of elements in the pattern patch in the x  – and y  – direction is 20 for both directions and the appearance constraint parameter$\gamma$ is 0.3.

Figure 7. (a) Design domain for a L-bracket. (b) Traditional compliance topology optimised design. (c) Drawn O-infill for the patch, (c) full O-infill pattern, (d) ROI, (e) converged interactive infill topology optimised design.

The traditional topology optimised L-bracket is shown in Figure 7(b). Based on their engineering intuition, a user identifies that the left upper member of the traditional optimised L-bracket would be the most likely to buckle under the loading shown in Figure 7(a). Therefore, the user draws the infill pattern shown in Figure 7(c,d) and ROI in Figure 7(e). The converged, interactively modified design is shown in Figure 7(f). As in the short cantilever example, the L-brackets are post-processed in Abaqus and yield the buckling shapes shown in Figure 8.

Figure 8. Lowest linear buckling eigenvalue and normalised displacements of the eigenmodes showing the failure location for the (a) traditional compliance optimised L-bracket and the (b) guided interactive infill topology optimised L-bracket, respectively.

Compared to the traditional optimised L-bracket shown in Figure 8(a), the interactive infill optimised L-bracket shown in Figure 8(b) has a 5% higher compliance value and 50% higher lowest linear buckling eigenvalue. The L-bracket example demonstrates that similar improvements in linear buckling resistance can be realised for more challenging design problems. However, the relative increase in the lowest buckling eigenvalue is smaller than the obtained increase in the linear buckling load for the simplistic short cantilever example shown in Figure 6.

5.2. Increase energy absorption

Introducing infill into a design is an established method of increasing the impact resistance and crashworthiness energy absorption of the structure, especially when subjected to blast and high rate impact loading. Increasing the energy absorption using infill is demonstrated on a 2D square design domain being fully clamped at the bottom edge and applied to pressure load on the top edge, shown in Figure 9(a). The number of elements for the design domain in both the x – and y – direction is 200, the number of elements in the pattern patch in both the x – and y – direction is 20, and the appearance constraint parameter γ is 0.3.

Figure 9. (a) Design domain for a short compressive square. (b) Traditional compliance topology optimised design.

The traditional compliance topology optimised design is shown in Figure 9(b) for a relative volume constraint $f=0.5$. Two interactive infill topology optimised designs are generated using two different infill patterns: an O and an X shape. These two infill patterns are chosen because they are symmetric, generating a continuous infill field, and vary significantly in their distribution of material in the small infill patch. An optimised design is generated first using an X-infill in a lower ROI as shown in Figure 10(a–d) and secondly using an O-infill in a similar ROI shown in Figure 10(e–h). The variations of infill patterns will determine the optimised compliance and energy absorption for crashworthiness and impact.

Figure 10. (a) Drawn X-infill for the patch, (b) full X-infill pattern, (c) lower ROI for the X-infill, (d) converged X-infill interactive infill topology optimised design. (e) Drawn O-infill for the patch, (f) full O-infill pattern, (g) lower ROI for the O-infill, (h) converged O-infill interactive infill topology optimised design.

The three compressive square designs are analyzed using the Abaqus solver and shell finite elements having elastoplastic constitutive steel material properties [78]. A compressive pressure load having a total of 70 kN is applied at the top edge of the structure creating force-displacement curves for the three optimised designs as shown in Figure 11, where the gray dotted line marks the maximum load.

Figure 11. Force-displacement responses and respective plastic strain magnitudes at maximum load for the traditional, X-infill pattern, and O-infill pattern optimised designs.

As demonstrated by the force-displacement responses and the corresponding plastic magnitude strain plots of the optimised designs at maximum load, the X-infill optimised design absorbs the most energy at maximum load, followed by the O-infill optimised design, and finally the traditional compliance topology optimised design. The X-infill optimised design has a 6.9% increase in compliance relative to traditional compliance topology optimisation. However, the design has 307.2% higher energy absorption. The O-infill optimised design increases the compliance by 3.0% and increases the energy absorption by 182.8%. Figure 11 demonstrates that the X-infill designs outperform the O-infill pattern for energy absorption, but both modified designs exhibit far more ductile behaviour than the traditional design. The authors explored multiple ROI locations, including locations closer to the upper edge where the load is applied, and found similar increases in the energy absorption. Unlike this example where the infill region is the location of plastic deformation, the alternative ROI locations similarly increased energy absorption even though the infill was not the direct mechanism of energy absorption. The guided patterns for interactive infill topology optimisation improve the external engineering performance metric for crashworthiness energy absorption without directly including an energy absorption measure in the objective function or constraints or applying non-linear finite element modelling. The improvement in performance is instead obtained by indirectly influencing the geometrical layout while still minimising compliance using linear finite element modelling subject to a 50% volume constraint.

6. Experimental example

The guided interactive infill topology optimisation allows one to incorporate customised infill patterns at selective regions of the optimised design, which is particularly relevant when considering manufacturing costs for additive manufacturing. Support structures are generated for the regions of the 3D printed designs where the material violates a given overhang angle, such as 45° compared to the baseplate for printing. One can reduce the overhang in the final optimised design using the interactive infill topology optimisation algorithm and reduce the amount of support structures and manufacturing costs for printing. The reduction of overhang support structures is demonstrated using the 2D MBB beam, shown in Figure 12(a) where half of the design domain is optimised. The 2D modelling is then extended to 2.5D modelling having the dimensions shown in Figure 12(b). Based on a rule of thumb established by previous, successful experimental three-point bending tests, the out of plane thickness of the beam is selected as 15 mm (10% of the beam length L = 150 mm) to prevent out of plane buckling. The dimensions shown in Figure 12(b) are used for all test specimens, therefore out of plane buckling has no impact on the experimental results and all beams exhibit flexural behaviour. The additive manufacturing uses a FDM printer where Polylactic Acid, or PLA, plastic filament is used to manufacture the optimised design and support structures. For this example, the relative volume constraint is $f=0.5$. The number of elements in the x – and y – direction is 300 and 100, respectively. The number of elements in the pattern patch in the x – and y – direction is 20 and 30, respectively, and the appearance constraint parameter $\gamma$ is varied.

Figure 12. (a) Half symmetrical 2.5D design domain for the MBB beam topology optimisation and (b) the full dimensions of the beam imposed to additively manufacture the 2.5D optimised designs.

Three beams are designed using the drawn infill shown in Figure 13(a,b), the three ROIs shown in Figure 13(c), and the appearance constraint parameter γ = 0.3, 0.5 and 0.7, respectively. By having the ROI and infill pattern constant for the three optimizations, these numerical and experimental studies explore the effect of the appearance constraint with respect to two design requirements: minimised mass of the support structures for reduced printing costs versus minimised compliance. The three guided interactive infill optimised designs are shown in Figure 14(c–e) and compared to the traditional topology compliance optimised design shown in Figure 14(a), as well as the optimised design subject to with an overhang constraint shown in Figure 14(b) [6]. Lower γ value yields a more rigid appearance constraint for replicating the X-infill in the three ROIs (Figure 14(c)), while a higher γ value (Figure 14(e)) loosely applies the patterns creating a design more similar to the traditional compliance topology optimised design.

Figure 13. (a) The user draws a symmetrical infill pattern, it is (b) replicated to the full design domain and then repeated across the interim 50-optimisation iteration design upon which (c) the user defines three ROIs for achieving the guided interactive infill optimised designs.

Figure 14. Topology optimised designs considering additive printing requirements printed with support structures (orange printed material), and post-manufacturing removing the support structures from the optimised designs obtained using (a) traditional compliance topology optimisation, (b) overhang constraint, and interactive infill enforced using (c) $\gamma =0.3$, (d) $\gamma =0.5$, and (e) $\gamma =0.7$.

Three specimens of each of the five designs were additively manufactured using natural colour Ultimaker PLA material on an Ultimaker S3 FDM printer using G-code generated by Cura for slicing the designs. The G-Code parameters include 100% infill, 70° support angle, and tree branch support, structures are printed using the orange Ultimaker PLA material. The 2.5D printed designs were tested on a 100 kN Instron 3-point bending testing machine as shown in Figure 15. The experimental compliance value was calculated from the slope of the force vs. displacement data in the beam’s elastic region, which is provided by the Instron testing machine.

Figure 15. Experimental set-up for the 3-point bending test using a 100 kN Instron machine.

The results of the numerical and experimental studies are shown in Table 1 and Figure 16. The numerical compliance predictions are within 3% of the experimental results, demonstrating excellent agreement between the loads and the boundary conditions used in the optimisation and the experimental 3-point bending test set-up. In general, as γ decreases and enforces a more rigid appearance constraint for the infill pattern, as expected, the compliance increases and the material for the support structures decreases. The correlation between the amount of material for the support structures and the compliance values suggest that a value of γ = 0.5 decreases the support material by 74% while maintaining stiffness performance within 3.1% of the traditional compliance topology optimised design. Additionally, the γ = 0.5 design outperforms the overhang constrained design by 0.7%, although the overhang design requires no support material. As previously mentioned, interactive infill topology optimisation provides a computationally cheap alternative method to improve performance metrics not included in the optimisation formulation. However, specialised numerical implementations will provide enhanced specific performance, as seen in the present example where the overhang has the largest reduction for the support material.

Figure 16. Comparison of compliance for the numerical results and physical experiments, as well as the mass required for support structures in the additive manufacturing process.

Table 1. Numerical and experimental compliance results for different topology optimised designs considering additive manufacturing requirements.

	Traditional	Overhang constraint	Infill pattern$\gamma =0.3$	Infill pattern$\gamma =0.5$	Infill pattern$\gamma =0.7$
Mass of supports (g)	2.53	0.00	0.57	0.67	2.43
Mass of printed beams without supports (g)	32.00	31.87	32.03	32.47	31.90
Numerical compliance (mm·kN)	1.53	1.61	1.61	1.58	1.54
Average experimental compliance (mm·kN)	1.51	1.57	1.60	1.55	1.55

7. Conclusions

The interactive infill topology optimisation enables the user to guide the design generation toward desired design requirements not directly included in the optimisation problem formulation. Instead, based upon engineering experience and design expertise, the present approach uses an appearance constraint to guide the material layout interactively. Established optimisation approaches can address challenging objectives and constraints, such as minimising stress [19–25] or maximising the lowest linear buckling eigenvalues [26–29], or maximising the crashworthiness energy absorption [30–33], and the designs of these approaches will have a better structural specific performance than those of the designs obtained using the guided interactive infill topology optimisation approach. However, the existing approaches often apply additional design parameters requiring tuning, as well as increased computational costs for both the finite element analysis and the number of optimisation iterations. The guided interactive infill topology optimisation is computationally cheaper as the modelling is static and linear elastic and is thereby approachable for novice users or those without computational experience in advanced modelling and optimisation but practical experience and know-how.

The present work first introduces a passive implementation for infill topology, allowing one to draw and prescribe a pattern in the final topology-optimised design. This control generates creative, unique aesthetic designs or creates custom void or solid regions for optimised parts nested within assemblies. Afterward, guided interactive infill topology optimisation is introduced and demonstrated on three numerical 2D optimisation benchmark examples and on one 2.5D example, where the optimised designs are examined both numerically and experimentally. It is worth noting that the proposed algorithm is formulated as a standalone differentiable function evaluated on a density field and added to the topology optimisation workflow. The guided interactive infill topology optimisation is demonstrated with compliance minimisation and a mass constraint for simplicity, but the method is, in principle, compatible with any other physical and geometrical objective or constraint found in the literature, as well as advanced boundary conditions such as kinematic couplings and multiple load case applications. It is also possible to apply the present interactive infill approach for multi-material topology optimisation. This could be achieved by the user drawing ROIs that selectively modify the Young’s Modulus throughout the domain, or by allowing infill patterns to have varying densities corresponding to the different materials in the design domain.

A limitation of the guided interactive infill topology optimisation is its difficulty in extending the approach to 3D applications. Although conceptually possible, future work will examine how to draw 3D infill patterns and intuitively select ROIs for 3D designs. Additionally, further exploration of the optimisation iteration for which the user should interactively modify the design requirements is needed, rather than the present default of using 50 optimisation iterations, as well as which $\gamma$ appearance values are best suited for various applications. Lastly, as with any interactive algorithm, the present approach relies upon knowledgeable decisions by the user for the infill patterns and the different ROI locations. The present modifications to the workflow enforcing a symmetric infill pattern and the filtering of the infill to enforce a minimum feature size requirement ensure compatible infill drawings but do not guarantee high structural performing designs.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Data availability statement

The authors confirm that the data supporting the findings of this study are available within the article.

Additional information

Funding

This material is based upon work supported by the National Science Foundation Graduate Research Fellowship Program under [grant number 2141064]. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation.