Combination of topology and shape optimization with finite element modeling in the case of an aerospace component produced by laser based additive manufacturing Open Access

Additive manufacturing (AM) has received considerable attention in the past two decades, although it had been present quite earlier, e.g., in the procedures to construct free-form topographical maps and photosculptures with 2D layers^1–3 even more than a century ago. In the 1960s, with the advent of computer aided design (CAD) and manufacturing (CAM), with processes such as powder fusion,⁴ photopolymerization,⁵ and sheet lamination,⁶ the groundwork for modern AM processes⁷ was set. Nowadays, AM processes are present in a wide range of sectors, like medical, space, and aeronautical.⁷ 

The wide usage of AM technologies relies on a set of advantages it provides over conventional manufacturing technologies (e.g., milling, drilling, or casting). It allows for increased geometrical complexity, especially important during light-weighting processes; a reduction in material waste due to the additive nature of the processes; and high customization, which allows for specialized, on-demand, printed parts. However, AM is not the choice for many applications due to its drawbacks and unsolved problems.⁸ These are, for instance, build size limitations, constrained by dimensions of the build chamber or low deposition rates; the production and hardware costs, which are higher than conventional manufacturing methods as the production scale increases; anisotropic mechanical performance due to layer-by-layer manufacturing; postprinting part distortion and cracking; and process complexity reliant on many parameters that impact the final part quality and mechanical performance.

Selective laser melting (SLM) is one of the most common AM processes and is featured in the current work. It involves the operation of a laser that selectively fuses layers of a powdered material along a predefined trajectory. As such, it is a powder-bed-based AM technology that allows the production of a part with a fine microstructure and similar properties to the ones obtained from conventional methods. However, the process parameters such as power, scanning speed, beam diameter, or cooling times are tightly interconnected and, thus, heavily influence the quality of the printed part.⁹ 

In order to maximize the quality of parts produced, a set of guidelines and a structured framework are necessary. These are conventionally named design for additive manufacturing (DfAM) guidelines and focus on the printing steps of a part, be it during its design and conceptualization stage, the printing process itself, or postprocessing procedures (e.g., support removal and heat treatments). The present work focuses on the application of structural optimization and lightweighting techniques, given the ability of AM technologies to print highly complex structures. Topology optimization (first introduced by Bendsøe and Kikuchi in 1988¹⁰) is the main optimization tool utilized. It allows us to further include AM-specific considerations (constraints) that impact a part’s performance and printability. These DfAM constraints can be subdivided into so-called directional and nondirectional categories. The first set is directly influenced by the print direction and includes anisotropy, support structures design and optimization, and the presence of distortion and residual stresses. The nondirectional constraints are not influenced by the build direction and include power enclosure and feature size control. To achieve a more complete implementation of DfAM guidelines, all constraints should be considered and implemented during the topology optimization designing phase. However, considering the application of the methodology in the context of the present work, with relevance to the aeronautical sector, where maximum weight reduction is sought-after, support structure optimization and anisotropy considerations are not included in the present formulation.

Alongside the structural optimization procedure, modeling of distortion and residual stress prediction is important, in order to achieve good quality in printing and to better understand and optimize the parameters that control the printing process.⁸ As such, a numerical AM process simulation finite element (FE) model is developed and applied to a highly complex geometry based on an uncoupled thermal-mechanical procedure, providing an extra design validation step.

A more comprehensive description of the current work can be found elsewhere.¹¹ 

The AMATHO project (Additive MAnufacturing Tilt-rotor HOusing) tackles directly the first goal set by Clear Sky 2, the European Union initiative for sustainable aviation. By taking advantage of the additive manufacturing potential, it aims at the design and manufacture of gearbox housing to be utilized by Leonardo Helicopter’s next generation civil tilt-rotor (NGCTR). The integration of AM in the design process allows for reduced weight and optimized performance, saving resources during both manufacturing and flight time.

The target problem of the present work is the optimum design of an aeronautical component to be employed by Leonardo Helicopters, utilizing SLM as the AM process.

