domain: spi chim phys

text:
Lander Vertonghen

François-Xavier Irisarri

Dimitri Bettebghor

Boris Desmorat

Gradient-based concurrent topology and anisotropy optimization for mechanical structures

Keywords: composite, laminate, density optimization, polar parameters

The present paper introduces a novel gradient-based optimization framework to obtain discrete topologies consisting of either none or complete presence of anisotropic material, modeled by means of the polar formalism. Current topology optimizations with the polar formalism are based on optimality criteria, and limited to performing compliance minimization for thermodynamically feasible materials. The proposed optimization approach uses a sequential approximations approach, based on the Methods of Moving Asymptotes. The material density, orientation and anisotropic modules are updated separately at each iteration, in parallel sub-problems constructed with different types of approximations and settings. The proposed optimization approach is successfully validated for compliance minimization against the Alternate Directions method for general orthotropic materials, defined by the thermodynamic bounds. The importance of the anisotropy initialization in the gradient-based approach is highlighted to obtain stiffer solutions. The gradient-based strategy is also extended to incorporate geometric bounds on the polar parameters, defining the domain of existence of composite laminates. Obtained results for compliance minimization with laminates are compared to published results using lamination parameters. Finally, new solutions are presented showing the improvement of the compliance with the expansion of the anisotropy design domains. This paper paves the way for strength-based topology optimization of parts made of orthotropic materials or composite laminates.

Introduction

In the field of transportation, the reduction of greenhouse gas emissions and costs implies to always reduce the amount of material used and thus the total mass of the vehicles while ensuring the mechanical resistance. During the design phases, structural mechanics engineers in design offices try to reach the best compromises through the use of structural optimization algorithms to find efficient shapes. Material mechanics engineers meanwhile seek the best mix of materials to take into account mechanical stresses and environment issues. These two design steps have traditionally been separated in design offices. This publication aims to consider both these steps simultaneously.

Early structural optimization methods involved size and shape optimization (Eschenauer and Olhoff [START_REF] Eschenauer | Topology optimization of continuum structures: A review*[END_REF], Rozvany [START_REF] Rozvany | A critical review of established methods of structural topology optimization[END_REF]), which is akin to a parametric study on predefined structure characteristics. Later, topology optimization was introduced: it has a more versatile application, by defining a structure without prior knowledge or assumption of underlying characteristics nor shape [START_REF] Bendsøe | Topology Optimization[END_REF]. First topology optimizations considered structures made of isotropic material [START_REF] Sigmund | On the Design of Compliant Mechanisms Using Topology Optimization[END_REF]. The topology is the design variable, and can be characterized with density [START_REF] Sigmund | A 99 line topology optimization code written in Matlab[END_REF] or boundary methods [START_REF] Dunning | Level-set topology optimization with many linear buckling constraints using an efficient and robust eigensolver: buckling constrained topology optimization[END_REF]. The former assigns a density representing the presence or absence of material, whereas the latter describes the boundary of the topology by means of mathematical functions. The anisotropy of the material elastic properties can be considered in such an optimization framework. Simple optimizations with a constant anisotropic stiffness over the domain are equivalent to isotropic ones. Nonetheless, more complex optimizations involving anisotropy have been studied.

Taking material anisotropy into account as an additional design variable allows to further increase the performances of the optimized parts. Anisotropy as used in topology optimization has been characterized in different ways. Some authors used microscale structures composed of void and isotropic material to obtain anisotropic homogenized elastic behaviors [START_REF] Wu | Topology optimization of multi-scale structures: a review[END_REF]. Homogenization techniques consists in defining the macroscopic elastic tensor of an underlying micro structure [START_REF] Allaire | Homogenization[END_REF]. Other authors chose to optimize the orientation of a transversely isotropic material, by rotating the stiffness tensor [START_REF] Schmidt | Structural topology optimization with smoothly varying fiber orientations[END_REF][START_REF] Nomura | General topology optimization method with continuous and discrete orientation design using isoparametric projection[END_REF][START_REF] Yan | Concurrent topology optimization of structures and orientation of anisotropic materials[END_REF][START_REF] Smith | Simultaneous material and topology optimization of composite laminates[END_REF][START_REF] Papapetrou | Stiffness-based optimization framework for the topology and fiber paths of continuous fiber composites[END_REF]. This is similar to considering anisotropy in 3D printing [START_REF] Zhou | Anisotropic Multicomponent Topology Optimization for Additive Manufacturing With Build Orientation Design and Stress-Constrained Interfaces[END_REF][START_REF] Dapogny | Shape and topology optimization considering anisotropic features induced by additive manufacturing processes[END_REF], all relating the anisotropy variables to the orientation and pattern of the physical fibers. These strategies can be classified as Continuous Fiber Angle Optimization (CFAO) [START_REF] Jiang | Continuous Fiber Angle Topology Optimization for Polymer Composite Deposition Additive Manufacturing Applications[END_REF]. Boddeti et al. [START_REF] Boddeti | Optimal design and manufacture of variable stiffness laminated continuous fiber reinforced composites[END_REF] used such CFAO framework to optimize the orientations in individual layers of a variable stiffness laminate, while using topology optimization to identify the presence or not of the laminate. Alternatively, in a layer by layer approach, the modules of anisotropic macroscopic materials can be optimized by considering a multi-material problem, where the stiffness properties are obtained as a combination of the stiffness of pre-selected materials such as in the Discrete Material Optimization (DMO) framework [START_REF] Stegmann | Discrete material optimization of general composite shell structures[END_REF]. DMO has later been extended to the Discrete Material and Thickness Optimization (DMTO) framework, to represent thickness variation by means of ply drops in laminate optimizations [START_REF] Sørensen | DMTO -a method for Discrete Material and Thickness Optimization of laminated composite structures[END_REF].

On the other hand, the components of the stiffness tensor in 2D or 3D can directly be the optimization variables as performed in the Free Material Optimization (FMO) framework [START_REF] Zowe | Free material optimization via mathematical programming[END_REF]. Nonetheless, FMO requires compatibility constraints on these components to ensure the stiffness tensor remains mechanically admissible. Thereafter, to identify microstructures obtained by FMO, inverse homogenization techniques are used [START_REF] Schury | Efficient Two-Scale Optimization of Manufacturable Graded Structures[END_REF]. In such a context, Hu et al. [START_REF] Hu | Cellular structure design based on free material optimization under connectivity control[END_REF] considered the connectivity of adjacent microstructures, whereas Tyburec et al. [START_REF] Tyburec | Modular-topology optimization of structures and mechanisms with free material design and clustering[END_REF] applied clustering to FMO results to seek modular designs for the microstructures. Instead, Weldeyesus and Stolpe [START_REF] Weldeyesus | Free material optimization for laminated plates and shells[END_REF] proposed a restriction of the FMO approach to represent laminates. FMO for membrane stiffness is considered for separate layers, and assembled taking into account shell kinematics.

