Issue 
Mechanics & Industry
Volume 21, Number 3, 2020



Article Number  304  
Number of page(s)  17  
DOI  https://doi.org/10.1051/meca/2020011  
Published online  07 April 2020 
Regular Article
Stressbased topology optimization of compliant mechanisms using nonlinear mechanics
^{1}
DMSM, ISAESUPAERO,
10 Avenue Edouard Belin,
Toulouse,
France
^{2}
Univ Toulouse, ISAE SupaeroINSAMines AlbiUPS, CNRS UMR5312, Institut Clément Ader,
Toulouse, France
^{3}
Airbus Operations SAS,
316 Route de Bayonn,
31300
Toulouse, France
^{*} email: gabriele.capasso@student.isaesupaero.fr
Received:
12
July
2019
Accepted:
29
January
2020
The present work demonstrates how a light structure can be easily designed through Topology Optimization even including complex analysis and sizing criteria such as hyperelastic NeoHookean materials for nonlinear analysis and aggregated stress constraints. The SIMP approach was adopted and two different strategies were analysed using an in house versatile MATLAB code. MMA was used as reference optimizer (in structural optimization) whereas a unified aggregation and relaxation method was adopted to deal with stress constraints. Feasibility was analyzed from the viewpoint of allowable stress verification. Two test cases are then studied: a morphing airfoil (for aeronautical applications) and a geometric inverter (for mechanics and biomedical applications). For both, a hyperelastic NeoHookean material was chosen. Finally a complementary study on the effects of constraints and the input force intensity is also presented.
Key words: Topology optimization / Nonlinear mechanics / Stressbased optimization / Morphing / Compliant mechanism
© G. Capasso et al., published by EDP Sciences 2020
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
The past twenty years have seen increasingly rapid advances in the field of structural optimization: in particular, Topology Optimization (TO) is becoming a key figure in the research panorama. This method leads to an optimal structural efficiency, through the removal of unnecessarily placed material, in order to accomplish the minimization of a given objective (such as the compliance of the system). Since the pioneering work of Bendsoe and Kikuchi (1988) [1] a number of approaches have been proposed, among which are: densitybased methods (the wellknown SIMP algorithm) [2,3], evolutionary strategies [4] and levelset based methods [5,6].
Moreover, several researchers have tried to implement and diffuse simple and efficient codes. One of the most wellknown and important results is due to Andreassen et al. [7], who created a very simple and efficient 88line code written in MATLAB, capable of optimizing an MBB beam through a SIMP algorithm. Compared to other available codes and the same 99line code (published by Sigmund et al. [8]), the former needs much less computational power, thanks to judicious use of vectorization and preallocation of most of the variables. The present work focuses on the SIMP approach and the code developed in our study is largely inspired by the top88.
In the literature, several promising applications of TO to classical aerospace structures can be found [9–13]. Recently, researchers have shown an increased interest in alternative architectures presenting wings with a variable geometry, better known as “morphing wings” [14–17]. Different types of morphing have been tested both numerically [18] and experimentally [19], revealing the advantages in terms of aerodynamic drag, noise and vibration control with respect to conventional aircraft configuration, based on flight control surfaces. Two main categories can be distinguished [20]: planform morphing, based on the modification of the plant itself of the wing, and performance morphing, which leads to variating the airfoil in shape (shape morphing) or in incidence (twist morphing). In this paper, a shape morphing design is analyzed.
Earlier work on TO was restricted to linear structural designs. Recent research explored the influence of nonlinear mechanics on the TO in several simple applications, including geometrical nonlinearity [21–24], material nonlinearity [25–30] and both geometrical and material nonlinearities simultaneously [31–33]. Work from [24] evidences the differences encountered in a TO problem, according to the choice of a linear or a nonlinear analysis, leading to totally different results.
The introduction of Nonlinear Mechanics into the TO of a morphing wing has already been treated by Bhattacharyya et al. [9]: in particular, geometric and material nonlinearities have been investigated. In [9], the optimization problem presents the objective of minimizing the actuating force, under constraints of desired volume fraction, compliance and trailing edge displacement. The obtained structure is, in fact, a compliant mechanism capable of getting a geometric advantage. The first of the two parts of the present work focuses on the same structure, whilst the objective is to design the lightest structure under stress constraints, in order to conceive a feasible structure.
In recent years other compliant mechanisms have frequently been proposed for aerospace systems due to their advantages over traditional rigidlink mechanisms. As compliant mechanisms can be manufactured using fewer parts (usually one single part), they are lighter than their traditional counterparts. They are also easy to manufacture with modern 3D printing techniques, and do not require assembly. The ability to 3D print mechanisms and tools on board the International Space Station (or future bases on the Moon and Mars) will save substantially on the cost and time involved in transporting these devices from Earth. To be prepared for deepspace exploration and habitation in the future, research and technology development concerning the utilization of additively construct structures from insitu materials is incorporated into the inspace manufacturing technology roadmap [34].
A compliant mechanism has been investigated by ConlanSmith and James [35]. Their study focuses on the optimal design via Topology Optimization of a bioinspired structure. The second part of the present work is inspired from [35], and aims at solving the same optimization problem by considering feasibility constraints, in terms of admissible stresses.
Stressbased Topology Optimization is an active research field. Several examples available in literature show the challenge posed by stress constraint in TO problems. Essentially, two problems arise [36]: first, a large number of constraints must be considered, since unlike stiffness, stress is a local quantity; second, stress is highly nonlinear with respect to design variables. Recent advances show possibilities to overcome these difficulties [37–42]. There are few recent examples in the literature where the stress distribution is considered as a constraint in TO problems with geometrical nonlinearities (as for instance [39]). The introduction of both stress constraint and nonlinear mechanics into a TO using nonlinear mechanics is the main challenge of the present work.
The present work aims at designing two compliant mechanisms: a morphing wing and a geometric inverter. The former is treated as an academic example conceived to validate the method, while the latter is a more realistic testcase. In both mentioned parts of this article, the optimal design is achieved using a nonlinear analysis of a hyperelastic material through TO, considering the stress as a constraint of the problem. Section 2 illustrates the theoretical background behind the optimization problem. Section 3 details the mathematical formulation and the newly introduced elements of the paper. Section 4 presents the different testcases, outlining the relative mathematical formalization. Optimized designs and numerical results are finally reported and discussed in Section 5.
2 Theoretical background
In this section the basic principles of nonlinear mechanics and Topology Optimization are explored.
2.1 Nonlinear mechanics
The present work focuses on the design of largedisplacement mechanisms, in which material doesn’t follow traditional linear mechanics. This topic will be just introduced in this section, although it is largely explained in [43].
Mechanical nonlinearities are distinguishable into four main categories [43]: (1) geometric nonlinearity (nonlinear straindisplacement relation), (2) material nonlinearity (nonlinear constitutive relation), (3) kinematic nonlinearity (nonconstant displacement BCs, contact) and (4) force nonlinearity (followup loads).
The resolution of a Nonlinear Finite Element Analysis (NFEA) is not straightforward, since the global system is nonlinear: in [43] several methods are illustrated. The most common algorithm is the one which goes under the name of NewtonRaphson, which is usually modified through the introduction of a stepforce control. This method guarantees better stability than the simple NewtonRaphson, preserving the quadratic convergence [43].
2.2 Topology Optimization (TO)
Topology Optimization methodology is a classical way to obtain an optimal structure by removing material where necessary. Currently, the four most common methods used to solve the optimization problem are:

evolutionary strategies [4];
The former is the main subject of the present article and will be explained in the following section.
2.2.1 SIMP algorithm
The densitybased approach, also known as Solid Isotropic Material with Penalization (SIMP), is the most popular TO method. SIMP has a solid mathematical foundation [1,47–50]. This method is capable of handling various objectives and constraints, and is relatively easy to implement within a finite element environment.
The basic idea behind NFEA based SIMP is that each finite element is associated with a fictitious pseudodensity variable 0 ≤ ρ ≤ 1, that essentially parameterizes the topology. The pseudodensities are optimized to reach the desired objective. The general algorithm is detailed in Algorithm 1.
In [51] the reader can find an exhaustive review on numerical instabilities typical of TO. They include essentially checkerboard patterns, meshdependencies and local minima. To deal with the first two issues, a meshindependency filter can beapplied [51]. A fundamental aspect of such filters is the distinction between the physical densities ρ and the merely numerical ones x: the latter category constitutes the design variable array updated in the optimization process at each iteration; once having updated numerical densities, these are filtered to obtain the physical distribution of the material. NFEA is based on physical densities ρ.
SIMP Algorithm
2.2.2 Optimizer: method of moving asympthotes
When the number of design variables is very large, gradient based optimization methods are the most efficient algorithms. They are preferred to gradientfree optimization methods since the latter typically require many evaluations.
As shown in [9], the MMA is one of the most suitable optimizers which could substitute the Optimality Criteria adopted in the 88line code by Sigmund [52], in order to accelerate convergence.
This method is based on a local convex approximation of objective and constraint function: this needs as input the local evaluation of the objective function, the constraints and all their derivatives with respect to the numerical densities x. The solution is said to have reached convergence when the first order Karush Kuhn Tucker (KKT) conditions are satisfied to an absolute tolerance of 2 × 10^{−3}.
3 New formulation
3.1 Material description
A nonlinear hyperelastic material was investigated: the same methods described in [9] were adopted, but modifying the proposed algorithm in order to achieve the same results as well as a fast convergence. Moreover, far too little attention has been paid to the feasibility of the structure itself: this can be observed by the absence of a constraint related to the maximum stress supportable by the material.
The chosen material is described by the potential density defined in [9,53] and reported in the formula below: (1)
Here J is the determinant of the deformation gradient F; the term tr(C) denotes the trace of the right CauchyGreen deformation tensor C = F^{T}F; finally,
constitute the Lamé constants of the material which are seen as functions of Poisson ratio ν and Young Modulus E. This latter is supposed to be related to the density ρ of each element (e), through one of the two following relations: (2) (3)
which will be used in order to make comparisons in our forthcoming NFEA. Following [7], equation (2) is referred to as calssical SIMP and equation (3) as modified SIMP.
3.2 Nonlinear finite elements analysis
We treat the NFEA as another minimization subproblem in the global optimization: the function to minimize is the total energy of the system, given by both internal and external forces, while the constraints are the displacements imposed in the boundary conditions. The domain is decomposed using triangular elements.
The definition of the objective function derives from the Principle of Virtual Works, where the external work is computed as the product between the force and the displacement of the point of application: (4)
Here, U represents the vector of the nodal displacements which has to be determined by Variational Principle. In order to make the Optimization Toolbox in MATLAB more efficient, we have to provide it with the gradient and the Hessian of the total energy. The gradient will be the residual R, computed as the difference between internal and external forces: (5)
The Hessian will be simply the tangent stiffness matrix K_{T}, computed as the second derivative of the compliance to the nodal displacements: (7)
The equilibrium condition is satisfied when: (8)
The MATLAB tool fmincon is a suitable instrument to perform NFEA, since it treats this subproblem as a constrained minimization problem. Moreover, it contains a stable implementation of the NewtonRaphson Method in the InteriorPoint Region [54]. This feature turns the fmincon process into a fast and reliable solver for nonlinear mechanical problems.
3.3 Stress constraint treatment
The Second PiolaKirchoff stress S relates forces in the reference configuration with areas also measured in the reference configuration (see [53]). It can be calculated as twice the derivative of the potential energy function with respect to the right CauchyGreen deformation tensor: (9)
The Cauchy stress, σ, which relates the forces in the deformed configuration with areas in the deformed configuration, can then be defined as: (10)
As a consequence, the microscopic (local) stress tensor will be evaluated as: (11)
where B is the left CauchyGreen tensor, which is defined as B = FF^{T}, and I is the secondorder identity tensor. Here, microscopic stress was considered by calculating the VonMises stress, under the hypothesis of plain strain. It is important to underline the fact that this does not depend directly on density distribution. In this work the unified approach proposed by Verbart et al. [40] was adopted to incorporate stress constraints in topology optimization. Assuming that the densities in SIMP represent a porous microstructure, one can distinguish the stress at macroscopic (effective) versus microscopic levels [55]. The former is the stress computed using the Young modulus of the SIMP model (c.f. Eqs. (2)–(3). The latter is computed considering a stress model that mimics the behavior of the “local stress” in a rank2 layered composite [55]. For this reason we identify the microscopic stress as: (12)
where S_{0} is the second PiolaKirchoff stress considering a material with full Young Modulus E = E_{0} (i.e. ρ = 1). To access singular optima classic of mathematical program with vanishing constraints, a relaxation method is adopted [40]. In particular, the relaxed constraint is defined as: (13)
where . In classical SIMP formulation (see Eq. (2)), we have ρ_{min} = 0.1 and ρ_{max} = 1 while in modified SIMP (see Eq. (3)) ρ_{min} = 0 and ρ_{max} = 1.
3.4 Adjoint sensitivity analysis
The method which is going to be illustrated has been proposed by Bhattacharyya [9]. Since the objective and constraint function depend both on density distribution and nodal displacements, an adjoint sensitivity analysis is to be performed. In particular, we use an augmented Lagrangian L[ρ, U_{F}(ρ)], which may be defined as follows: (14)
Here f indicates a general function, indifferently the objective or the constraint; λ is a Lagrange multiplier ensuring that the residual of forces R vanishes as required in equation (8). Deriving this expression and using the chain rule, we obtain: (15)
Collecting all the implicit terms indicated by , we obtain: (16)
In order to eliminate the implicit dependence of free nodal displacements on density distribution, the lagrangian multiplier can be computed, remembering the definition of the derivative of the residual related to the nodal displacements, as follows: (17)
One can observe that the evaluation of gradients only requires a linear system of equation resolution, extremely convenient if compared with direct or finite difference approaches.
3.5 Density filtering
As pointed out in Section 2.2, a meshindependency filter was applied to deal with certain numerical instabilities, namely checkerboard patterns and meshdependent solutions [51].
Here the method suggested by Bruns and Tortorelli [56] was implemented, where: (19)
Here w_{ij} is the weight associated with the ith element within the prescribed neighborhood of the jth one, determined as: (20)
The parameter R_{min}, called filter radius, is the radius of thespecified neighborhood and r_{ij} the distancebetween the centroids of the two elements. This formulation proposed in [56] is valid under the assumption of uniform mesh. This filter makes the results independent from the adopted mesh, improving the resolution as R_{min} decreases. However, the maximal resolution is linked to the mesh itself, since the minimum applicable value R_{min} cannot be smaller than the size of the largest element.
When computing derivatives with respect to the numerical densities x (necessary to gradientbased optimization methods), two steps are followed: firstly, derivatives with respect to physical densities ρ have to be computed; then, chain rule is applied to obtain the desired gradient.
3.6 Constraint aggregation
Treating the relaxed stress distribution, there would be an excess of constraints. So the maximum may be taken into account to represent the whole domain. Given the fact that the max function is nonderivable, the lower bound KreisselmeierSteinhauser function [36,57] was employed: (21)
The larger the factor P, the better the approximation of the max, but also the more expensive the computational cost of the overall optimization process. In fact the nonlinearity of the optimization problem is enhanced by the use of higher values of P. Therefore, P has to be chosen as an appropriate compromise between computational burden and accuracy in stress control.
4 Applications
4.1 First optimization problem: NACA airfoil
As an academic example, this part focuses on the minimization of the weight of a shape morphing wing.
Using the classical SIMP approach (see Eq. (2)), the whole problem can be formalized as follows: (22)
Here V represents the volume fraction, U_{e} the displacement at the trailing edge and σ_{lim} is the highest tolerable stress in the structure (to avoid failure). The constraint relative to the σ in each element is due to the desire to design a feasible structure while not reaching the weakest one: in particular, we will consider the VonMises stress, under the hypothesis of plain strain. Finally, the constraint on U_{e} derives fromthe necessity of designing an efficient shape morphing wing.
We shall also perform further analyses based on the modified SIMP approach (see Eq. (3)), for comparisons. The new optimization problem is then: (23)
Here another constraint has to be added, since the volume fraction minimization under stress constraint, would generate gray regions, preventing the creation of a continuous chain between force application point and trailing edge. Such constraint will be an inferior limit on the compliance C of the trailing edge region [35]: this value is arbitrarily defined as the initial compliance.
Fig. 1 Initialdomain: airfoil optimization problem. 
4.2 Second optimization problem: geometric inverter
The second optimization problem is a more realistic test case, inspired by the work in [35]. Figure 2 illustrates the initial domain.
The idea is to design a compliant mechanism which would invert the input displacement, creating as output another displacement on the opposite side of the structure. The actual input is a force. The optimized mechanism is defined as the one capable of reaching the highest ratio between output and input; as the two displacements have opposite signs, this is translated into a minimization problem: (24)
The second pseudodensity parameterization of the SIMP Algorithm is adopted. Before the main loop, the topology is initialized by imposing ρ = V_{0} in every element of the domain. The constraints are quite similar to the ones adopted in the previous section. We have thus:

volume fraction: the structure has to be as light as possible, to reduce costs and weight;

stress: the whole mechanism has to support the applied load without breaking, or reaching the yield point;

compliance: a continuous chain between input and output is ensured, by imposing that the output region has to keep the same compliance at the first iteration loop (where ρ = V_{0} in every finite element) [35].
Another aspect which appeared fundamental is the symmetric nature of the domain: to accelerate computation, only half of the whole domain has been simulated.
Fig. 2 Geometric inverter scheme. 
5 Results and discussion
In this section, numerical results are presented and discussed. Simulations have been performed using a unique MATLAB code, adopting the original MMA code provided by Svanberg [52] and without altering its default parameters. Given the fact that the constraint is approximated (as explained in Sects. 3.3 and 3.6) stress distribution postprocessing on the optimal structure is necessary. The advantage of our methodology is that stress distribution is already computed in our optimization loop (“Constraints evaluation” in Algorithm 1.
5.1 NACA airfoil
A first analysis is performed on the two problems 22 and 23, considering the same domain (represented in Fig. 1). The parameters adopted in this section are reported in Table 1.
Default parameters used in the Airfoil optimization problem.
5.1.1 First formulation results
In this subsection, the results to the problem (22) are presented. The base structure can be observed in Figure 3a. This can be decomposed into four main sectors:

Front sector:This is responsible for stress redistribution, as illustrated in Figure 3d. The stresses use this path in order to create an equilibrating moment. The middle part, situated between the front chain and the application point of the external force, is empty: in fact, the volume fraction constrained minimization produces this result, where the stresses have to be collocated in an intermediate location. Moreover, this sector connects the two series of points where the homogeneous Boundary Conditions are applied.

Central chain:This sector is the transmitter of forces between the application point of the external force and the rear chain (discussed later), as shown in the macroscopic stress distribution (Fig. 3d). It is directly linked to all the other parts of the optimized structure.

Inferior chain:It is the main link between the central chain and the inferior BCs. In the base structure, It is also directly linked to the front part: this is due to the fact that it starts exactly at the left limit of the inferior boundary points affected by the application of homogeneous BCs.

Rear chain:This is the only part not interested in any redistribution of stresses, as shown in Figure 3c and d. Its main role is transmitting the commanded displacement at the trailing edge, thanks to an input force provided by the central chain. This sector affects the volume fraction negatively, given that it doesn’t have to carry any stress.
The optimized structure in deformed configuration is illustrated in Figure 3b. The reader can observe the logic at the base of the mechanism. By introducing a single input force at the node in green, the effort is transmitted to the rear chain through the central sector. Simultaneously, the inferior chain and the front sector connect the central chain to the nodes in blue (where homogeneous Dirichlet Boundary Conditions are applied). Finally, the rear chain converts the input force into a displacement at the trailing edge.
The distribution of the relaxed constraint proves the approximated nature of the G_{KS} function: this is not always guaranteed, but it is the product of a compromise between precision, stability and computational cost. However, the only elements where the relaxed constraints are not satisfied are, effectively, “empty”.
Fig. 3 Results relative to the Airfoil problem (first formulation). 
5.1.2 Second formulation results
When using the second formulation to parameterize the topology, a totally different structure is obtained. Results from problem (23) are presented in this subsection. The differences are evident by comparing Figure 3 and Figure 4.
As was visible in the previous numerical experiment, the aggregation method does not always guarantee the respect of the stress limits: this is a consequence of the G_{KS} approximation. This time, an overshoot can be found on the filled structure. A lower limit on the stress should have been introduced in order to adapt to the real constraint.
The logic behind the optimal structure presented in Figure 4 (solution of the problem (22)) is unaltered compared with to the first one in Figure 3 (solution of the problem (23)).
However, this formulation leads to a more efficient application of the stress constraint, since all empty elements do not contribute to the macroscopic stress distribution. An immediate consideration is the removal of most of the Front sector and the Rear chain from the previous results: globally the obtained volume fraction would be inferior.
As a consequence, the modified SIMP approach (see Eq. (3)) is more suitable for topology optimization problems presenting a stress constraint. This consideration is applied in the second part of the present work.
Fig. 4 Results relative to the Airfoil problem (second formulation). 
5.2 Geometric inverter
Out of the previous academic examples, the second formulation has been identified as a better compromise between severity of the penalizationon stress constraint and reality of the final structure, leading to a much lighter final structure. However, a lower limit on the stress has to be introduced.
In the second test case of the present article, modified SIMP was adopted. Moreover, only half of the domain has been considered in the simulation to accelerate computation (symmetric properties have been exploited). All the parameters characterizing the Inverter problem are reported in Table 2.
In Figure 5 the final structure is presented. In the deformed configuration (Fig. 5b), the gain in terms of displacement inversion is easily visible.
Each half of the compliant mechanism can be decomposed into four main sectors:

Central Bulb: It absorbs the applied force and transmits it to the rest of the structure.

Support Net: It guarantees an efficient redistribution of stresses all over the structure: the ramifications create the necessary moments to equilibrate the whole mechanism. In fact, as the available volume fraction constrained, the resultingpositioning produces the best compromise.

Angular chain: This connects the Support Net with the clamped point, guaranteeing the force equilibrium of the structure.

Terminal chain: It make the morphing nature of the whole mechanism possible, generating the displacement in output.
The relaxed stress constraint and macroscopic stress distributions present some differences: in fact, the former reports a stress constraint violation in the elements of the terminal chain of the actuator, where the displacements are higher, due to the aggregation;this violation disappears in the macroscopic stress distribution. On the contrary, the maximum effective stress is in the areas of higher curvature.
Default parameters used in the Inverter optimization problem (symmetrical problem).
Fig. 5 Results relative to the Inverter problem. 
5.2.1 Influence of stress constraint
Further numerical analyses were performed, in order to investigate the effects of the stress constraint. In particular, both increases and decreases on the stress limits were considered.
If the stress limit decreases to 0.2 MPa, the obtained structure is totally altered (results reported in Fig. 6). However, the global decomposition follows the same logic adopted in the previous section. Four main differences may be noticed:

the Central Bulb is more compact and does not present any hole;

the Support Net presents only two ramifications;

the Angular Chain and the Support Net cannot be distinguished;

the Terminal Chain is much longer than before.
The reader can easily understand the necessity of introducing an additional buckling constraint if the stress limit is exaggeratedly low.
An analogue analysis was performed for an augmented allowed limit stress of the material (results reported in Fig. 7). The differences are less evident than an analogue reduction on stress limits. However, the Support Net is more ramified, allowing for a different allocation of the efforts.
Another fundamental aspect is the imperfect exploitation of the material: in fact, the stress constraint isn’t active, leading to a structure where the maximum macroscopic stress isn’t reached (refer to Fig. 7d).
In Table 3 the numerical results are shown for each value of σ_{lim}. In particular, it has been demonstrated that higher limit in terms of allowed stress doesn’t necessarily imply an augmentation in final performances. Moreover, the number of iterations to convergence augments as σ_{lim} decreases.
Fig. 6 Results relative to the inverter problem (σ = 0.2 MPa). 
Fig. 7 Results relative to the inverter problem (σ = 0.4 MPa). 
5.2.2 Influence of the amplitude of the applied load
The applied load shows little differences with respect to the base structure. In particular, the Support Net is less ramified for higher loads: this reflects again the same effect of a very small reduction in acceptable limit stress.
Moreover, as the amplitude of the force increases, an evident overshoot on the macroscopic stress can be observed: this leads to the idea that another approximation for the max is needed.
In the Appendix, different designs of the optimized structure are presented.
In Table 4 the numerical results are presented for each value of F_{app}. By increasing the amplitude of the applied load, the geometric advantage improves: however, the gain isn’t proportional to the relative augmentation in terms of applied force. On the other hand, a general rule governing the number of iterations cannot be exploited.
6 Conclusion
Topology Optimization was used to design two different compliant mechanisms: a morphing wing and a geometric inverter. A continuous material distribution problem following the SIMP approach was formulated to find the optimal layout throughout a design domain. For our purpose, a hyperelastic Neohookean material was chosen.
We considered two TO problems. The first one consists in the minimization of the volume fraction of a NACA 0033 airfoil, which must be conceived as capable of reaching a fixed displacement at its trailing edge. The second TO problem aims at maximizing the geometric gain in a structure with given volume fraction, in order to design a geometric inverter starting from a rectangular initial domain. The stress distribution was considered as the main constraint in all the problems presented inthis paper. In particular, a relaxation and aggregation method was adopted to treat this constraint optimally.
The problem of morphing wing design was solved using two different models of the relation between Young Modulus E and density ρ of each element, namely classical SIMP and modified SIMP (see Eqs. (2)–(3)). The objective is always to minimize the volume fraction V, under constraints of performance (in terms of trailing edge displacement) and of feasibility (in terms of stress). Results showed a substantial gain in terms of V when we adopted the modified SIMP rather than classical SIMP: we pass from a volume fraction of 51.64% to 19.62%. The second structure is also easier to manufacture using a 3D printing technique or additive manufacturing, which are the natural outputs of Topology Optimization.
Modified SIMP approach was applied to parameterize the domain in the geometric inverter problem. The first results show a ratio of 318% between output and input displacement (respectively U_{out} and U_{in} in Fig. 2). Several analyses were performed, by modifying the stress limits or the applied load. One can observe that increasing the stress constraint does not imply a significant gain in terms of structure performances. By analyzing the effects of the applied load, the geometric gain increases as the force increases, but the relation is not linear.
Looking at the number of iterations, by reducing the allowed stress, convergence results to be slower; on the other hand, a general dependence between number of iterations and applied load cannot be proved.
Globally, the present work demonstrated how a light structure can be easily designed through TO even with hyperelastic material responses and stress constraints. Nonlinear stress constraints were treated trough aggregation based on the lower bound KreisselmeierSteinhauser function: effects of different aggregation functions could be further investigated.
The position of the input force was considered as an input of each optimization problem presented in this paper, but could be treated as a new design variable as well.
Moreover, the two optimization problems could be extended to 3D analysis, allowing also the introduction of more realistic Boundary Conditions.
Future works will focus on the application of stressbased TO on the design of a morphing winglet subject to very high deformations.
Looking forward, further attempts with different optimization problems could prove beneficial to literature on the subject.
Appendix A
Fig. A.1 Results relative to the inverter problem: applied force F_{app} = 0.15 N/mm. 
Fig. A.2 Results relative to the inverter problem: Variation of the applied force F_{app} =0.25 N/mm. 
Fig. A.4 Evolution of optimization constraints (expressed in negative null form) for problem (24): variations in σ_{lim}. 
Fig. A.6 Evolution of optimization constraints (expressed in negative null form) for problem (24): variations in F_{app}. 
References
 M.P. Bendsøe, N. Kikuchi, Generating optimal topologies in structural design using a homogenization method, Comput. Methods Appl. Mech. Eng. 71, 197–224 (1988) [Google Scholar]
 M.P. Bendsøe, O. Sigmund, Material interpolation schemes in topology optimization, Arch. Appl. Mech. 69, 635–654 (1999) [Google Scholar]
 M. Zhou, G. Rozvany, The COC algorithm, Part II: topological, geometrical and generalized shape optimization, Comp. Methods Appl. Mech. Eng. 89, 309–336 (1991) [Google Scholar]
 Y.M. Xie, G.P. Steven, A simple evolutionary procedure for structural optimization, Comput. Struct. 49, 885–896 (1993) [Google Scholar]
 G. Allaire, F. Jouve, A.M. Toader, Structural optimization using sensitivity analysis and a levelset method, J. Comput. Phys. 194, 363–393 (2004) [Google Scholar]
 M.Y. Wang, X. Wang, D. Guo, A level set method for structural topology optimization, Comput. Methods Appl. Mech. Eng. 192, 227–246 (2003) [Google Scholar]
 E. Andreassen, A. Clausen, M. Schevenels, B.S. Lazarov, O. Sigmund, Efficient topology optimization in MATLAB using 88 lines of code, Struct. Multidiscip. Optim. 43, 1–16 (2010) [Google Scholar]
 O. Sigmund, A 99 line topology optimization code written in Matlab, Struct. Multidiscipl. Optim. 21, 120–127 (2001) [Google Scholar]
 A. Bhattacharyya, C. ConlanSmith, K.A. James, Topology optimization of a bistable airfoil using nonlinear elasticity. In: 18th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference. American Institute of Aeronautics and Astronautics (2017) [Google Scholar]
 G. Capasso, R. Amargier, S. Coniglio, J. Morlier, C. Mabru, M. Di Sciuva, Structural optimization for propulsion aiframe, M.Sc. thesis, Politecnico di Torino, 2019 [Google Scholar]
 G. Capasso, S. Coniglio, M. Charlotte, J. Morlier, Optimisation topologique de structures adaptatives (bistables) en mécanique nonlinéaire, in 14éme Colloque National en Calcul des Structures (2019) [Google Scholar]
 S. Coniglio, C. Gogu, R. Amargier, J. Morlier, Pylon and engine mounts performance driven structural topology optimization, in World Congress of Structural and Multidisciplinary Optimisation. Springer (2017), pp. 1349–1363 [Google Scholar]
 J.H. Zhu, W.H. Zhang, L. Xia, Topology optimization in aircraft and aerospace structures design, Arch. Comput. Methods Eng. 23, 595–622 (2016) [Google Scholar]
 P. Bettini, A. Airoldi, G. Sala, L.D. Landro, M. Ruzzene, A. Spadoni, Composite chiral structures for morphing airfoils: numerical analyses and development of a manufacturing process. Compos. Part B: Eng. 41, 133–147 (2010) [Google Scholar]
 P.R. Budarapu, Y.B. Sudhir Sastry, R. Natarajan, Design concepts of an aircraft wing: composite and morphing airfoil with auxetic structures, Front. Struct. Civil Eng. 10, 394–408 (2016) [Google Scholar]
 K.C.W.Cheung, Digital cellular solids : reconfigurable composite materials. PhD thesis, Massachusetts Institute of Technology (2012) [Google Scholar]
 A. Airoldi, M. Crespi, G. Quaranti, G. Sala, Design of a morphing airfoil with composite chiral structure, J. Aircraft 49, 1008–1019 (2012) [Google Scholar]
 S. Joshi, Z. Tidwell, W. Crossley, S. Ramakrishnan, Comparison of morphing wing stategies based upon aircraftperformance impacts, in 45th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics Materials Conference. American Institute of Aeronautics and Astronautics (2004) [Google Scholar]
 J.M. Martinez, D. Scopelliti, C. Bil, R. Carrese, P. Marzocca, Cestino, G. Frulla, Design, analysis and experimental testing of a morphing wing, in 25th AIAA/AHS Adaptive Structures Conference. American Institute of Aeronautics and Astronautics (2017) [Google Scholar]
 D. Wagg, I.P.W. Bond, M. Friswell, Adaptive structures – engineering applications (2007) [Google Scholar]
 T. Buhl, C. Pedersen, O. Sigmund, Stiffness design of geometrically nonlinear structures using topology optimization, Struct. Multidiscip. Optim. 19, 93–104 (2000) [Google Scholar]
 H.C. Gea, J. Luo, Topology optimization of structures with geometrical nonlinearities, Comput. Struct. 79, 1977–1985 (2001) [Google Scholar]
 S.M. Han, S.I. Kim, Y.Y. Kim, Topology optimization of planar linkage mechanisms for path generation without prescribed timing, Struct. Multidiscipl. Optim. 56, 501–517 (2017) [Google Scholar]
 C.B.W. Pedersen, T. Buhl, O. Sigmund, Topology synthesis of largedisplacement compliant mechanisms, Int. J. Numer. Methods Eng. 50, 2683–2705 (2001) [Google Scholar]
 M.P. Bendsøe, J.M. Guedes, S. Plaxton, J.E. Taylor, Optimization of structure & material properties for solids composed of softening material, in IUTAM Symposium on Optimization of Mechanical Systems. Springer, Netherlands (1996), pp. 17–24 [Google Scholar]
 M. Bogomolny, O. Amir, Conceptual design of reinforced concrete structures using topology optimization with elastoplastic material modeling, Int. J. Numer. Methods En. 90, 1578–1597 (2012) [Google Scholar]
 K. Maute, S. Schwarz, E. Ramm, Adaptive topology optimization of elastoplastic structures, Struct. Optim. 15, 81–91 (1998) [Google Scholar]
 G.H. Yoon, Y.Y. Kim, Topology optimization of materialnonlinear continuum structures by the element connectivity parameterization, Int. J. Numer. Methods Eng. 69, 2196–2218 (2007) [Google Scholar]
 K. Yuge, N. Iwai, N. Kikuchi, Optimization of 2d structures subjected to nonlinear deformations using the homogenization method, Struct. Optim. 17, 286–299 (1999) [Google Scholar]
 K. Yuge, N. Kikuchi, Optimization of a frame structure subjected to a plastic deformation, Struct. Optim. 10, 197–208 (1995) [Google Scholar]
 X. Huang, Y. Xie, Topology optimization of nonlinear structures under displacement loading, Eng. Struct. 30, 2057–2068 (2008) [Google Scholar]
 X. Huang, Y.M. Xie, G. Lu, Topology optimization of energyabsorbing structures, Int. J. Crashworthiness 12, 663–675 (2007) [Google Scholar]
 D. Jung, H.C. Gea, Topology optimization of nonlinear structures, Finite Elem. Anal. Des. 40, 1417–1427 (2004) [Google Scholar]
 M.J. Werkheiser, J. Dunn, M.P. Snyder, J. Edmunson, K. Cooper, M.M. Johnston, 3d printing in zerog ISS technology demonstration, in AIAA SPACE 2014 Conference and Exposition. American Institute of Aeronautics and Astronautics (2014) [Google Scholar]
 C. ConlanSmith, A. Bhattacharyya, K.A. James, Optimal design of compliant mechanisms using functionally graded materials, Struct. Multidiscipl. Optim. 57, 197–212 (2017) [Google Scholar]
 R. Yang, C. Chen, Stressbased topology optimization, Struct. Optim. 12, 98–105 (1996) [Google Scholar]
 G. Da Silva, E. Cardoso, Stressbased topology optimization of continuum structures under uncertainties, Comp. Methods Appl. Mech. Eng. 313, 647–672 (2017) [Google Scholar]
 C. Kiyono, S. Vatanabe, E. Silva, J. Reddy, A new multipnorm formulation approach for stressbased topology optimization design. Comp. Struct. 156, 10–19 (2016) [Google Scholar]
 S.J. Moon, G.H. Yoon, A newly developed QPrelaxation method for element connectivity parameterization to achieve stressbased topology optimization for geometrically nonlinear structures, Comp. Methods Appl. Mech. Eng. 265, 226–241 (2013) [Google Scholar]
 A. Verbart, M. Langelaar, F. van Keulen, A unified aggregation and relaxation approach for stressconstrained topology optimization, Struct. Multidiscipl. Optim. 55, 663–679 (2016) [Google Scholar]
 C. Le, J. Norato, T. Bruns, C. Ha, D. Tortorelli, Stressbased topology optimization for continua, Struct. Multidiscipl. Optim. 41, 605–620 (2010) [Google Scholar]
 M. Bruggi, P. Duysinx, Topology optimization for minimum weight with compliance and stress constraints. Struct. Multidiscipl. Optim. 46, 369–384 (2012) [Google Scholar]
 N.H. Kim, Introduction to Nonlinear Finite Element Analysis. Springer US (2014) [Google Scholar]
 S. Coniglio, J. Morlier, C. Gogu, R. Amargier, Generalized geometry projection: a unified approach for geometric feature based topology optimization, Arch. Comput. Methods Eng. 1–38 (2019) [Google Scholar]
 W. Zhang, W. Yang, J. Zhou, D. Li, X. Guo, Structural topology optimization through explicit boundary evolution, J. Appl. Mech. 84, 011011 (2016) [Google Scholar]
 W. Zhang, J. Yuan, J. Zhang, X. Guo, A new topology optimization approach based on moving morphable components (MMC) and the ersatz material model, Struct. Multidiscipl. Optim. 53, 1243–1260 (2015) [Google Scholar]
 G.I. Rozvany, Aims, scope, methods, history and unified terminology of computeraided topology optimization in structural mechanics, Struct. Multidiscip. Optim. 21, 90–108 (2001) [Google Scholar]
 G.I. Rozvany, M. Zhou, T. Birker, Generalized shape optimization without homogenization, Struct. Optim. 4, 250–252 (1992) [Google Scholar]
 M. Zhou, G. Rozvany, The COC algorithm, Part II: Topological, geometrical and generalized shape optimization, Comp. Methods Appl. Mech. Eng. 89, 309–336 (1991) [Google Scholar]
 T. Bruns, A reevaluation of the simp method with filtering and an alternative formulation for solid–void topology optimization, Struct. Multidiscip. Optim. 30, 428–436 (2005) [Google Scholar]
 O. Sigmund, J. Petersson, Numerical instabilities in topology optimization: a survey on procedures dealing with checkerboards, meshdependencies and local minima, Struct. Optim. 16, 68–75 (1998) [Google Scholar]
 K. Svanberg, The method of moving asymptotes—a new method for structural optimization, Int. J. Numer. Methods Eng. 24, 359–373 (1987) [Google Scholar]
 R.W. Ogden, NonLinear Elastic Deformations. Courier Corporation (2013) [Google Scholar]
 P.Y. Papalambros, D.J. Wilde, Principles of Optimal Design. Cambridge University Press (2000) [Google Scholar]
 P. Duysinx, M.P. Bendsøe, Topology optimization of continuum structures with local stress constraints, Int. J. Numer.Methods Eng. 43, 1453–1478 (1998) [Google Scholar]
 T. Bruns, D. Tortorelli, Topology optimization of nonlinear elastic structures and compliant mechanisms. Comput. Methods Appl. Mech. Eng. 190, 3443–3459 (2001) [Google Scholar]
 G. Kreisselmeier, R. Steinhauser, Systematic control design by optimizing a vector performance index, in Computer aided design of control systems. Elsevier (1980), pp. 113–117 [Google Scholar]
Cite this article as: G. Capasso, J. Morlier, M. Charlotte, S. Coniglio, Stressbased topology optimization of compliant mechanisms using nonlinear mechanics, Mechanics & Industry 21, 304 (2020)
All Tables
Default parameters used in the Inverter optimization problem (symmetrical problem).
All Figures
Fig. 1 Initialdomain: airfoil optimization problem. 

In the text 
Fig. 2 Geometric inverter scheme. 

In the text 
Fig. 3 Results relative to the Airfoil problem (first formulation). 

In the text 
Fig. 4 Results relative to the Airfoil problem (second formulation). 

In the text 
Fig. 5 Results relative to the Inverter problem. 

In the text 
Fig. 6 Results relative to the inverter problem (σ = 0.2 MPa). 

In the text 
Fig. 7 Results relative to the inverter problem (σ = 0.4 MPa). 

In the text 
Fig. A.1 Results relative to the inverter problem: applied force F_{app} = 0.15 N/mm. 

In the text 
Fig. A.2 Results relative to the inverter problem: Variation of the applied force F_{app} =0.25 N/mm. 

In the text 
Fig. A.3 Convergence curves for problem (24): variations in σ_{lim}. 

In the text 
Fig. A.4 Evolution of optimization constraints (expressed in negative null form) for problem (24): variations in σ_{lim}. 

In the text 
Fig. A.5 Convergence curves for problem (24): variations in F_{app}. 

In the text 
Fig. A.6 Evolution of optimization constraints (expressed in negative null form) for problem (24): variations in F_{app}. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.