The currently used component assembly is made of two individual parts, and it can be seen in Fig. 1, highlighting an upper and lower manifold.

Three servo-actuators are connected to the part through 18 mounting bosses, six on each of the three faces. Both load application points (⁠ $A F S i$⁠) and servo-actuator centers of mass (⁠ $C G i$⁠) are provided. Static strength of the component is defined under two loading conditions (axial and bending) and three load intensities for each. Under proof stress loading intensity, the part shall not experience any yielding. While not directly assessed in this work, fatigue loading is also used to justify the definition of two additional bending load cases, in order to achieve topology optimization solution symmetry. The servo-actuator is subjected to a minimum axial stiffness requirement in the z-direction, whereas the mounting bosses themselves are subjected to both a relative and normal stiffness requirement. Finally, a minimum first natural frequency is also defined.

To the above-mentioned problem, a DfAM-guided approach is employed to obtain a design compliant with requirements defined by Leonardo Helicopters, with the goal of obtaining a matching or preferably lower weight, as compared to the currently used design. First, a topology optimization task is carried out based on the limiting geometry, design requirements, and loading conditions. The interpretation of the results is carried out by using polyNURBs in Altair Inspire. Following design interpretation, such part is subjected to same loading conditions and modal analysis to verify the compliance with constraints. If the result is not satisfactory, the topology optimization parameters are adjusted, until a design is chosen. Such resulting part is then subjected to a shape optimization procedure to address static strength requirements and obtain stress homogeneity (with added durability in the case of fatigue loading). Identically, the shape optimization process is carried out until a satisfactory result is obtained and could be composed of several iterations, with the goal of minimizing stress peaks, thus extending fatigue life of the part. Finally, a thermo-mechanical FE simulation is carried out to highlight the part’s performance under the complex thermal-mechanical process it undergoes, namely, regarding distortion and residual stresses. Figure 2 is the proposed design workflow for this project.

Structural optimization has the objective of finding the best way an assembly of materials may be set, to withstand a set of loading and boundary conditions. However, the definition of the goal may vary depending on the project’s objectives. As an example, one may seek a minimum weight approach. Therefore, the definition of the structural optimization problem involves the search of the assembly of material that, while respecting a given set of requirements (stiffness, weight, eigenvalues, displacement, etc.), optimally sets the objective of interest, which could be, for example, the stiffness or weight of the said structure.

A general structural optimization problem is defined by an objective function [ $f(x,y)$⁠, representing the objective to be optimized], a design variable (⁠ $x$⁠, which describes the design and is typically the geometry of the structure), and the state variable [ $y(x)$⁠, which represents the structure’s response for a given design variable].

In order to obtain a design with well-defined voids and solid regions (1-0 design), intermediate density values are to be steered toward the two extremes. This is achieved through the use of interpolation methods, with the most commonly used being the SIMP method. However, SIMULIA’s recently introduced mass interpolation material penalization (MIMP) method¹³ is of interest for the problem at hands and is implemented in this work.

The MIMP algorithm is based on a combination between SIMP for the constitutive material interpolation and a SIMULIA-developed physical density material interpolation. Its suggested use is in strength optimization with applied stress, either as an objective or as a constraint. Additionally, it reportedly improves results with mass-related responses, such as natural frequency, since previously the algorithm would by default change the interpolation scheme to PEDE (the so-called Niels–Pedersen-based approach, also proprietary to SIMULIA¹⁵).

The shape optimization applied in this work is based on a nonparametric approach.^16,17 It is determined on implicit parameters, defined from the set of surface design nodes from the FE model, namely, the displacement vector. The minimization of the maximum objective formulation is considered, in order to achieve overall higher robustness.¹⁵ 

The AM simulation problem is defined by coupling between three main fields, the thermal, mechanical, and metallurgical. The metallurgical field, while important, for example, due to the dependence of distortion and residual stress on solid phase transformations that occur, will have its implementation limited in the context of the present work, with the focus being on thermal and mechanical problems. It can be stated that the plastic deformation-induced friction heat is negligible, when compared to the thermal energy present in the system, leading to a weakly coupled or uncoupled model that is currently employed in most analyses.^18–20 