Topology optimizations of shell structures have also been conducted by means of lamination parameters, which represents the homogenized stiffness of a complete laminate [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF][START_REF] Hu | A hybrid multilevel method for simultaneous optimization design of topology and discrete fiber orientation[END_REF][START_REF] Tong | Optimization of Combining Fiber Orientation and Topology for Constant-Stiffness Composite Laminated Plates[END_REF][START_REF] Bohrer | Concurrent topology and stacking sequence optimization of composite laminate plates using lamination parameters[END_REF]. In these works, the laminate macroscopic stiffness tensor is optimized with respect to a fixed local frame. Both types of variables, the material elastic modules and the local frame, can be conveniently handled using the polar formalism [START_REF] Vannucci | Plane Anisotropy by the Polar Method[END_REF], an invariant-based representation of the elastic tensors. The polar formalism can be applied not only to composite laminates, but to any 2D material. The domain of existence of the polar parameters for any material is delimited by the thermodynamic bounds. These thermodynamic bounds are a simple condition on the polar parameters, and are the equivalent to the constraints on the components of the stiffness tensor in FMO. Laminated structures constituted from a known base ply can also be defined conveniently as a subset of the thermodynamic domain on the polar parameters, the geometric domain [START_REF] Vannucci | A Note on the Elastic and Geometric Bounds for Composite Laminates[END_REF]. The corresponding geometric bounds have only been used thus far in thickness optimization problems [START_REF] Izzi | Strength and mass optimisation of variable-stiffness composites in the polar parameters space[END_REF][START_REF] Montemurro | A general multi-scale two-level optimisation strategy for designing composite stiffened panels[END_REF]. Nonetheless, the polar parameters have been used as design variables in a topology optimization considering the thermodynamic bounds in [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF].

Compliance minimization problems have been optimized with a variety of algorithm types. Sigmund [START_REF] Sigmund | A 99 line topology optimization code written in Matlab[END_REF] used a heuristic optimality criterion update scheme for isotropic materials, whereas the optimality criterion used in the Alternate Directions (AD) of Allaire et al. [START_REF] Allaire | The homogenization method for topology and shape optimization. Single and multiple loads case[END_REF] has a proof of convergence, for both isotropic and anisotropic material [START_REF] Ranaivomiarana | Simultaneous optimization of topology and material anisotropy for aeronautic structures[END_REF]. The AD algorithm is used for topology optimization of anisotropic materials in [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF] with the restriction of being only applicable to a compliance minimization problem within the thermodynamic bounds. Other algorithms to solving compliance minimization include gradient-based [START_REF] Sigmund | Topology optimization approaches[END_REF] and metaheuristic ones [START_REF] Wu | Topology optimization of structures using modified binary differential evolution[END_REF]. The latter two algorithm classes have the advantage of being able to consider more versatile mechanical responses as constraint or objective, which is important for actual engineering usage. But as was noted by Sigmund [START_REF] Sigmund | On the usefulness of non-gradient approaches in topology optimization[END_REF], non-gradient metaheuristic algorithms are computationally more expensive than gradient-based methods. Gradient-based topology optimization for isotropic material with different mechanical considerations have been carried out considering material stress [START_REF] Gao | Improving the overall performance of continuum structures: A topology optimization model considering stiffness, strength and stability[END_REF][START_REF] Holmberg | Stress constrained topology optimization[END_REF][START_REF] Le | Stress-based topology optimization for continua[END_REF][START_REF] Senhora | Topology optimization with local stress constraints: a stress aggregation-free approach[END_REF][START_REF] Capasso | Stress-based topology optimization of compliant mechanisms using nonlinear mechanics[END_REF], linear buckling [START_REF] Gao | Improving the overall performance of continuum structures: A topology optimization model considering stiffness, strength and stability[END_REF][START_REF] Ferrari | Revisiting topology optimization with buckling constraints[END_REF], thermo-elastic buckling [START_REF] Dunning | Level-set topology optimization with many linear buckling constraints using an efficient and robust eigensolver: buckling constrained topology optimization[END_REF], vibrations [START_REF] Zargham | Topology optimization: a review for structural designs under vibration problems[END_REF], compliant mechanisms [START_REF] Sigmund | On the Design of Compliant Mechanisms Using Topology Optimization[END_REF][START_REF] Capasso | Stress-based topology optimization of compliant mechanisms using nonlinear mechanics[END_REF] and manufacturing constraints [START_REF] Zhou | Gradient-based multi-component topology optimization for stamped sheet metal assemblies (MTO-S)[END_REF][START_REF] Zuo | Manufacturing-and machining-based topology optimization[END_REF]. Among the gradient methods, one can mention the Augmented Lagrangian method as used by [START_REF] Senhora | Topology optimization with local stress constraints: a stress aggregation-free approach[END_REF], and equally the broadly used Method of Moving Asymptote (MMA) [START_REF] Svanberg | MMA and GCMMA -two methods for nonlinear optimization[END_REF] for a range of mechanical responses [START_REF] Gao | Improving the overall performance of continuum structures: A topology optimization model considering stiffness, strength and stability[END_REF][START_REF] Holmberg | Stress constrained topology optimization[END_REF][START_REF] Le | Stress-based topology optimization for continua[END_REF][START_REF] Ferrari | Revisiting topology optimization with buckling constraints[END_REF].

Additional mechanical design criteria can also be considered in topology optimization problems with anisotropic materials using gradient-based approaches. Mirzendehdel et al. [START_REF] Mirzendehdel | Strength-based topology optimization for anisotropic parts[END_REF] considered anisotropic strength allowables in 3D printing, yet with isotropic material stiffness. Ma et al. [START_REF] Ma | Strength-constrained simultaneous optimization of topology and fiber orientation of fiberreinforced composite structures for additive manufacturing[END_REF] instead incorpo-rated a Tsai-Wu failure criterion in a CFAO framework with Uni-Directional fibers. In case of laminates, Lund used the DMO framework for buckling maximization [START_REF] Lund | Buckling topology optimization of laminated multi-material composite shell structures[END_REF] or considered failure criteria in DMTO [START_REF] Lund | Discrete Material and Thickness Optimization of laminated composite structures including failure criteria[END_REF]. Alternatively, the FMO framework for topology optimizations has been used with strength constraints in 2D [START_REF] Kočvara | Free material optimization for stress constraints[END_REF][START_REF] Kočvara | Free material optimization: recent progress[END_REF] and 3D [START_REF] Haslinger | Multidisciplinary Free Material Optimization[END_REF], or fundamental eigenfrequency [START_REF] Stingl | Free Material Optimization with Fundamental Eigenfrequency Constraints[END_REF]. Weldeyesus [START_REF] Weldeyesus | Free material optimization with local stress constraints for laminated structures[END_REF] also considered strength criteria, in the dedicated FMO framework for laminated plates and shells.