The accurate modeling of the SLM process requires the definition and consideration of a set of process elements, these being toolpath/scan path information, the progressive material deposition information, the thermal-structural material properties, and the laser/heat input properties. Software dedicated to AM simulation already exists and it is used. However, it is commonly associated with specific machines, resulting in license and data format limitations. In order to overcome this difficulty, both Abaqus with its built-in user subroutines and mesh intersection tools and Autodesk Netfabb were utilized to define the scan path based on AM machine printing parameters.

The representation of both the laser source and the recoater roller is done through a toolpath. The main geometries of interest are the point and infinite line, respectively, for the laser source and the recoater roller. However, others can be considered for different applications, such as the box geometry, useful to characterize material deposition in polymer extrusion or wire-feed methods.

The event series is a functionality in Abaqus that define time- and space-varying fields, independently of both the mesh and the simulation time increment. The field definition is done through a step function representation, linearly discretized in both space and time. Any event series has a fixed format set for both mandatory time and space definitions, to which a maximum of six user-defined fields or dependent variables can be added.

Two event series are of particular interest in the simulation of the powder bed AM method, namely, for the definition of the recoater roller motion, and the laser power, both defined in both time and space. The laser event series is used to define the laser power from which the heat flux can be calculated. On the other hand, the roller event series is used to progressively activate elements and assumes a boolean dependent variable equal to 1 when the roller is passing and 0 otherwise.

Progressive element activation in Abaqus is employed to simulate the material deposition event that occurs in powder bed-type AM processes due to layer-by-layer raw material deposition from each recoater roller passage, and it had been used previously,^22,23 e.g., in the case of laser cladding.^24,25

The path taken by the roller is represented by an infinite line moving along the same direction as the roller, defining a plane parallel to the new layer and perpendicular to the printing direction. Such infinite line intersects the finite element mesh, either partially or fully activating said elements.

The newly deposited layer of the raw material is fused due to the passage of the heat source, characteristic of the AM process of interest (laser, electron beam, etc.). This makes the modeling of the heat source essential to properly define the AM process, mainly through its power, scan path, and heat distribution type.

User subroutines are implemented to increase the capabilities of Abaqus to simulate the AM process.¹⁵ Two particular subroutines are utilized to aid the simulation of both the moving heat flux and the progressive element activation process. These are UMDFLUX and UEPActivationVol and the respective accompanying data setup subroutines UMDFLUXSetup and UEPActivationSetup.

UMDFLUXSetup provides required data structures to the mesh intersection tool, such as mesh and event series details related to the moving heat source at the start of the step. At the beginning of the increment, it is once again called to compute the heat flux per element based on intersection between the event series laser path and the mesh.

UMDFLUX is responsible for the definition of heat flux distribution of multiple moving heat sources, continuously with respect to time and space. The subroutine accounts for different heat source models, such as the Goldak heat source model²⁶ and the concentrated point energy source model. It allows for a continuous definition of one or more flux sources that do not depend on the mesh and instead are based on the event series information. For each element, the subroutine obtains the total flux for each event provided by the mesh intersection tool during UMDFLUXSetup.

UEPActivationSetup, at the beginning of the step, creates and provides to the mesh intersection tool the required data structures to perform progressive element activation, with accordance to both the mesh and the recoater roller event series. On the other hand, at the start of the increment, it triggers the mesh intersection tool to compute the intersection between the recoater roller event series and the mesh for the specific increment time interval. The volume fractions are stored by the mesh intersection tool.

UEPActivationVol obtains the volume fraction previously computed by UEPActivationSetup at the start of the increment for each element of interest. The volume fraction increase defines the partial element activation process.

The diagram in Fig. 3 highlights the order of interaction between user subroutine calls and the mesh intersection tool during the Abaqus analysis.

The AM modeler plug-in is developed and made available by SIMULIA, with the main goal of aiding the definition of the necessary keywords that coordinate extensive data needed by user subroutines. It requires the previous definition of both thermal and structural models, both the heat source path and recoater roller event series, and simulation parameters.

The geometry preprocessing process is divided into three major steps: printing support design and generation, voxelization, and part slicing and scan path information.

With the CAD geometry of the part to be printed, first the supports need to be designed and generated alongside the part. This is done utilizing Autodesk Netfabb’s automated tools. Following this procedure, the new mesh file is rasterized in order to create a voxelized representation of the part and supports. The result of the procedure is a regular hexahedral element-based mesh, constructed along the printing direction. Finally, the resulting mesh can be exported as a STL file and imported again in a slicing software to obtain the scan part information.

Voxelization of a continuous geometry consists of the discrete subdivision of such a geometry into cells or voxels. In the present case, each voxel is defined by the same geometry of a regular hexahedral finite element and the volume fraction property. The volume fraction property is a scalar discrete field that characterizes the intersection between a continuous geometry (reference geometry) and a finite element mesh (intersection grid). Therefore, the parts needed to define the voxelization process are as follows:

• Reference Geometry: A CAD file of the part to be printed, which will be transformed to a voxel-based geometry; and

• Intersection Grid: A hexahedral element-based mesh that bounds the reference geometry. The elements are oriented with respect to the printing direction to mimic progressive layer deposition and facilitate both progressive element activation and heating of elements from the moving heat flux.

From the description of the process, it is apparent that the quality of geometrical approximation is dependent on dimensions of the intersection grid. However, a higher resolution intersection grid warrants more computational resources, so a balance must be found such that important details are accurately captured, and the simulation is feasible to be ran in reasonable time according to available computational power.

Despite further comparisons between conventional meshing with tetrahedral elements and voxel-based meshing with hexahedral elements being necessary, key features of both methods hold the following items apparent:

• The voxel-based mesh showcases a uniform meshing grid that is built based on the progressive layer deposition characteristic of a powder bed AM process. Since the voxel-based mesh is constructed along the printing direction, each element is activated at a similar time and rate compared to its neighboring elements in the layer plane, potentially improving convergence.

• The conventional mesh does a better job at representing the surface of the part, whereas the voxel-based one exhibits a jagged representation of curved surfaces, which can lead to stress concentration spots.

• The transformation to an adequately sized voxel-based mesh means less complexity and computational resources consumed.

• The likelihood of bad elements being present is reduced with the use of a voxel-based mesh, in spite of the geometrical complexity of the original part.

With the aim of achieving a generalized workflow applicable to any machine parameters and complex geometry, it is first studied how to take advantage of Autodesk Netfabb’s geometry slicing and support structures creation capabilities as well as the AM process setup, and import this information into FEM simulation, namely, in Abaqus.

The slicing process depends on the geometry itself, part orientation, build direction, and printing layer size, which are parameters Autodesk Netfabb allow to adjust. As these are set, the software proceeds to cut the tessellated triangles such that each layer is made up of closed polygons that constitute the cut surface. These closed perimeters are defined as walls or contours, whereas the inside region is the infill. The hatch filling geometry is defined through the set of parameters made available in Autodesk Netfabb and is automatically added to each layer.

Finally, with the entirety of the scan path defined geometrically, the speed of the laser and its power remain to be set. The parameters can be defined separately for both the contour and the infill, allowing to define the position of the laser in both space and time. The output of the slicing process is chosen to be the laser vector LSR file format.

As a first approach and in order to achieve a computationally efficient simulation, a lower resolution simulation is proposed using a coarser mesh (element size larger than a unitary printing layer thickness), a larger time increment, and a concentrated moving heat source model. The larger time increment results in the lumping of the printing events, such that the near-field action is not accurately captured, but instead its effect is averaged, while the far-field region is still described.