The present work aims at combining the advantages of the polar parametrization of the elasticity tensors with the versatility of gradient-based optimization. The use of a gradient-based framework enables to perform the 2D material anisotropy optimization either in the thermodynamic domain corresponding to the orthotropic materials or in the geometric domains of the laminated composite materials. It also opens the way to the incorporation of additional engineering constraints in topology optimization whilst considering material anisotropy. Thus, in the following, the algorithm to solving the optimizations is switched from the AD algorithm used in [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF][START_REF] Ranaivomiarana | Simultaneous optimization of topology and material anisotropy for aeronautic structures[END_REF] to MMA. Such a sequential approximation approach has been applied over large scale optimization problems to incorporate different objectives and constraints including composite structures with fixed topology, such as [START_REF] Bruyneel | Composite structures optimization using sequential convex programming[END_REF], as well as for the largest topology optimization problems, see [START_REF] Mukherjee | Accelerating Large-scale Topology Optimization: State-of-the-Art and Challenges[END_REF] for a review. A sequential approximation approach is therefore well-suited to incorporate the additional design variables used to describe the material anisotropy.

This articles therefore proposes a MMA based optimization strategy for compliance minimization with respect to both topology and anisotropy, represented by means of the polar parameters. The method is first validated against the AD algorithm as benchmark with thermodynamic bounds for orthotropic materials. Thereafter, the strategy is extended to considering composite laminates only, by means of the geometric bounds on the polar parameters [START_REF] Vannucci | A Note on the Elastic and Geometric Bounds for Composite Laminates[END_REF]. The remainder of this paper is setup as follows: the background of the topology and anisotropy modeling is first presented in Section 2. Section 3 details the optimization strategy, followed by the gradient calculations in Section 4. The results of the optimizations are presented and discussed in Section 5 and finally this paper is concluded in Section 6.

Methodology

This work deals with a two-dimensional compliance minimization problem. The corresponding mathematical formulation reads as follows: min

ρ,γ C(ρ, γ) s.t. V (ρ) ≤ V 0 ρ ∈ [ρ min , 1] g(γ) < 0 (1)
where C = U T KU is the compliance. ρ represents the topology distributed optimization variable, and the material anisotropy distributed variables are regrouped under the γ variables. The volume V (ρ) is bound by a prescribed volume V 0 constraint, while the anisotropy variables are subjected to physical existence constraints g(γ). The mechanical problem is solved by means of the Finite Element Method (FEM), where U is the solution to KU = F, K being the global stiffness matrix and F the force vector.

Topology parametrization

The topology is parametrized by means of a density approach, where each element e is assigned a scalar density value ρ e . The densities take value in [ρ min , 1], where ρ min is set to 10 -3 to avoid having a near singularity stiffness matrix and subsequent numerical problems with the FEM resolution. To obtain distinct black an white solutions, an exponent p is used to penalize intermediate densities and converge to the presence (ρ = 1) or absence (ρ = ρ min ) of material. This method is similar to the classical SIMP (Solid Isotropic Material with Penalization) approach, although it is applied here to the case of anisotropic materials. Thus, each element e is assigned a stiffness tensor Q e used for the mechanical analysis:

Q e = ρ p e Q (0) e , (2)
where Q (0) e is the pristine elastic tensor of the element. In order to have a length scale control and avoid numerical artifacts such as checkerboard, a density filter is used [START_REF] Sigmund | Numerical instabilities in topology optimization: A survey on procedures dealing with checkerboards, mesh-dependencies and local minima[END_REF]. The filtered densities ρ e are the physical variables used to setup the FEM analysis with Equation 2 and compute the volume constraint, while the densities ρ i are the optimization variables:

ρ e = i∈Ωe w ei ρ i A i i∈Ωe ρ i A i . ( 3
)
A i is the element area and Ω e the set of elements i whose centroid x i is within the filter radius R of the centroid x e of element e. The associated filter weights w ei are obtained as

w ei = R -∥x i -x e ∥ . (4)

Material parametrization

Material in-plane anisotropy is represented by means of the polar formalism [START_REF] Verchery | Les invariants des tenseurs d'ordre 4 du type de l'élasticité[END_REF]. The elastic behavior of the to be designed material is imposed to be orthotropic, in which case the polar formalism expresses its fourth order elasticity tensor by means of 5 invariants [START_REF] Vannucci | Plane Anisotropy by the Polar Method[END_REF]: the four moduli T 0 , T 1 , R 0 and R 1 and the orthotropy shape parameter K. The parameter K takes a discrete value of 0 or 1. Equation 5shows the relation between the Cartesian components and the polar ones for an orthotropic stiffness tensor, where the angle ϕ 1 represents the direction of orthotropy:

Q (0) 1111 = T 0 + 2T 1 + (-1) K R 0 cos 4ϕ 1 + 4R 1 cos 2ϕ 1 Q (0) 1122 = -T 0 + 2T 1 -(-1) K R 0 cos 4ϕ 1 Q (0) 1112 = (-1) K R 0 sin 4ϕ 1 + 2R 1 sin 2ϕ 1 Q (0) 2222 = T 0 + 2T 1 + (-1) K R 0 cos 4ϕ 1 -4R 1 cos 2ϕ 1 Q (0) 2212 = -(-1) K R 0 cos 4ϕ 1 + 4R 1 sin 2ϕ 1 Q (0) 1212 = T 0 -(-1) K R 0 cos 4ϕ 1 . (5)
T 0 and T 1 dictate the spherical behavior of the stiffness tensor and must both be strictly positive. They are kept constant during the optimization, as the trivial stiffest solution would be for T 0 and T 1 to have infinite values. Thus, the material optimization is performed with respect to the parameters that influence the anisotropic terms in Equation 5: the polar invariants R 0 , R 1 and K and the orientation ϕ 1 . The normalized anisotropic material modules η 0 and η 1 are introduced, defined as:

η 0 = (-1) K R 0 T 0 , η 1 = R 1 √ T 0 T 1 . (6)
The normalization not only allows to reduce the amount of variables, by regrouping R 0 and K into η 0 , but foremost includes the discrete variable K into a continuous one. The continuous variables facilitate the gradient-based optimization process, which does not perform well with discrete variables. As the isotropic moduli T 0 and T 1 remain constant during the optimization, they do not require normalization. The elasticity tensor of Equation 5 is now defined as follows:

Q (0) 1111 = T 0 + 2T 1 + η 0 T 0 cos 4ϕ 1 + 4η 1 √ T 0 T 1 cos 2ϕ 1 Q (0) 1122 = -T 0 + 2T 1 -η 0 T 0 cos 4ϕ 1 Q (0) 1112 = η 0 T 0 sin 4ϕ 1 + 2η 1 √ T 0 T 1 sin 2ϕ 1 Q (0) 2222 = T 0 + 2T 1 + η 0 T 0 cos 4ϕ 1 -4η 1 √ T 0 T 1 cos 2ϕ 1 Q (0) 2212 = -η 0 T 0 cos 4ϕ 1 + 4η 1 √ T 0 T 1 sin 2ϕ 1 Q (0) 1212 = T 0 -η 0 T 0 cos 4ϕ 1 . (7)

Thermodynamic bounds and geometric bounds

In the most general case, for 2D orthotropic materials, η 0 and η 1 have to take value within the thermodynamic bounds that ensure the positive definiteness of the elasticity tensor Q (0) [START_REF] Vannucci | Plane Anisotropy by the Polar Method[END_REF]:

2η 2 1 -1 < η 0 < 1 . (8)
The domain corresponding to all possible laminates constructed from a single base layer, the so-called geometric domain [START_REF] Vannucci | A Note on the Elastic and Geometric Bounds for Composite Laminates[END_REF], forms a subset of the thermodynamically admissible domain. With η L 0 and η L 1 the normalized anisotropic material properties of the base orthotropic layer, the geometric bounds for an orthotropic laminate reads as:

2 η 1 η L 1 2 -1 ≤ η 0 η L 0 ≤ 1 . (9)
Interestingly, Equation 8and Equation 9share the same form. Thus, by taking a base ply whose normalized anisotropic properties tends towards 1, i.e. η L 0 = η L 1 = 1 -ε, with ε an infinitely small positive number, the geometric domain tends to overlap with the thermodynamic domain and both domains can be described using Equation 9.

Remapping

The feasible domain for the normalized anisotropic variables η 0 and η 1 is remapped using a change of variables. A similar approach is proposed in [START_REF] Izzi | Strength and mass optimisation of variable-stiffness composites in the polar parameters space[END_REF], with a different variable change. The relation between the new design variables, α and β, and the anisotropic modules, η 0 and η 1 , is the following:

η 0 η L 0 = 1 -8α(1 -α)β η 1 η L 1 = 2α -1 , (10)
From Equation 9, the inequality 2

η1 η L 1 2 -1 ≤ 1 implies that -1 ≤ η1 η L 1 ≤ 1, therefore Equation 10 implies 0 ≤ α ≤ 1. When considering 0 < α < 1, the inequality 1 -η0 η L 0 ≥ 0 reads 8α(1 -α)β so that β ≥ 0. Also, the inequality η0 η L 0 -(2 η1 η L 1 2 -1) ≥ 0 reads 8α(1 -α)(1 -β) ≥ 0 so that β ≤ 1.
When considering α = 0 or 1, meaning η0 η L 0 = 1, β is undetermined and can be chosen in [0, 1]. By using the proposed remapping, the non-linear feasible domain boundary given in Equation 9can thus be replaced by α, β ∈ [0, 1] × [0, 1], and the Problem 1 can be rewritten as:

min ρ,ϕ1,α,β C(ρ, ϕ 1 , α, β) s.t. V ≤ V 0 ρ ∈ [ρ min , 1] ϕ 1 ∈ [-π/2, π/2] α, β ∈ [0, 1] × [0, 1] . ( 11
)
The change of optimization constraint contributes to significantly reduce the computational cost, but even more importantly, it ensures that the optimized material always satisfies the thermodynamic or geometric bounds, thus avoiding failure of the optimization process due to material singularities in the FEM computation.

Optimization strategy

The MMA optimization algorithm [START_REF] Svanberg | MMA and GCMMA -two methods for nonlinear optimization[END_REF] is used to solve Problem 11. The MMA algorithm is based on an iterative optimization process that consists in solving a succession of approximated convex sub-problems. At each iteration, a sub-problem is constructed based on the values and gradients at the current iteration point of both the objective and constraint functions with respect to the optimization variables.

Several approximation types can be used to construct the sub-problem, and have been found to influence the convergence and the progress of the optimization. Since the regularity of the design responses taken into account in the problem differs with respect to the nature of the variable, a strategy called SplitMMA is suggested to solve Problem 11. This strategy uses approximations which are already better tailored to a variable's influence, to facilitate convergence as was hinted by Bruyneel et al. [START_REF] Bruyneel | A family of MMA approximations for structural optimization[END_REF].

The general functioning of MMA is retained in the suggested strategy. But instead of constructing a unique sub-problem taking into account all the variables at each iteration, three separate sub-problems are devised depending on the type of variable (density, orientations and anisotropic modules) and their derivatives. These sub-problems are then solved in parallel at each iteration, meaning once per FEM analysis. All gradients can then be computed from this common mechanical analysis. The old variables and asymptote update history is however kept independently per sub-problem throughout the optimization. This process parallelization leads to neglecting the cross-influences of the sub-sets of variables, which has little influence on the resolution of Problem 11. Indeed, the volume and domain of existence constraints are separate and do not depend on the same variable. Optimizing with respect to each type of variable separately with SplitMMA allows to construct more simple sub-problems to be solved, with independent settings for the solver. Those settings can be individually tailored to the behavior of the optimization in junction with an adequate selection of the approximation type to improve convergence. The density variables are optimized with the standard monotonous MMA approximation, shown to behave well in literature [START_REF] Bendsøe | Topology Optimization[END_REF][START_REF] Bruyneel | A family of MMA approximations for structural optimization[END_REF], and its own set of algorithm settings dictating the asymptotes update. Equally, to update the anisotropic modules η 0 and η 1 , the optimization with respect to α and β is performed with a separate call of MMA and its monotonous approximation, giving good convergence. Having this separate MMA call allows to use a different set of settings to guide the optimization of the anisotropic components. Finally, the orientations are optimized with the convex but non-monotonous Global Convergent MMA (GCMMA) without inner iterations. This helps to mitigate the highly non linear behavior and response surface of the compliance with respect to the orientation. This GCMMA call has again individual settings. Using GCMMA to optimize with respect to orientations is similar to the works performed by Kiyono et al. [START_REF] Kiyono | A novel fiber optimization method based on normal distribution function with continuously varying fiber path[END_REF]. An overview of the optimization strategy is given in Figure 1.

Sensitivity analysis

The gradient of the volume with respect to the density variables can be obtained straightforwardly as follows, where V i is the volume of each element:

∂V ∂ ρi = V i . (12)
In the case of the compliance sensitivity, the gradient formulation is self adjoint [START_REF] Bendsøe | Topology Optimization[END_REF]. Therefore, the gradient with respect to density ρi of element i is given by:

∂C ∂ ρi = -d T i ∂k i ∂ ρi d i = -d T i pρ p-1 i k i d i , (13)
where d i is the displacement vector of the element's nodes at the current iteration, whilst k i is the element's local stiffness matrix.

The sensitivity with respect to the direction of anisotropy ϕ 1 can be obtained analogously to the sensitivities with respect to the density [START_REF] Yan | Concurrent topology optimization of structures and orientation of anisotropic materials[END_REF], as:

∂C ∂ϕ 1i = -d T i ∂k i ∂ϕ 1i d i . (14)
where the local stiffness matrix is given by:

k i = Si B T ρp i Q (0) BdS i , (15)
with B the strain-displacement matrix of the element, and S is the surface over which the element is integrated.

As the stiffness tensor does not depend on the coordinates of the element, the partial derivative of Equation 14 reads:

∂k i ∂ϕ 1i = Si B T ρp i ∂Q (0) i ∂ϕ 1i BdS i . (16)
The gradient of the stiffness tensor Q (0) i can be obtained straightforwardly from Equation 7, and the numerical integration is performed in the same way as in the original FEM analysis. The sensitivities with respect to η 0 and η 1 are obtained similarly.

To obtain the gradient with respect to the design variables α and β, the following chain rule must be applied:

∂C ∂α = ∂C ∂η 0 ∂η 0 ∂α + ∂C ∂η 1 ∂η 1 ∂α ∂C ∂β = ∂C ∂η 0 ∂η 0 ∂β + ∂C ∂η 1 ∂η 1 ∂β . ( 17
)
Finally, the influence of the density filter on the sensitivities is taken into account using a another chain rule, as follows:

∂f ∂ρ i = e∈Ωi ∂f ∂ ρ e ∂ ρ e ∂ρ i , (18)
with:

∂ ρ e ∂ρ i = w ei j∈Ωe w ej , (19)
where f can be either the objective or a constraint, such as the compliance or volume [START_REF] Gao | Improving the overall performance of continuum structures: A topology optimization model considering stiffness, strength and stability[END_REF].

Results and discussions

The proposed strategy is applied to the three test cases shown in Figure 2. The first considered test case (#1) is the cantilever beam problem of aspect ratio of 2:1. The prescribed volume fraction v f is 50%, with a mesh size of 0.5 × 0.5 mm 2 and a 1.1 mm filter radius R. A 200 N downwards load is distributed over the 5 middle nodes on the right hand size, while the left hand is clamped. The second test case (#2) is the suspended bridge problem with a prescribed volume fraction v f of 20%, a mesh size of 1.875 × 1.85 mm 2 and a 5.75 mm filter radius R. A 900 N downwards load is distributed over the deck, while some elements at the bottom are clamped. The deck is 2 elements thick. The clamping areas are symmetric, as defined in Figure 2. The black elements only have their anisotropy optimized, their densities remaining 1. Both test cases (#1) and (#2) are used in [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF]. The third test case (#3) is the cantilever beam problem of aspect ratio 1:1 used in [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF]. The prescribed volume fraction v f is 60%, with a mesh size of 5 × 5 mm 2 . A 9.5 mm filter radius R is assumed here, while in the original work, an implicit filter is used. A downwards distributed load totaling 10 N is applied across the 2 right most bottom nodes, while the left hand is clamped. A regular mesh is used in all cases, as is the most convenient for topology optimization. This however does not prohibit the proposed optimization strategy from being applied on a free mesh. For test case (#1) and (#2), the filter radius is taken according to the common practise as to have two to three elements within the radius. On the other hand, the radius in test case (#3) is taken as to match the topology of Peeters et al. [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF], since their filter is implicitly implemented in their approach.

The material properties used for test cases (#1) and (#2) are given in Table 1, while the material properties used for test case (#3) are presented in Table 2. The results presented in Section 5.1 are obtained using the thermodynamic bounds of Equation 8. These strict inequalities are enforced by using Equation 9with η L 0 = η L 1 = 1-ε, where a ε = 0.05 offset is imposed. The isotropic modules T 0 and T 1 as set equal to T L 0 and T L 1 respectively. In Section 5.2, results are obtained using the geometric bounds defined by Equation 9, in which case the materials given in Table 1 and Table 2 correspond to the base plies of composite laminates. The isotropic modules T 0 and T 1 are still equal to T L 0 and T L 1 . The continuation strategy for the optimizations is given in Table 3. All optimizations are carried out for a fixed number of iterations. The filter reduction phase consists in gradually decreasing the filter radius every 15 iterations until its value is smaller than the smallest element dimension. In this step, the SIMP exponent is still equal to 5. The filter reduction is used to minimize the amount of intermediate density values, measured by the measure of non-discreteness (M nd ) [START_REF] Sigmund | Morphology-based black and white filters for topology optimization[END_REF]: The algorithm settings for the SplitMMA strategy are given in Table 4. A Python implementation for MMA and GCMMA1 was used. The initial density and anisotropy distributions are uniform. The density is set to the imposed volume fraction of the optimization, whereas the anisotropy is initialized with ϕ 1 = 0 • , η 1 = 0.001.η 0 is set to its maximum possible value. For the thermodynamic bounds, with a ε = 0.05 offset to satisfy the strict inequality, η 0 is thus initialized as 0.95. When the geometric bounds are used, the anisotropy is initialized with ϕ 1 = 0 • , η 1 = 0.001 and η 0 = η L 0 . The AD optimizations that are used as reference in Section 5.1, are initialized with all its anisotropic components equal to 0.

M nd = n i=1 4 ρ i (1 -ρ i ) n × 100% . (20)
In the following, for the sake of clarity, the absolute value of the variable η 1 is shown, since a change of sign of η 1 corresponds to a 90 • rotation of the material orthotropy axes. The orientations ϕ 1 are also corrected accordingly, for easier visualization of the stiffest direction. The sequential approach is first presented to investigate the influence of the topology and anisotropy separately. It is applied to both test case (#1) and (#2). The first step consists in obtaining a topology with an isotropic material. This step, done with only the density variables in the SplitMMA strategy and AD algorithm (η 0 and η 1 are set to 0), results in the topology solutions of the sequential optimizations in Figure 3 and Figure 4. The anisotropic components, optimized in a separate stage, of both these sequential optimizations are equally shown in Figure 3 and Figure 4. Lastly, the convergence graph for all the iterations together is also shown. In case of test case (#1), the topology and compliance are similar between both AD and SplitMMA, be it with a higher M nd % for the AD algorithm. For test case (#2), the topology remains similar, be it with a higher compliance for SplitMMA. Furthermore, the local topology at the junction is different between SplitMMA and AD, with small holes emerging. Nonetheless, the disposition of the bars and junctions are coherent between both strategies.

Simultaneous optimization

This part presents the simultaneous optimization with respect to the topology and anisotropy for compliance minimization. The results of these concurrent optimizations with the AD and SplitMMA strategy for both test cases are given in Figure 3 and Figure 5 respectively. In case of the second SplitMMA optimization (anisotropy initialization) of test case (#2), the 35 initial iterations are performed only with respect to ϕ 1 , η 0 and η 1 , with all element densities kept at ρ = v f and SIMP p = 1. However, the total volume is greater than v f in this step, due to the presence of the imposed full elements. When the SIMP exponent p is increased after iteration 35, the compliance is minimized with respect to all variables as in the other concurrent optimizations. The volume constraint is then satisfied.

The results for test case (#1) have similar density and anisotropy distributions, for the same compliance yet small differences in M nd value. On the other hand, the SplitMMA solutions are more compliant than the AD solution for test case (#2). Nonetheless, the stiffest of the SplitMMA solutions, the one with the anisotropy initialization, resembles closely the AD solution. Furthermore, it is stiffer than any of the sequential optimization, with either the AD algorithm or SplitMMA strategy. Finally, Figure 8 shows the distribution of the anisotropic modules for the SplitMMA optimization within the thermodynamic domain for test case (#1).

Discussion

The non-discreteness measure (M nd ) is important in the comparison of the solutions with AD and SplitMMA. Indeed, despite the same filter radius being used in both the SplitMMA and AD optimizations, they are not the same filter. The AD algorithm in [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF] is programmed with an energy filter [START_REF] Desmorat | Structural rigidity optimization with frictionless unilateral contact[END_REF], which filters the deformation energy instead of only the densities. The most notable difference from this comes as the energy filter with a given filter radius will result in a topology with little intermediate densities. On the contrary, the density filter as used with the SplitMMA strategy with the same active filter radius will have intermediate densities due to the averaging and blurring effect on the boundary of the topology. Therefore the filter removal in SplitMMA lessens the blurring, and allows the topology to converge to more distinct density values, with lower M nd %. This has a twofold advantage. First, the distinct topology makes it easier to define the boundary of the structure for later post-processing steps in a design loop. Secondly, the intermediate densities are penalized the most by the SIMP approach. This therefore deteriorates the compliance for a similar topology with blurred boundaries. Having the intermediate densities removed and similar levels of M nd facilitates the comparison of the compliance value between solutions.

Thereafter, there is an equivalence between the isotropic topology with the AD and MMA algorithm for both test cases (the first step of the sequential optimizations). This was also reported by Fanni et al. [START_REF] Fanni | A Comparison between Different Topology Optimization Methods[END_REF] for optimization with isotropic material with optimality criteria, the basis of the AD algorithm, and MMA. On the other hand, both test case (#1) and (#2) show that the anisotropy changes the topology.

The comparison of the overall anisotropic distribution is in good agreement between AD and SplitMMA for all cases. Looking more in detail at η 0 of the SplitMMA solutions, its distribution is locally different along the jagged edge of the oblique bars between algorithms for the sequential approach of test case (#1) in Figure 3. Figure 6 shows the difference in optimized values of the polar parameters after having applied the optimality criteria on the final SplitMMA solution. η 0 values on the jagged edge differ significantly from the optimal values, also having a different sign. This different sign means the shape of orthotropy is different, and the overall elastic behavior changes. This leads to equally suboptimal values for ϕ 1 and η 1 .

Local negative η 0 values are also present in test case (#2), but in neither test case (#1) or (#2) does it have a major impact on the global distribution and compliance value.

The AD solutions furthermore show that η 0 is constant at its upper bound. This is true regardless of the stress state [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF], and the reason for the initialization choice of η 0 as the upper limit. On the other hand, the η 1 distribution is in good agreement, also when the optimality criteria are applied to the SplitMMA solution as given in Figure 6 and Figure 7. This clearly highlights that η 1 is well optimized with respect to the local stress state, be it with some noise present on the distributions. As a possible solution to the noise for the η 0 and η 1 distribution and negative η 0 locations, further improvements could look at including a filter on these anisotropic modules to smooth the variation.

Afterwards, the ϕ 1 distribution also corresponds well between algorithms. The final orientations in SplitMMA greatly coincide with the theoretical optimality conditions used in the AD algorithm, being aligned with the maximal of the absolute value of the principal stresses. This is visible in Figure 6 and7. A note should be made about the periodicity of the orientations: despite neighbouring elements having 90 • and -90 • orientations, the corresponding mechanical properties are identical. The discontinuity in the orientation field is not a source of error, but merely a periodicity effect. Therefore, to avoid 180 • orientation differences being averaged out, a prospective filter should address this issue, such as the strategy proposed by Schmidt et al. [START_REF] Schmidt | Structural topology optimization with smoothly varying fiber orientations[END_REF]. Nonetheless, all these elements of analysis show that the anisotropy is well optimized and conforming to the topology, and hence well taken into account in the optimization. Both test cases also highlight the importance of the initialization for the concurrent optimization with SplitMMA compared to the AD algorithm. This is because the optimality criteria in the AD can change the orientations and modules abruptly, as it only depends on the stress state of the current iteration. On the other hand, the variable values in the gradient approach can only evolve from the variable values of the previous iteration by the maximum stepsize allowed in MMA. The gradient-based approach seeks the best stepsize of the variables to improve the objective. However, the variable update is dictated by the validity of the gradient and its approximation, hence the variable change between iterations is usually less than optimality criteria changes. Figure 8: Anisotropic modules distribution within the thermodynamic bounds for the concurrent topology and anisotropy optimization with SplitMMA for test case (#1). Modules shown for elements with ρ ≥ 0.9.

For the uniform initialization, the SIMP p = 1 optimization step with respect to all variables is required. It allows for the interaction of the anisotropy and topology to converge to an intermediate configuration, even with these smaller steps. The intermediate convergence is due to the lack of penalization hindering any excess favouritism of any variable on the compliance. If the SIMP exponent p were to be raised too early, or the optimization started with SIMP p > 1, the synergy between anisotropy and topology would not have time to take place, as the higher density penalization would lead the optimization. In the case of SIMP p > 1, the smaller anisotropy variables steps and slower changes would be conforming to the already defined and predominant topology variables.

On the other hand, test case (#2) is started with an anisotropy only initialization. The optimization is carried out only with respect to the anisotropy in the first 35 iterations. This allows to have an anisotropy distribution already conforming to the load path over the design domain at the start of the topology optimization. Then the optimization with respect to all variables is begun with SIMP p = 3, meaning the influence of the anisotropy is directly taken into account by the topology variables. This is actually equivalent to the AD algorithm workings. Although the AD algorithm starts directly with SIMP p = 3, its local optimization starts by minimizing the compliance with respect to the anisotropy. Thereafter, the densities are updated with the optimal anisotropy but by keeping the initial stress distribution of the iteration, as one global iteration corresponds to one FEM analysis. But the anisotropy update is a discrete change, based only on the stress field. This is the same effect occurring in the first 35 iterations in SplitMMA, where only the anisotropy is optimized. Only, as stated earlier, the changes can not be abrupt with the gradient method, but rather with small steps from the previous iteration variable values to obtain the anisotropic field best align with the stress field. Adding to the small possible variable variation between iterations with SplitMMA, allowing both positive and negative η 1 in the optimization permits to quickly switch the direction of highest stiffness between two iterations, even though the orientation ϕ 1 can only vary by a few degrees. This property is useful when the direction of principal stress flips due to bar-like features being created during the optimization. In case of test case (#1), the anisotropic initialization is not required, only the SIMP p = 1 optimization step with all variables. This is allocate to the nature of the problem with many similar local minima, as highlighted by the AD optimizations showing both the sequential and concurrent compliance as equivalent.

This approach of having first an anisotropy initialization leads thus to a different topology, which for the bridge test case is less bulky, but foremost also a stiffer solution. From a physical perspective, the arch is extremely stiff in bending in case of the right solution in Figure 5 compared to middle solution in Figure 5, Furthermore, with a lower height, more material can be used to obtain an additional stiffening effect by increasing the inertia of the arch and having an added support.

Finally, the difference between the best AD and SplitMMA solution is because the gradient optimization is prone to getting stuck in local minima. Nonetheless, the difference is admissible, as the AD is specific to compliance minimization, but a wider range of optimization problems can be considered with the SplitMMA framework. Furthermore, SplitMMA upholds the concurrent approach yielding significant compliance gains over the sequential approach, by about 8% for the bridge case with the anisotropic initialization. This corresponds to the conclusion of Ranaivomiarana et al. [START_REF] Ranaivomiarana | Concurrent optimization of material spatial distribution and material anisotropy repartition for two-dimensional structures[END_REF]. All these observations allow to validate the suggested SplitMMA strategy, which is therefore well suited to incorporating anisotropy in topology optimization with thermodynamic materials. The solutions with SplitMMA converge to a similar solution to that of the benchmark AD algorithm, which is a powerful means for validation. These numerical experiments show that neglecting the cross-influences in SplitMMA does not hinder the strategy for compliance minimization. Indeed, the cross-influences are considered in the closed form solutions for the AD algorithm.

Optimization with laminates considered by means of the geometric bounds

Comparison with optimizations using lamination parameters

The topology optimization with geometric bounds on the polar parameters is validated by comparing the results obtained to similar optimization results with lamination parameters on test case (#3), available in [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF]. In the case of an orthotropic stiffness tensor, the relation between the lamination parameters and the polar parameters given by Panettieri et al. [START_REF] Panettieri | Blending constraints for composite laminates in polar parameters space[END_REF] can be expressed with normalized modules as follows:

V 1 = η 1 η L 1 cos 2ϕ 1 , V 2 = η 1 η L 1 sin 2ϕ 1 , V 3 = η 0 η L 0 cos 4ϕ 1 , V 4 = η 0 η L 0 sin 4ϕ 1 . (21)
To reproduce the same orthotropy conditions of a balanced laminate as Peeters et al. [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF], ϕ 1 is fixed to 0 • and not optimized. Only η 0 and η 1 are optimized concurrent to the topology, with a uniform initialization as used in the initial work. The optimization is performed in two steps with p = 1 and p = 3. The results of this balanced laminate optimization with polar parameters and geometric bounds is shown in Figure 9, where the polar parameters have been converted to lamination parameters through Equation 21.

The results of the optimization with geometric bounds in Figure 9, where the polar parameters have been converted to lamination parameters, are in good agreement for both topology and anisotropy distribution as compared to Peeters et al. [START_REF] Peeters | Combining topology and lamination parameter optimisation[END_REF]. This shows that the SplitMMA strategy is equally well suited to considering geometric bounds in a simultaneous topology and anisotropy optimization.

Topology optimizations with different types of anisotropy

Lastly, SplitMMA can be used to simulate different degrees of freedom of the anisotropy, as defined by the constraints placed on the anisotropy design variables, and therefore the type of anisotropy that is obtained. Figure 10 compares solutions obtained for different material for the same volume and mass with test case (#1). The isotropic case is optimized only with respect to the density, where T 0 = T L 0 , T 1 = T L 1 , η 0 = 0 and η 1 = 0. This corresponds to the first step of the sequential approach. The steered fiber case is optimized with respect to the density and orientation ϕ 1 only, with

T 0 = T L 0 , T 1 = T L 1 , η 0 = η L 0 and η 1 = η L 1 .
The geometric case optimizes with respect to all variables, limited to the geometric bounds based on η L 0 and η L 1 . The optimization solution with thermodynamic bounds is the one from Section 5.1.2. The visual representation of the geometric feasible domain within the thermodynamic domain is shown in Figure 11, for the properties as listed in Table 1. The optimizations with steered fiber and geometric bounds use the anisotropy initialization strategy.

Analyzing the solutions of the topology optimizations with different anisotropic degrees of freedom shows two important facts. First of all, all the topologies are different, demonstrating once again the importance of considering the anisotropy concurrently to the topology. Secondly, the more anisotropic freedom is given, the less compliant the solution becomes. The laminate anisotropy by means of the geometric bounds is stiffer than just steered fibers on a structure, although less than the thermodynamic bounds. Figure 11 shows the distribution of the anisotropic modules within the geometric domain. The reduced domain for the anisotropic modules with geometric bounds compared to the thermodynamic ones is the reason for the more compliant result: the restriction means less than ideal anisotropy can be used from a theoretical point of view. Even so, the optimized orthotropy direction ϕ 1 with geometric bounds still coincides with the direction of the maximum of the absolute value of the principal stresses, as shown in Figure 12.

The solution with the geometric bounds holds nonetheless an advantage over its counterpart with thermodynamic bounds. Anisotropy defined by the thermodynamic bounds is valid mathematically, but part of the domain has no known corresponding material. This is for example visible where η 0 tends to 1, meaning no shear stiffness. This is where the geometric bounds are superior, as they represent a feasible stacking sequence and therefore a material which can be manufactured. However, this stacking sequence retrieval, performed in subsequent steps after the optimization, is not part of the current topic. The geometric bounds had already been used for laminates with thickness optimization before [START_REF] Izzi | Strength and mass optimisation of variable-stiffness composites in the polar parameters space[END_REF][START_REF] Montemurro | A general multi-scale two-level optimisation strategy for designing composite stiffened panels[END_REF], but not yet incorporated in a topology optimization routine. The current framework offers the prospect of topology and anisotropy optimization with additional optimization constraints such as strength or buckling, for future research.

Conclusion

This paper presents an optimization strategy for the simultaneous optimization of topology and material anisotropy combining the Method of Moving Asymptotes (MMA) algorithms and the polar method, an invariant-based representation of the elasticity tensors. The density-based topology formulation allows to converge towards distinct designs (i.e., black and white solutions made of void (ρ = ρ min ) or anisotropic material (ρ = 1). The anisotropic modules are normalized, resulting at the same time in a set of continuous design variables, for easier use in the optimization. The considered materials are either general orthotropic materials or composite laminates. A change of variables is used to implicitly satisfy either of the corresponding existence constraints during the optimization. A strategy called SplitMMA is suggested to solve the compliance minimization problem, based on the MMA algorithm class. The density, orientation and anisotropic modules are updated separately at each iteration, in parallel sub-problems. Each sub-problem is constructed with a different type of approximation and settings, selected to tailor best the regularity of the problem with respect to the different types of variables. The proposed SplitMMA method is first compared with respect to the Alternate Direction (AD) algorithm in the case of general orthotropy in 2D. For both sequential and simultaneous topology and anisotropy optimizations, the SplitMMA solutions are a close match to their AD counterparts, both in objective value and variable distributions. In particular, the SplitMMA strategy preserves the improvement in compliance due the simultaneous consideration of topology and anisotropy. Moreover, the importance of the anisotropy initialization in SplitMMA is highlighted to achieve the improvements. The approach is further validated in the case of composite laminates by comparing the results obtained with published results using density and lamination parameters as design variables. Finally, the influence of different types of material anisotropy in topology optimization are presented with novel results for composite laminates. With the same prescribed mass and the same isotropic part of the material, the compliance of the solutions decreases with increasing anisotropy design domains, while the topologies vary. Further work will deal with the incorporation of stress criteria in the proposed SplitMMA framework, for strength-based topology optimization of 2D parts made of anisotropic materials.

Figure 1 :

 1 Figure 1: Overview of the optimization strategy SplitMMA.

Figure 2 :

 2 Figure 2: Test case representation.

Figure 3 :

 3 Figure 3: Results of the cantilever beam test case (#1) with thermodynamic bounds, for the sequential and simultaneous optimizations of the topology and anisotropy with the AD algorithm and SplitMMA strategy. Anisotropy variables are shown where ρ ≥ 0.9.

Figure 4 :

 4 Figure 4: Results of the bridge test case (#2) with thermodynamic bounds, for the sequential optimizations of the topology and anisotropy with the AD algorithm and SplitMMA strategy. Anisotropy variables are shown where ρ ≥ 0.9.

Figure 5 :

 5 Figure 5: Results of the bridge test case (#2) with thermodynamic bounds, for the simultaneous optimizations of the topology and anisotropy with the AD algorithm and SplitMMA strategy. Anisotropy variables are shown where ρ ≥ 0.9.

Figure 6 :

 6 Figure 6: Difference between the anisotropic values obtained with SplitMMA and the optimal values obtained by applying the optimality criteria to SplitMMA solution for the sequential approach of test case (#1).

Figure 7 :

 7 Figure 7: Difference between the anisotropic values obtained with SplitMMA and the optimal values obtained by applying the optimality criteria to SplitMMA solution for the concurrent approach of test case (#1).

Figure 9 :

 9 Figure 9: Optimization results of test case (#3) with geometric bounds: (a) the topology, (b) V 1 and (c) V 3 .

Figure 10 :

 10 Figure 10: Topology and anisotropy distributions for simultaneous optimizations with different degrees of anisotropy: isotropic material design, steered-fiber design, laminated solution (geometric bounds) and general orthotropic materials (thermodynamic bounds).

Figure 11 :

 11 Figure11: Anisotropic modules distribution within the geometric bounds for the concurrent topology and anisotropy optimization with SplitMMA for test case (#1). Modules included where ρ ≥ 0.9.

Figure 12 :

 12 Figure 12: Optimized anisotropy orientation and ϕ 1 for the geometric bounds optimization of Figure 10. Modules shown for elements with ρ ≥ 0.9.

Table 1 :

 1 Optimization and material properties for the optimizations of test case (#1) and (#2).

	Module	Value	Polar components Value
	Longitudinal Young Modulus [MPa] 181.0e3 Transverse Young Modulus [MPa] 10.3e3 Shear Modulus [MPa] 7.2e3 ν LT [-] 0.28	T L 0 [MPa] T L 1 [MPa] R L 0 [MPa] R L 1 [MPa] K L [-]	26.9e3 24.7e3 19.7e3 21.4e3 0

Table 2 :

 2 Optimization and material properties for the optimization of test case (#3).

	Module	Value	Polar components Value
	Longitudinal Young Modulus [MPa] 177.0e3 Transverse Young Modulus [MPa] 10.8e3 Shear Modulus [MPa] 7.6e3 ν LT [-] 0.27	T L 0 [MPa] T L 1 [MPa] R L 0 [MPa] R L 1 [MPa] K L [-]	26.6e3 24.3e3 19.0e3 20.9e3 0
	8

Table 3 :

 3 Continuation strategy of the optimization parameters in number of iterations for the different strategies.

		Sequential	Sequential	Concurrent	Concurrent
		SplitMMA	AD	SplitMMA	AD
	SIMP p = 1	-	-	35	35
	SIMP p = 3	90	90	90	90
	SIMP p = 5	35	35	35	35
	Filter reduction	90, where R new = 0.8 × R old every 15 iterations
	Anisotropy optimization	35	35	-	-
	Total number of iterations	250	250	250	250

Table 4 :

 4 MMA and GCMMA algorithm settings for the SplitMMA.

		MMA ρ GCMMA ϕ 1 MMA α & β
	epsimin	10 -10	10 -10	10 -10
	raa0	10 -5	0.01	10 -4
	raa0eps	-	10 -6	-
	move	0.5	-	0.2
	albefa	0.9	0.985	0.965
	asyinit	0.8	0.7	0.7
	asyincr	1.2	1.2	1.2
	asydecr	0.8	0.6	0.6
	5.1. Optimization with thermodynamic bounds
	5.1.1. Sequential optimization

https://github.com/arjendeetman/GCMMA-MMA-Python, accessed on the 17/11/2020

Acknowledgments

The authors are thankful to Airbus Atlantic for the financial support for this research, and specifically to Attilio Chiappini for the exchanges and discussions on the topic.

Declaration of competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Replication of results

We state that the results presented in this paper can be reproduced by the implementation details provided herein.