The annealing effect is a feature allowing to simulate stress relaxation and creep processes that occur typically at temperatures higher than a so-called relaxation temperature $T r e l a x$ to be determined. In this case, the stress and strain components are set to zero if the element temperature exceeds $T r e l a x$⁠. This is important because elements are activated at high temperatures. In fact, as the time step of the low-resolution simulation is not enough to capture high-temperature peaks at the melting spot, an initial temperature is applied to each newly activated element, such that it defines initial thermal contraction during structural simulation, and so a proper constitutive description at high temperature is needed. As such, the value chosen for the present work for Ti-6Al-4V, $T r e l a x=690 °$C, is taken from literature.²⁷ 

The uncoupled thermal-mechanical simulation is started by the definition of the thermal history of the process. The approach taken in this work focuses on a transient thermal analysis, with a duration longer than the one of the printing of the part itself, to allow for part cooling.

Progressive cooling is present throughout the entire building process, namely, through both convection and radiation effects. Heat loss due to radiation is defined by the ambient temperature, $T a m b$⁠, inside the printing chamber as well as the emissivity of the exposed surface of the part, $ε$⁠. On the other hand, heat loss through convection also takes places, with the key-parameter heat transfer or film coefficient, $h$⁠. It is defined according to the inert gas atmosphere present in the printing chamber (typically the inert gas Argon).

Partial element activation adds heat content to the system and also plays a role on cooling effects due to both convection and radiation, as it allows to change the exposed surfaces of the elements cut by the activation plane.

Both initial and boundary conditions are essential to accurately represent temperature evolution of the part during printing. In order to define substrate process heating, a fixed temperature boundary condition is applied to the bottom of the build plate part $T s u b$⁠. Separately, the powder bed preheating effect is defined by a predefined field applied to the entire build plate in the initial step, named $T p r e$⁠.

The initial temperature of the deposited powder from each printing layer, $T p a r t$⁠, is also set through a predefined field applied to the printing part.

The mechanical field is driven by the results of thermal simulation. The same laser-path and recoater roller event series INP files are needed, in order to properly activate the elements, and the temperature field information stored in the thermal simulation ODB file is also required. Similarly, for numerical simplicity, the same mesh geometry is shared between thermal and structural simulations, with the sole difference being the transformation of heat transfer elements into general 3D solid elements.

Under low-resolution assumption, there is a need to define the initial temperature of elements as they are activated, as the relaxation temperature, $T r e l a x$⁠, during the structural analysis.

To prevent rigid body motion, the displacement constraints are applied to the building plate. The building part itself is then defined through the tie constraint connection with the building plate. Finally, the results of the thermal simulation are introduced through a predefined field, which makes use of the ODB file generated by the thermal analysis.

For the topology optimization procedure, the load application points are connected to each of the six respective mounting bosses through structural distributing coupling constraints, such that relative motion between each is possible. Identically, the centers of mass are connected to the same mounting bosses through a continuum distributing coupling constraint (see Fig. 4).

The bolt holes at the bottom region are used to apply boundary conditions to the component. The bolt hole displacement constrained in x and y directions [Fig. 5(a)], the top and underside bolt hole regions displacement constrained in the z direction [Fig. 5(b)], and the fillet feature on the conic element displacement constrained in x and y directions [Fig. 5(c)].

The load cases considered are purely related to proof stress loading, such that Ti6Al4V ELI titanium alloy linear-elastic properties can be used throughout all topology optimization simulations.

Leonardo Helicopters provided both the design envelope and fixed elements’ CAD files. The design envelope defines the volume subject to optimization and changes, whereas fixed elements represent mounting regions to other helicopter components and parts and shall remain unchanged. However, Leonardo accepted changes that were carried out to the fixed elements’ region to mend broken or missing regions (since the provided file was obtained directly from the current component that has bolts connecting the two parts), to include the central hydraulic region (as its optimization is beyond the scope of this work and, therefore, remains identical), and to simplify the geometry to reduce the need for a fine mesh [see Fig. 6(a), where red regions highlight the removed material and green regions highlights the added material]. Such changes are thereafter used to update the design envelope [see Fig. 6(b)]. The assembly of the parts is meshed with quadratic tetrahedron elements (C3D10) and an element size of 7 mm.

A volume minimization approach is considered with displacement and frequency constraints, utilizing both SIMP and MIMP interpolation methods. The stiffness requirements are approximated to displacement limits based on the loading at each point (axial stiffness and normal stiffness) or pair of points (relative stiffness). Several trial cases are defined to better understand the implementation of topology optimization to the part geometry and the impact of the parameters on the final design solution. More precisely, SIMP and MIMP methods are compared in trial cases 1 to 4, whereas trial case 5 marks the introduction of numerous relative and normal stiffness constraints. Trial case 1 [Fig. 7(a)] is defined by a volume minimization topology optimization task, with axial displacement constraints, utilizing SIMP, whereas trial case 2 [Fig. 7(b)] utilizes MIMP. Both solutions are identical, since no mass-related design response is implemented, with trial case 1 showcasing a lower number of gray elements (with relative density between 20% and 80%).

Trial case 3 [Fig. 8(a)] is the same as trial case 1, with the added inclusion of the minimum frequency constraint. Similarly, the trial case 4 [Fig. 8(b)] is the same as trial case 2, with the added inclusion of the minimum frequency constraint. The introduction of the frequency design response drives both trial cases to different solutions. In fact, trial case 3 fails to achieve convergence, while trial case 4 converges to a well-connected design that also respects the frequency constraint.

A smoothing process is done by Tosca Structure.smooth, so the STL resulting file can be interpreted in Altair Inspire using polyNURBs and turned into a solid part to be validated in Abaqus under the same loading and boundary conditions. A summary of the relative mass and frequency can be seen in Table I, revealing the inability of the SIMP method to comply with the frequency requirement. All interpreted designs have a mass lower than the current geometry. Additionally, all trial cases respect axial stiffness requirements.

An additional trial case is introduced with the inclusion of relative and normal stiffness constraints in the form of displacement design responses. Since the trial case 4 is the only one to satisfy the frequency constraint, it is chosen to be the trial case to expand on with the new constraint. In fact, the newly obtained design is vastly different, with a frame truss-like structure connecting the outermost mounting bosses in each face. The interpreted final design can be seen meshed in Fig. 9. It respects both the axial stiffness and frequency constraints, while having a vastly superior performance with regard to the relative and normal stiffness constraints when compared to the trial case 4 design. It also has a similar mass at 0.9771 of the current design. Therefore, it is the chosen topology optimized design.

One of the main aims of the present study is to minimize and homogenize stress hotspots that occur during the static analysis. Additionally, the proof stress loading case is considered, and not the ultimate stress, as the behavior of the component above fracture failure is not accurately modeled with the type of elements that were considered. The shape optimization task utilizes the controller-based optimality criteria to perform Von Mises stress maximum minimization, considering the same load cases as in the topology optimization task. The optimization constraint is defined as a constant volume. The results of the shape optimization process can be seen in Fig. 10, where the surface displacement vector is represented. Regions in red expand outwardly, and blue regions inward.

Based on said surface vector results, the polyNURBs file referring to trial case 4 is modified. The first iteration results in a 1% mass decrease and smoothing out of the Von Mises stress, as seen in Figs. 11 referring to the red path highlighted in Fig. 9.

An additional iteration is performed to further reduce the weight of the component, while maintaining the stress homogeneity and ensuring it complies with the design requirements. In fact, a final mass reduction of 5.12% is achieved compared to the current component design, while respecting the proposed design requirements.

The resulting geometry from the applied structural optimization process is seen in Figs. 12 and 13, where the stress and displacement magnitude fields are shown, respectively, and compared to the current component performance counterpart. The results are obtained for the proof stress loading case under bending. It should also be noted that the scale utilized is entirely identical between the two components (omitted).

A low-resolution AM process simulation is implemented on a scaled version of the part to verify the printability of the proposed design, following the implementation of a benchmarking simulation that validates the chosen method, utilizing a model defined by the National Institute of Standards and Technology (NIST) (see  Appendix). Additionally, for consistency with said simulation, the same temperature field boundary conditions are applied. The substrate is preheated to a temperature of 80  $°$C, and the newly deposited layer is at a temperature of 40  $°$C, the same as the surrounding environment temperature. The first step involves the definition of support structures in Autodesk Netfabb utilizing its automated creation tools. These depend on the printing process itself. A regular hexagonal bar support was chosen, and it is shown in Fig. 14(a).

The voxelization process, in addition to the surface mesh describing the part to voxelize, requires the creation of the intersection grid, a hexahedral element-based mesh “box” that encapsulates the part. A grid with element size 0.4 mm is defined, such that each element contains, approximately, seven physical layers. Finally, the part intersection translates to the definition of a volume fraction in grid elements. The ones with a value higher than 0 (intersection took place) are extracted and translated into an INP file. The INP file characterizes the node and element sets that allow for the definition of a solid mesh in Abaqus [see Fig. 14(b)].

The temperature field drives the mechanical analysis in the proposed AM simulation. The snapshots from the transient thermal simulation represent the moment the new layer deposition occurs (described in Subsec. IV E). Examples of these can be seen in Fig. 15(a) (representing the first physical layer), Fig. 15(b) (representing a halfway through the physical layer), Fig. 15(c) (representing the last physical layer), and Fig. 15(d) (representing the moment after the postprinting cooling interval). The peak temperatures are not as high as expected from a fine time increment simulation. However, the top part (last printed layer) converges to approximately 40  $°$C, the temperature of the surrounding environment, whereas the bottom part (first printed layer) converges to approximately 80  $°$C, the temperature of the substrate.

The temperature plot with respect to time can be seen in Fig. 16, showcasing the nodal temperature (NT11) of chosen nodes, calculated as moving averages (Eulerian reference frame), from the first printed layer to the last.

The results of the mechanical simulation reveal high distortion in a region with overhanging material and no support structure, as it is expected (see Fig. 17).

The representation of residual stress is done with Von Mises stress distribution (depicted in Fig. 18). However, the low-resolution nature of the simulation has its drawbacks revealed in stress distribution, namely, the presence of stress concentration points at the surface of the part, as a result of the voxelization process, that causes sharp corners at the part’s surface. The region is, as expected, subjected to larger temperature gradients due to the close contact to cooler elements, such as the inert gas. While displacement is less susceptible to mesh changes, the stress field showcases larger sensitivity, leading to observed stress peaks at single nodes. Nevertheless, the stress distribution is, for the most part, fairly homogeneous.

It has been shown that the MIMP interpolation method, while based on the SIMP method, achieves both convergence and better results for the given design constraints. In fact, a mass-dependent analysis benefits from the use of MIMP, which is in fact verifiable. Comparing the results from trial cases 1 and 2, both volume minimization problems with displacement constraints, they exhibit the same geometry result, except for the increase in gray elements in trial case 2 (from 1.32% with SIMP, to 3.97 with MIMP). Considering that the design interpretation is at the discretion of the user, the proposed solution being similar at a grand scale, it means that both can be designed the same way. However, the introduction of the frequency constraint in trial cases 3 and 4 highlights the differences between both interpolation methods. As SIMP fails to achieve convergence and the frequency constraint oscillates throughout the optimization process, MIMP converges with fair swiftness. By comparing all results, MIMP is chosen to be the interpolation method to be used for remaining trial cases.

The introduction of further design constraints in the optimization task results in definite geometrical changes, namely, in trial case 5, which features a frame surrounding each face, in order to comply with relative displacement constraints between mounting bolts. Considering this and the final mass savings of 2.29% of the original component, trial case 4 is deemed the best design and the one to undergo further analysis. Shape optimization compliments the design process by addressing another set of design requirements, namely, the ones pertaining to static strength. The final design represents a 5.12% mass decrease compared to the current component. Such reduction is largely attributed to DfAM methodology employed, as it allows for a single component design, not limited to typical manufacturing constraints that do not allow for organic-looking designs as the one presented.

The obtained thermal field in the aerospace part follows the predicted behavior based on the two driving effects, these being the fact that as the printed part height increases, the influence of heat conduction to the substrate decreases, leading to a tendency to increase the temperature, seen between layers 1, 167, and 334, and the increasing effect of convection and radiation cooling effects. The superposition of the two effects leads to the initial increase in the temperature per layer, followed by a plateau and eventual decrease. Moreover, the voxelization process with such a complex geometry implies initialization of layers with vastly different sections and number of elements one after the other, leading to discontinuities present in Fig. 16, as the heat stored for each step can be vastly different.

The results obtained from AM simulation reveal distortion that is in fact higher in regions with an unsupported overhanging material, as expected, reaching a maximum of 0.6277 mm.

Combination of topology and shape optimization with finite element modeling can be a powerful tool in the design of aerospace components produced by laser-based additive manufacturing. This can be achieved, namely, by combining mass reduction and stress/distortion minimization capabilities of optimization algorithms, with the predictive power of the finite element method to describe thermal and residual stress fields, in what concerns the production of such components by means of laser-based additive manufacturing.

This paper is in memory of Professor Jyoti Mazumder, a pioneer in Laser Deposition Methods and Additive Manufacturing of Metals, as well as countless contributions to the fields of Laser Materials Processing and Laser Based Diagnostics. The project leading to this application has received funding from the Clean Sky 2 Joint Undertaking (JU) under Grant Agreement No. 717194 (AMATHO). The JU receives support from the European Union’s Horizon 2020 research and innovation programme and the Clean Sky 2 JU members other than the Union. A. M. Grande and G. Sala are grateful to S. Sartori and F. Montagna (Leonardo Helicopters) for the provided design data and for helpful discussions. J. M. Guedes was partially supported by LAETA Project No. UIDB/50022/2020. A. M. Deus was partially supported by CeFEMA Project No. UIDB/04540/2020. Both projects are sponsored by FCT, the Portuguese Agency for Science and Technology.

The authors have no conflicts to disclose.

R. J. Ferreira: Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Visualization (equal); Writing – original draft (equal). A. M. Grande: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). J. M. Guedes: Methodology (equal); Project administration (equal); Supervision (equal); Writing – review & editing (equal). A. M. Deus: Methodology (equal); Supervision (equal); Validation (equal); Writing – review & editing (equal). G. Sala: Funding acquisition (equal); Methodology (equal); Project administration (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from Leonardo Helicopters. Restrictions apply to the availability of these data, which were used under license for this study. Data are available from the authors upon reasonable request and with the permission of Leonardo Helicopters.

The AM simulation is implemented in Abaqus by addressing AMB2018-01 additive manufacturing benchmark problem defined by NIST,²⁸ such that the methodology implemented here is validated, to be later used in the part featured in the present study (see Sec. VII). Built-in user subroutines and mesh intersection tools are utilized to define both the moving heat flux /laser source and progressive material deposition performed by the recoater roller in powder bed-type AM processes.

The benchmarking simulation geometry is made available as an STL file in the NIST website,²⁸ as seen in Fig. 19.

The residual elastic strain and distortion results are compared against the x-ray diffraction measurements provided by NIST (2018).²⁸ Namely, the z and x direction residual elastic strains are obtained at the $y=2.5$ mm middle cross section of the printed part. Figures 20(a) and 21(a) present the residual elastic strain in x and z-directions, respectively. The residual elastic strain in x and z-directions provided by NIST (2018)²⁸ can be seen in Figs. 20(b) and 21(b), respectively.

Additionally, elastic strain in z-direction path data at $z=2.75$ mm [Fig. 22(a)] and $z=10.75$ mm [Fig. 22(b)], both at $y=2.5$ mm, are compared against the NIST benchmark test data. The obtained results show a good correlation with the tests at both build heights.

Finally, distortion is measured at eleven top ridges found on the part at $z=12.5$ mm in the central plane at $y=2.5$ mm. The plot in Fig. 23 shows both benchmark test results and simulation results obtained, again suggesting a reasonably accurate prediction of the part behavior under AM printing conditions.