T Hirschler
email: thibaut.hirschler@insa-lyon.fr

R Bouclier

• A Duval

• T Elguedj

• J Morlier

Isogeometric sizing and shape optimization of thin structures with a solid-shell approach

Keywords: Shape Optimization, Sizing Optimization, Isogeometric Analysis, Stiffened structures, Solid-shell, Kirchhoff-Love

This work explores the use of solid-shell elements in the the framework of isogeometric shape optimization of shells. The main difference of these elements with respect to pure shell ones is their volumetric nature which can provide recognized benefits to analyze, for example, structures with non-linear behaviors. From the design point of view, we show that this geometric representation of the thickness is also of great interest since it offers new possibilities: continuous sizing variations can be imposed by modifying the distance between the control points of the outer surfaces. In other words, shape and sizing optimization can be performed in an identical manner. Firstly, we carry out a range of numerical experiments in order to carefully compare the results with the commonly adopted technique based on the Kirchhoff-Love formulation. These studies reveal that both solid-shell and Kirchhoff-Love strategies lead to very similar optimal shapes. Then we apply a bi-step strategy to integrate shape and sizing optimization. We highlight the potential of the proposed approach on a stiffened cylinder where the cross-section along the stiffener is optimized leading to a final design with smooth thickness variations. Finally, we combine the benefits of both Kirchhoff-Love and solid-shell formulations by setting up a multi-model optimization process to efficiently design a roof.

Introduction

Shape optimization is one of the great application of Isogeometric Analysis (IGA) since the latter integrates Computer Aided Design (CAD) and analysis [START_REF] Hughes | Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement[END_REF][START_REF] Cottrell | Isogeometric Analysis: Toward Integration of CAD and FEA, 1st edn[END_REF]. Structural optimization requires a suitable mix of an accurate geometric description and an efficient analysis model. The classical Finite Element Method (FEM) faces few difficulties originated from the geometric approximation inherent to the finite element mesh. Two approaches can be distinguished in the FEM framework: the CAD-based [START_REF] Imam | Three-dimensional shape optimization[END_REF][START_REF] Braibant | Shape optimal design using b-splines[END_REF][START_REF] Ramm | Shape optimization of shell structures[END_REF][START_REF] Lund | Finite element based design sensitivity analysis and optimization[END_REF]) and the FE-based [START_REF] Bletzinger | Optimal shapes of mechanically motivated surfaces[END_REF][START_REF] Le | A gradient-based, parameter-free approach to shape optimization[END_REF][START_REF] Firl | Regularization of shape optimization problems using FE-based parametrization[END_REF][START_REF] Bletzinger | A consistent frame for sensitivity filtering and the vertex assigned morphing of optimal shape[END_REF] approach. The CAD-based approach uses a CAD model to define the shape design and the corresponding shape parameters. The main problem is that the finite element mesh has to be generated at every design update. The FE-based approach consists in using the spatial location of the nodes as the design variables. It leads to a high number of design variables and irregular shapes with possible distortion meshes. Expensive mesh regularization techniques are performed to overcome the shape design representation [START_REF] Firl | Regularization of shape optimization problems using FE-based parametrization[END_REF][START_REF] Bletzinger | A consistent frame for sensitivity filtering and the vertex assigned morphing of optimal shape[END_REF]. Beyond this, the IGA concept fills the needs of embedding efficient geometric and analysis models by discretizing the structure with its intrinsic, computer-aided geometric definition. This study relies on a Non-Uniform Rational B-Splines (NURBS) based framework [START_REF] Cohen | A flexible approach for coupling nurbs patches in rotationless isogeometric analysis of kirchhofflove shells[END_REF][START_REF] Piegl | The NURBS Book (2Nd Ed[END_REF][START_REF] Farin | Curves and Surfaces for CAGD: A Practical Guide[END_REF]. In fact, before the development of IGA, NURBS have been used in shape optimization due to their ability to describe smooth design, and because they offer an attractive way to impose the design updates [START_REF] Imam | Three-dimensional shape optimization[END_REF][START_REF] Braibant | Shape optimal design using b-splines[END_REF][START_REF] Ramm | Shape optimization of shell structures[END_REF][START_REF] Lund | Finite element based design sensitivity analysis and optimization[END_REF]. Finally, since CAD uses the boundary representation (B-rep), IGA is especially suitable for analyz-ing structures whose geometry is easily derived from a surface, as is the case for shells [START_REF] Kiendl | Isogeometric shell analysis with Kirchhoff-Love elements[END_REF][START_REF] Benson | Isogeometric shell analysis: The Reissner-Mindlin shell[END_REF][START_REF] Echter | A hierarchic family of isogeometric shell finite elements[END_REF][START_REF] Dornisch | Isogeometric Reissner-Mindlin shell analysis with exactly calculated director vectors[END_REF][START_REF] Bouclier | Efficient isogeometric NURBS-based solid-shell elements: Mixed formulation and Bmethod[END_REF][START_REF] Caseiro | On the Assumed Natural Strain method to alleviate locking in solid-shell NURBS-based finite elements[END_REF][START_REF] Cardoso | Blending moving least squares techniques with NURBS basis functions for nonlinear isogeometric analysis[END_REF]Bouclier et al 2015b,a;[START_REF] Caseiro | Assumed natural strain NURBS-based solid-shell element for the analysis of large deformation elasto-plastic thinshell structures[END_REF]. Therefore, isogeometric shape optimization of shells is a promising field whose effectiveness has already been highlighted in some recent works [START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF][START_REF] Bandara | Isogeometric shape optimisation of shell structures using multiresolution subdivision surfaces[END_REF].

Isogeometric shape optimization has been successfully applied to a wide range of applications (see for example [START_REF] Wall | Isogeometric structural shape optimization[END_REF][START_REF] Qian | Full analytical sensitivities in NURBS based isogeometric shape optimization[END_REF][START_REF] Nagy | Isogeometric sizing and shape optimisation of beam structures[END_REF][START_REF] Nagy | Isogeometric design of elastic arches for maximum fundamental frequency[END_REF][START_REF] Nguyen | Isogeometric Shape Optimization for Electromagnetic Scattering Problems[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF][START_REF] Taheri | Simultaneous isogeometrical shape and material design of functionally graded structures for optimal eigenfrequencies[END_REF][START_REF] Fußeder | Fundamental aspects of shape optimization in the context of isogeometric analysis[END_REF]Wang and Turteltaub 2015;[START_REF] Ding | Accurate analysis and thickness optimization of tailor rolled blanks based on isogeometric analysis[END_REF][START_REF] Bandara | Isogeometric shape optimisation of shell structures using multiresolution subdivision surfaces[END_REF][START_REF] Herrema | A framework for parametric design optimization using isogeometric analysis[END_REF]Wang et al 2017c,b)). A general procedure, which has been improved over the years [START_REF] Daxini | Parametric shape optimization techniques based on Meshless methods: A review[END_REF], is commonly adopted. It is based on a multilevel design concept which consists in choosing different refinement level of the same NURBS-based geometry to define both optimization and analysis spaces [START_REF] Nagy | Isogeometric sizing and shape optimisation of beam structures[END_REF][START_REF] Nagy | Isogeometric design of elastic arches for maximum fundamental frequency[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF]Wang and Turteltaub 2015). Shape updates are represented by altering the spatial location of the control points, and in some case the weights, on the coarse level. The finer level defines the analysis model and is set to ensure good quality of the solution. The optimization and analysis refinement levels are independently chosen which provides a problem-adapted choice of the spaces. Gradient-based optimization algorithms are generally called upon to solve the isogeometric optimization problems cited therein. In fact, isogeometric shape optimization can provide accurate sensitivities, especially when geometric properties such as curvature are involved. Several works tackle the issue of suitable IGA-based sensitivities which has led to improve the optimization procedures by reducing the discretizationdependency during the shapes updates (Wang et al 2017a; [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF]. From this overview, it can be seen that main researches focus on modeling and optimization aspects leading today to an accurate general procedure of resolution. However, discussion from the analysis point of view is somehow missing. We believe that this is a crucial issue especially regarding the shape optimization of shells.

In the framework of isogeometric shape optimization of shells, it seems that Kirchhoff-Love NURBS elements have been mainly considered [START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF][START_REF] Bandara | Isogeometric shape optimisation of shell structures using multiresolution subdivision surfaces[END_REF], while little research works adopt a Reissner-Mindlin NURBS shell formulation [START_REF] Kang | Isogeometric shape optimization of trimmed shell structures[END_REF]. However, the influence of the shell formulations on the optimization results is not well-known, and, unfortunately, the choice of one shell formulation rather than another is usually dictated by numerical needs. Thus, the Kirchhoff-Love NURBS element is attractive due to its low computational cost, but it should be borne in mind that it is restricted to very thin shells. In order to provide a broader framework that is also suitable for thicker shells, we propose to consider in this work solid-shell NURBS-based element as well. Indeed, we believe that the solid-shell element can bring interesting features in the shell shape optimization framework from both the analysis and the design points of view. First and foremost, it has been demonstrated that such elements are particularly adapted to compute shells with variable shell thicknesses and with complex behaviors (geometric, material and contact nonlinearities) [START_REF] Bouclier | Efficient isogeometric NURBS-based solid-shell elements: Mixed formulation and Bmethod[END_REF][START_REF] Caseiro | On the Assumed Natural Strain method to alleviate locking in solid-shell NURBS-based finite elements[END_REF][START_REF] Cardoso | Blending moving least squares techniques with NURBS basis functions for nonlinear isogeometric analysis[END_REF]Bouclier et al 2015b,a;[START_REF] Caseiro | Assumed natural strain NURBS-based solid-shell element for the analysis of large deformation elasto-plastic thinshell structures[END_REF].

Then, regarding more precisely the design aspect, it appears that the geometric representation of the thickness offers new possibilities. Indeed, the volume representation of the structure enables to optimize the thickness profile continuously by modifying the control point coordinates. In this respect, [START_REF] Ding | Accurate analysis and thickness optimization of tailor rolled blanks based on isogeometric analysis[END_REF] highlight the benefits of making use of a NURBS solid element to accurately represent and analyze the thickness variations encountered in tailor rolled blanks. Thickness optimization becomes thus natural with the solid-shell approach and its inherent continuous transition between the thin and the thick case appears attractive in a large range of applications. For example, automotive industry already uses rolling processes to enable weight reduction and a better use of materials [START_REF] Kopp | Flexibly Rolled Sheet Metal and its use in Sheet Metal Forming[END_REF][START_REF] Merklein | A review on tailored blanks-Production, applications and evaluation[END_REF][START_REF] Hirt | Selected Processes and Modeling Techniques for Rolled Products[END_REF]. Blended composite laminate structures also presents continuous thickness transition which has a fundamental role to play in this context [START_REF] Adams | Genetic algorithm optimization and blending of composite laminates by locally reducing laminate thickness[END_REF]. Regarding the optimization, the number of design variables is reduced in comparison to more classic methods using discrete sizing parameters at the element level. Applying the thickness variations with the NURBS geometry tends to regularize the problem and it acts as a filter. Therefore, no additional filtering techniques are required to avoid undesired phenomenon as, for example, checkerboard pattern experienced with discrete thickness approaches [START_REF] Nha | An evolutionary structural optimization method for sizing problems with discrete design variables[END_REF][START_REF] Dfvlr Lam | Performance characteristics of resizing algorithms for thickness optimization of plate structures[END_REF]. The benefits of solid-shells must, however, be qualified: the automatic generation of analysis-suitable volumetric NURBS models is still a challenging task [START_REF] Liu | Volumetric T-spline construction using Boolean operations[END_REF][START_REF] Akhras | Towards an automatic isogeometric analysis suitable trivariate models generation-Application to geometric parametric analysis[END_REF]. All the same, for slender structures the task is affordable since the geometry is derived from a surface.

Therefore, this study has two main objectives. Firstly, we seek to extend the existing methodology to be able to consider solid-shell elements. Comparisons with Kirchhoff-Love elements are performed in order to validate the proposed solid-shell strategy. Then, we study the benefits of the solid-shell formulation for the shape optimization of variable thickness structures, and its interest in carrying out continuous sizing optimization of shells. Hence, this pa-per is organized as follows. In section 2, we review the concept of isogeometric shape optimization, and we give practical aspects on the design parametrization of shells. In section 3, we remind the basics of the solid-shell and the Kirchhoff-Love NURBS formulations and specify their treatment to be used for isogeometric shape optimization. Section 4 provides optimization examples on which comparisons between the optimal results obtained with both shell formulations are made. In section 5, we apply the proposed solid-shell strategy for thickness optimization of shells. We perform an integrated shape and sizing optimization on a stiffened cylinder. Then, we combine both benefits of Kirchhoff-Love and solid-shells into a multi-model approach. The global optimal shape is obtained through a first optimization by using Kirchhoff-Love shells. Then, the solid-shell formulation is used to further improve the structures by optimizing the thickness variation. Finally, section 6 concludes on this work by summarizing our most important points and motivating future research in this direction.

Isogeometric shape optimization of shells

The framework of this study draws on research dealing with isogeometric shape optimization (see for example [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF][START_REF] Fußeder | Fundamental aspects of shape optimization in the context of isogeometric analysis[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Bandara | Isogeometric shape optimisation of shell structures using multiresolution subdivision surfaces[END_REF]). This section reviews the principal insights of the approach by highlighting main theoretical points and also practical aspects. In particular, we discuss the multilevel optimization process and the definition of design variations in the context of isogeometric shape optimization of shells.

NURBS and Isogeometric Analysis

NURBS are a generalization of B-splines and standard in CAD and computer graphics for geometry modeling (see [START_REF] Cohen | A flexible approach for coupling nurbs patches in rotationless isogeometric analysis of kirchhofflove shells[END_REF][START_REF] Piegl | The NURBS Book (2Nd Ed[END_REF][START_REF] Farin | Curves and Surfaces for CAGD: A Practical Guide[END_REF]). Only the fundamentals are given in the following. For futher details, the interested read is referred to the references cited therein. A general expression for a NURBS geometry with parameter ξ ∈ R d (d being the dimension of the space) is written as:

S S S(ξ ) = n ∑ I=1 R I (ξ )P I , ( 1
)
where n is the number of control points P I and R I are the NURBS basis functions. The multivariate NURBS basis functions R I are obtained as weighted B-splines N I by the following rational definition:

R I (ξ ) = N I (ξ )w I W (ξ ) , with W (ξ ) = n ∑ k=1 N k (ξ )w k , (2)
where w I is the weight of the I th control point. The multivariate B-spline functions N I are defined as the tensor product of univariate B-spline functions. Finally, the 1D B-spline basis functions are piecewise polynomials defined by the polynomial degree p and a set of parametric coordinates ξ i collected into a knot vector Ξ. These B-spline basis functions are obtained recursively using the Cox-de Boor recursion formula (see [START_REF] Cohen | A flexible approach for coupling nurbs patches in rotationless isogeometric analysis of kirchhofflove shells[END_REF]). The concept behind IGA [START_REF] Hughes | Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement[END_REF][START_REF] Cottrell | Isogeometric Analysis: Toward Integration of CAD and FEA, 1st edn[END_REF] consists in discretizing the displacement field using the NURBS parametrization already introduced to describe the geometry:

u u u(ξ ) = n ∑ I=1 R I (ξ )U I , (3)
where U I is the displacement corresponding to the I th control point. By substituting the NURBS approximations in the weak formulation of a boundary value problem, a linear system to be solved is obtained as in standard FEM. An interesting feature of NURBS is their high degree of continuity. If m is the multiplicity of a given knot, the functions are C p-m continuous at that location. This is attractive from both design and analysis points of view. In particular, it allows to define smooth free form shapes with high continuity. Numerical errors in the analysis are reduced because the geometry can be exactly preserved. The higher continuity within a NURBS patch increases the accuracy of the analysis in comparison to standard FEM where C 0 regularity is encountered. Furthermore, NURBS present efficient refinement procedures which allows to enhance the design space without changing the geometry nor the parametrization. In particular, k-refinement, in which order and continuity of the basis functions are simultaneously increased, can be performed. In practice and especially in the structural optimization framework, it is efficient to use the matrix representation of refinement methods of NURBS. For more details on refinement strategies of NURBS and their matrix representation, reference is made to [START_REF] Piegl | The NURBS Book (2Nd Ed[END_REF][START_REF] Lee | Degree elevation of nurbs curves by weighted blossom[END_REF].

Optimization flowchart

Figure 1 describes the main steps of the optimization process. It gets an initial CAD model as a main input. CAD model has here to be understood as an analysis suitable NURBS geometry. As an output the optimal shape is given and, an interesting point to notice is that this shape is a NURBS geometry too. Therefore, the result is directly usable by designers and, thus, for production. There is no need to re-design the final geometry in a CAD environment as it would be the case with FE-based optimization. The optimal NURBS geometry is obtained from an iterative process with three main steps. Each iteration consists in:

1. updating the current shape in a hopefully better one, 2. running the analysis for this new geometry and then inferring mechanical and/or geometrical properties of interest, 3. computing new design variations to further improve the structure.

A more detailed description of each step of the process is given in the following paragraphs.

Multilevel Design

A major asset of isogeometric shape optimization is the possibility to properly choose both optimization and analysis spaces ( It is clear that, the finer the analysis model, the better the results. The choice of the refinement level used for the analysis model is principally dictated by the need for realistic computational cost. When it comes to the optimization model, the situation is quite different. In this case, the choice of the refinement level has an impact on the complexity of the optimal shape reached. A coarse optimization model may provide a simpler optimal shape than a fine optimization model. Therefore, few questions can be asked in order to accurately define the refinement levels for the analysis and optimization models :

• How different from the initial geometry the optimal shape is allowed to be ? • Depending on the intended level of shape variation complexity, how fine should the analysis model be to ensure reliable results ?

These points are illustrated and discussed in more details in section 4 of this paper.

Design variables

The isogeometric framework provides a proper and accurate way to vary the geometry. It has long been established that shape control of NURBS is very comfortable and leads to smooth optimal design [START_REF] Imam | Three-dimensional shape optimization[END_REF][START_REF] Braibant | Shape optimal design using b-splines[END_REF][START_REF] Ramm | Shape optimization of shell structures[END_REF]. Indeed, a straightforward and natural way of modifying the shape of the structure is to move the control points. NURBS also offer an other interesting way to apply shape variation by modifying the control point weights, but this point is not investigated here. As pointed out by [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF], for the optimization of free form shapes it is generally sufficient to use only the control point coordinates. Nevertheless, note that taking the weights as design variables may improve the optimal shape by adding proper shape adjustments if the optimization is performed on a very coarse NURBS model [START_REF] Wall | Isogeometric structural shape optimization[END_REF][START_REF] Nagy | Isogeometric sizing and shape optimisation of beam structures[END_REF][START_REF] Qian | Full analytical sensitivities in NURBS based isogeometric shape optimization[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF].

Although attractive, taking the control point coordinates as design variables requires attention. Problems are similar to those encountered with FE-based optimization where the FE nodal points are directly used as design variables. In fact, attention must be paid to mesh distortion in order to avoid local non-injectivities due to fold-overs. If the control points are allowed to independently move in every spatial direction then the optimization process constantly needs to check the mesh quality and often requires cumbersome mesh regularization techniques [START_REF] Bletzinger | Optimal shapes of mechanically motivated surfaces[END_REF][START_REF] Choi | A mesh regularization scheme to update internal control points for isogeometric shape design optimization[END_REF][START_REF] Le | A gradient-based, parameter-free approach to shape optimization[END_REF][START_REF] Nguyen | Isogeometric Shape Optimization for Electromagnetic Scattering Problems[END_REF][START_REF] Firl | Regularization of shape optimization problems using FE-based parametrization[END_REF]. In the context of isogeometric optimization, links between design variables and move directions are usually set up to avoid these challenges [START_REF] Lee | Dangers inherited in shells optimized with linear assumptions[END_REF][START_REF] Wall | Isogeometric structural shape optimization[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF][START_REF] Taheri | Simultaneous isogeometrical shape and material design of functionally graded structures for optimal eigenfrequencies[END_REF]Wang et al 2017c,b). Figure 3 shows this idea in the case of a hemisphere. By letting only each design control point moving in the radial direc-tion, we avoid mesh regularity problems. Moreover, control points located on patch boundaries are often (but not always) coupled with their neighbors to ensure a G 1 geometric continuity between patches [START_REF] Kiendl | Isogeometric shell analysis with Kirchhoff-Love elements[END_REF]. This prevents the emergence of kinks during the optimization and it finally leads to smoother optimal shapes. Linking design variables may also be relevant to preserve initial shape properties as reflectional, rotational and translational symmetries. It also guarantees a reduction in the number of design variables which leads to lower computational cost and again smoother optimal shapes.

Shape gradient and Normalization

A typical shape optimization problem can be formulated as the minimization of a given objective function carried out over a set of admissible domains. In this work, we limit ourselves to the case of minimizing the compliance. Other objective function could have been considered as the mass, the displacement of a given material point, the maximal stress, or the critical buckling load and so on. The design space is limited by lower bound s l and upper bound s u which surround the n s design variables collected in vector s = (s 1 , s 2 , . . . , s n s ). This unconstrained optimization problem can be defined mathematically as:

min s∈Ω f (s) where: Ω = {s ∈ R n s | s l ≤ s ≤ s u } f (s) = 1 2 F • U = 1 2 U • K • U. (4)
F denotes the force vector, K the stiffness matrix, and U the displacement vector which satisfies the equilibrium equations. Moreover, note that constraints, as for example a maximal mass, can be added to this optimization problem. An illustration regarding a constrained optimization problem is provided in 5.1.

In order to use a gradient-based optimization algorithm, the derivatives of the compliance with respect to the design variables need to be computed. A simple chain rule leads to the following result (see, e.g., [START_REF] Nagy | Isogeometric sizing and shape optimisation of beam structures[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF]; [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF] for more details):

d f ds = U • dF ds - 1 2 U • dK ds • U. ( 5
)
The derivatives of the force vector and stiffness matrix with respect to the design variables are approximated by firstorder finite differences

dF ds ≈ F(s + ∆ s) -F(s) ∆ s , dK ds ≈ K(s + ∆ s) -K(s) ∆ s . ( 6
)
Combining equation ( 6) with equation 5 gives semianalytical sensitivities which are easy to compute and offer practical advantages [START_REF] Lund | Finite element based design sensitivity analysis and optimization[END_REF][START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF]. In fact, the computation of semi-analytical sensitivities does not require more development than those necessary to build the stiffness matrix and the force vector. In addition, whatever the shell formulation, computation of the sensitivities is identical. Note that analytical sensitivities may also be available, depending on the shell formulation, but require more numerical development [START_REF] Cho | Isogeometric shape design optimization: exact geometry and enhanced sensitivity[END_REF][START_REF] Qian | Full analytical sensitivities in NURBS based isogeometric shape optimization[END_REF][START_REF] Nagy | Isogeometric design of anisotropic shells: Optimal form and material distribution[END_REF][START_REF] Taheri | Simultaneous isogeometrical shape and material design of functionally graded structures for optimal eigenfrequencies[END_REF][START_REF] Ha | Generalized isogeometric shape sensitivity analysis in curvilinear coordinate system and shape optimization of shell structures[END_REF].

In practice, the derivatives of the compliance with respect to the design variables are not directly used as the search direction. Indeed, it has been shown [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF] Wang et al 2017a) that normalization approaches help to reduce discretization-dependency. As a results, we use here the diagonally-lumped mapping matrix normalization, also called sensitivity weighting, which leads to the following search directions:

d i = - 1 M ii d f ds i , with M ii = ∑ k∈Ω i D R k dD . (7)
Geometric factor M ii can be interpreted as the fraction of surface or volume linked to the i th design variable. It integrates in the physical space D all NURBS basis functions related to the i th design variable (Ω i lists the indexes of the control points governed by the i th design variable). Finally, this normalized semi-analytical gradient can be used in gradient based optimization algorithms as, for example, those provided in the SciPy package [START_REF] Jones | SciPy: Open source scientific tools for Python[END_REF]. In this work, the SLSQP method is used [START_REF] Kraft | A Software Package for Sequential Quadratic Programming[END_REF]).

Isogeometric shell analysis

In this work, two isogeometric shell formulations are envisaged to perform the optimization. On the one hand, we extend the shell shape optimization framework in order to be able to use solid-shell elements for the modeling. On the other hand, we consider the already established optimization process based on the Kirchhoff-Love elements in order to compare its results with those provided by the solid-shell element-based procedure. The two shell formulations, along with their treatment to be used for isogeometric shape optimization, are introduced in this section.

3.1 NURBS-based solid-shell element

Basics

The solid-shell NURBS element used in this work is the classical one introduced by [START_REF] Bouclier | Efficient isogeometric NURBS-based solid-shell elements: Mixed formulation and Bmethod[END_REF]. It is based on a 3D continuum mechanics-type formulation. The idea is to take a 3D solid continuum element and, in order to limit the computational costs, discretize the shell using a single layer of elements through the thickness. This element is used in its raw form and no particular treatment (strainprojection, reduced integration, etc.) is put in place here. To alleviate locking, the basic strategy of increasing the polynomial degree of the analysis model is performed. Despite the simplicity of the formulation, this element has shown similar behaviours to those of Kirchhoff-Love or of Reissner-Mindlin NURBS shell elements through several shell benchmarks (Bouclier et Since based on classical 3D continuum elements, the standard solid-shell formulation is well-know. We only remind the basics in order to underline differences with the Kirchhoff-Love NURBS shell formulation. No particular kinematic assumptions are made on the displacement field u u u. Therefore, strain and stress fields have six components using Voigt notation. In the context of small perturbations, the components of the infinitesimal strain tensor ε ε ε are:

ε ε ε i j = 1 2 ∂ u i ∂ x j + ∂ u j ∂ x i ∀ i, j = {1, . . . , 3} 2 , ( 8
)
where u i denotes the i th component of the displacement at point M located at (x 1 , x 2 , x 3 ) in the Cartesian frame. In IGA, the displacement field is approximated using the NURBS functions which enables to build the geometry of the structure. This implies an approximation ε ε ε h of the strains as fol- lows:

ε ε ε h = ∑ A B A U A with B A =            R A,1 0 0 0 R A,2 0 0 0 R A,3 0 R A,3 R A,2 R A,3 0 R A,1 R A,2 R A,1 0            , (9)
where U A is the vector of the control variables at control point A, and R A,i denotes the derivative of the NURBS function at control point A with respect to coordinate x i . Introducing the approximation fields into the classical weak formulation leads to the following linear system to be solved:

K • U = F where K ab = D B a • C • B b dV, F a = D R a f f f dV + Γ F R a F F F dγ, (10)
and with C the Hooke tensor, f f f as body force in D, and F F F as surface load over the Neumann boundary Γ F . It has been shown that for solid-shell element it is unnecessary to increase the degree of approximation of the functions through the thickness beyond 2. When speaking of solid-shell elements of degree p, the reader should understand that degree of p is used in the main directions (commonly referred to as ξ and η) and that degree 2 is selected trough the thickness of the shell. Finally, this basic solid-shell NURBS element is nothing other than the classical IGA solid element. Consequently, being able to perform shape optimization of shell with this element offers another attractive interest: with the simplest IGA code in hand, one can run optimization without any further development.

Imposing shape variations

With solid-shell elements, the structure is described by its volume and not by its mid-surface. Thus, the standard strategy depicted in section 2.4 needs to be further improved: the shape control has to deal with volume. A solid-shell is characterized by two outer skins. A natural idea is to define one of the outer skin as the master surface and the other as the slave surface. This is carried out in this work in two steps. A first set of design variables updates the master surface as it is commonly done with shell structure and as already discussed. Then, a second set of design variables moves the slave surface starting from the master surface. This idea is illustrated in figure 4 which shows the shape parametrization of a tube.

Each control point from the master surface has its equivalent on the slave surface. If a variable var1 moves the master point from initial position P 0 m to P m in direction n 1 and a second variable var2 further updates the neighbour slave point form P 0 s to P s in direction n 2 , then the design update follows

P m = P 0 m + var1 n 1 , [step1] P s = P 0 s + var1 n 1 + var2 n 2 . [step2] (11)
By setting additional bounds to var2 we avoid undesirable shapes as for example overlapping of the outer skins. Generally, in shape optimization of shells, the thickness is kept constant. For this purpose, we apply a simple approach where the second set of design variables is not included in the optimization. In other words, same shape variations are applied on both outer skins. If the optimization model has a higher degree through the thickness than 1, then the intermediate control points are interpolated between both outer skin. Finally, the use of solid-shell NURBS element offers the possibility to optimize through the thickness. The thickness optimization issue seems to be promising [START_REF] Ding | Accurate analysis and thickness optimization of tailor rolled blanks based on isogeometric analysis[END_REF], and it is investigated in the present work (see section 5).

3.2 Kirchhoff-Love NURBS shell formulation

Basics

At the beginning of NURBS-based shell formulations, structural model approach was firstly developed (see, e.g., the initial works of [START_REF] Kiendl | Isogeometric shell analysis with Kirchhoff-Love elements[END_REF] In the Kirchhoff-Love NURBS shell formulation, the structure is defined by its mid-surface D and its thickness t. In the IGA concept, the mid-surface is described by a NURBS function S providing at each material point x ∈ D a straightforward parametrization in terms of a system of curvilinear coordinates as

x(ξ , η, ζ ) = S(ξ , η) + ζ a 3 (ξ , η), - t 2 ≤ ζ ≤ t 2 . ( 12
)
The unit normal a 3 and the standard covariant vectors of the mid-surface a α are given by

a 1 = ∂ S ∂ ξ = S ,ξ , a 2 = ∂ S ∂ η = S ,η , a 3 = a 1 × a 2 |a 1 × a 2 | . (13)
The displacement at each material point on the whole structure is directly linked to the displacement on the mid-surface ū u u. This is the result of the well-know Kirchhoff kinematic assumptions: straight lines normal to the mid-surface are characterized by rigid-body displacements and remain normal to the mid-surface after deformation. It follows that the linearized strain tensor of the shell is found to be of the form

ε ε ε αβ = e e e αβ + ζ κ κ κ αβ , (14)
where membrane strains e e e and bending strains κ κ κ are given by

e e e αβ = 1 2 ū u u ,α • a β + ū u u ,β • a α , (15)
κ κ κ αβ = -ū u u ,αβ • a 3 (16) + 1 |a 1 × a 2 | ū u u ,ξ • a α,β × a 2 + ū u u ,η • a 1 × a α,β + a 3 • a α,β |a 1 × a 2 | ū u u ,ξ • a 2 × a 3 + ū u u ,η • a 3 × a 1 .
Greek indices take the values 1 and 2. In fact, the kinematic assumptions vanish the transverse shear strains ε ε ε α3 = 0. As for the mid-surface, the displacement field is approximated using the NURBS basis functions. Introducing the approximation ūh into (15) and ( 16) gives the membrane and bending strains as follows:

e e e( ūh

) = n ∑ i=1 B m i Ūi , and κ κ κ( ūh ) = n ∑ i=1 B f i Ūi , (17)
where the matrices B m i and B f i can be found in A.1. The equilibrium configurations of the shell follow from the principle of minimum potential energy which, once expressed in weak form, leads to the linear system K • Ū = F where:

Kab = D t B m a • H • B m b + t 3 12 B f a • H • B f b dA, Fa = D R a p p p dA + Γ t R a t t t dγ, ( 18
)
and with p p p as distributed loads per unit area of D, t t t as axial forces per unit length of border Γ t , and H the material ten- sor describing the linear elastic behaviour of the shell. Its expression is given in A.2.

Multi-patch analysis

By comparing the basic solid-shell NURBS element with the Kirchhoff-Love NURBS shell element, it is clear that the first is simpler. In addition to this practical aspect, Kirchhoff-Love element faces some difficulties when it comes to enforce the G 1 -continuity between patches, fix the angle between surface folds, enforce symmetry conditions and prescribe rotational Dirichlet boundary conditions. Computational formulations to overcome these difficulties are usually either expensive or difficult to implement and this issue is still of interest [START_REF] Apostolatos | Domain decomposition methods and kirchhoff-love shell multipatch coupling in isogeometric analysis[END_REF][START_REF] Guo | A parameter-free variational coupling approach for trimmed isogeometric thin shells[END_REF][START_REF] Duong | A new rotation-free isogeometric thin shell formulation and a corresponding continuity constraint for patch boundaries[END_REF]Coox et al 2017).

In this paper the bending-strip approach (Kiendl et al 2010;

Goyal and Simeon 2017) is used to handle multi-patch discretization and thus, in particular, to perform the coupling between stiffeners and panels when a stiffened structure is investigated. The idea is to introduce additional patches of fictitious material, namely the bending strips, where the patches are joined with C 0 continuity. The material has zero membrane stiffness and non-zero bending stiffness only in the direction transverse to the strip. In practice, adding the bending strips consists in adding terms to the global stiffness matrix introduced in 18. The stiffness matrix of a bending strip is the same as the one of a classical Kirchhoff-Love NURBS patch, except the material tensors: The membrane material tensor is null and the bending material tensor H bnd-strip is given by

Kbnd-strip ab = D ✘ ✘ ✘ ✘ ✘ ✘ t B m a • H • B m b + t 3 12 B f a • H bnd-strip • B f b dA.
H bnd-strip = P T • H bnd-strip •P, H bnd-strip =   E 0 0 0 0 0 0 0 0   (19)
and with E as the directional bending stiffness. Matrix H bnd-strip describes the stress-strain relationship in the local coordinate system defined by:

e 1 = a 1 |a 1 | , e 2 = a 2 -(a 2 • e 1 ) e 1 |a 2 -(a 2 • e 1 ) e 1 | , e 3 = a 3 . ( 20
)
Finally, matrix P maps the components of the strain tensor with respect to the contravariant basis vectors to those with respect to local orthonormal basis e i .

Investigation on preliminary optimization examples

Description of the numerical examples

In order to investigate the use of isogeometric solid-shell elements for shape optimization of shell structures, two optimization examples depicted in figure 5 are first studied. The motivation is to ensure compatibility of solid-shell NURBS elements with isogeometric shape optimization. The effectiveness of the isogeometric Kirchhoff-Love shell formulation for shape optimization of thin structures no longer needs to be proven. Therefore, there is a great interest in comparing the optimum results obtained by solid-shell elements to those obtained by the Kirchhoff-Love formulation.

Both the pinched cylinder and the pinched hemisphere are known from the shell obstacle course (Belytschko et al 1985). These optimization problems are non-convex and are numerically challenging. In fact, the pinched hemisphere and the pinched cylinder are commonly studied to test highefficient shell elements. The hemisphere is challenging in terms of locking: it exhibits almost no membrane strains. The pinched cylinder is even more severe. The loading and the displacement boundary conditions lead to a highly localized strain state. The local loading associated with the curved geometry make the problem challenging numerically in terms of both locking and representation of complex membrane states. Optimizing this set of problems automatically requires accurate analysis results. Approximation error of sensitivity analysis can lead to different local minimum and unfortunately lead to poor optimal geometry [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF].

Pinched cylinder

We first deal with the test case of the pinched cylinder (see figure 5(a)). The optimization model is chosen relatively coarse in order to keep the shape simple. Because of the complexity of the strain state generated by the local loading, a fine optimization model would lead to a too complex final shape. Due to the symmetry of the problem, one eighth of the cylinder is considered, hence only one NURBS patch is needed. In the framework of the results presented in figure 6, the optimization model is built with uniform knot vectors which discretize the patch in 3 and 4 elements in the circumferential and the axial direction respectively. Degree 2 is chosen in both directions. With such an optimization model in hand, 6 design variables update the geometry as depicted in figure 6(a). Variables named x i move some control points in the x-direction, whereas variables denoted r i move some control points in the radial direction. At the patch boundaries, groups of two or four control points are set in order to preserve a smooth G 1 -continuous shape as explained in figure 3. All variables are constrained by bounds of [0,60]. Analysis is performed with cubic elements. Refinement level of 4 and 3 are set in the circumferential and the axial direction respectively, leading to an analysis model of 1536 elements.

We run the shape optimization with Kirchhoff-Love NURBS shell elements and with Solid-shell NURBS elements. Final shapes are identical as shown in figure 6. The difference between both results is not visible when looking at the design modifications. To further compare the results, the distance between the mid-surfaces is presented. The maximum gap between the optimal surface obtained by Kirchhoff-Love elements and the final mid-surface obtained by Solid-shell elements is about 0.5 which is, by comparison with the range of allowed design variation, very small. The maximum gap normalized by 60 is lower than 1%. Ta- The small difference in final design when using either Kirchhoff-Love elements or Solid-shell elements may be a result of the difference in converged solutions. It has be shown that the converged solution using NURBS solidshell elements for the pinched cylinder problem is a little softer than the reference solution [START_REF] Bouclier | Efficient isogeometric NURBS-based solid-shell elements: Mixed formulation and Bmethod[END_REF][START_REF] Hughes | Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement[END_REF]. This little difference in solution may influence the shape updates and finally lead to a slightly different optimal design. Beyond that, the shape optimization of the pinched cylinder reveals very similar behavior for both shell formulations. During the optimization process, the shape updates are similar and it has been observed that the convergence is reached after the same amount of iterations.

Pinched hemisphere

The second example of the pinched hemisphere helps us to further investigate the use of solid-shell NURBS elements for shape optimization of shells. Note that the same test case has been conducted by [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF] in the context of isogeometric shape optimization. The symmetry of the problem allows to consider only one quarter of the structure. Hence, only one NURBS patch is needed. The optimization model is built with uniform knot vectors which discretize the geometry in 5 quadratic elements per parametric direction. 25 design variables depicted in figure 7(a) are defined in order to prescribe the shape modifications in radial direction. Once again, the boundary control points are coupled in order the preserve a smooth G 1 -continuous shape. All variables are constrained by bounds of [-1, 1]. Further k-refinement is applied to perform the analysis, leading to an analysis model with 400 cubic elements. Tests showed that this refinement level was adequate to converge to the good looking optimal design.

Once again, we run the optimization with the Solidshell NURBS elements and compare the results with those provided by Kirchhoff-Love NURBS elements. Figure 7(b) presents the optimal design obtained by both formulations. Despite the complexity of this optimization problem, the final shapes are very similar: no visible difference can be observed. Similar hollows and bumps are located on both optimal designs. To compare more in detail the results, the distance between the two mid-surfaces is computed (see figure 7(c)). The maximum gap in shape is about 6.34 • 10 -3 which is very small in comparison with the allowed range of design modification. Normalized by 2, the maximum relative distance between the optimal mid-surfaces is lower than 0.5%. Final values of the cost function also demonstrate the similarity between the optimal design obtained with the solid-shell NURBS elements and the one obtained with the Kirchhoff-Love NURBS elements. The final compliance is 1.34 • 10 -3 for the solid-shell case and 1.30 • 10 -3 for the Kirchhoff-Love case. However, the final shape obtained by Kirchhoff-Love NURBS elements is not better. Note that in order to correctly compare the final shapes, the compliance has to be evaluated on the same analysis model [START_REF] Kiendl | Isogeometric shape optimization of shells using semi-analytical sensitivity analysis and sensitivity weighting[END_REF]. Strain energy is evaluated for the optimized shapes in a post-processing step using a fine mesh with quartic Kirchhoff-Love NURBS elements. This postprocessing procedure reveals that both shapes give a final compliance of 1.34 • 10 -3 which is 1.5% of the initial value 9.24 • 10 -2 .

Beyond the similarity in final shape, the whole optimization histories look alike. Figure 8 shows the convergence of the objective function in both cases either with solidshell NURBS elements or with Kirchhoff-Love NURBS elements. Optimization processes converge at the same speed and stop after an equivalent number of iterations within a given tolerance on function value. Some of the intermediary shapes are also shown in figure 8. Shape updates are very similar and the different hollows and bumps appear at the same iteration during the optimization.

Application to shape optimization of structures with variable thickness

The volume representation of the structure offers the possibility to optimize the thickness profile continuously. This interesting feature is investigated through two examples. First, simultaneous shape and thickness optimization is performed on a stiffened cylinder. Secondly, we present a multi-model approach to combine benefits of both Kirchhoff-Love and solid-shell formulations.

Stiffened cylinder: shape and sizing optimization

Fig. 9 The stiffened cylinder: two stiffeners are added to the initial pinched cylinder test case. The left side depicts a global 3D view of the problem. The right side shows a top view of one quarter of the stiffener and its cross-section for several angles (Notes: thickness is magnified by a factor of 4 for better visualization; this detailed view of the stiffener will be useful to depict the optimal designs in figure 11).

The test case of the stiffened cylinder on which both sizing and shape optimization will be performed is derived from the already discussed pinched cylinder problem (see figures 5 and 6). Two stiffeners are added to the cylinder as depicted on figure 9. The stiffeners are located at heights 150 and 450. Initially, these stiffeners have a constant rectangular cross-section of height 30 and thickness 3. The other geometric and material parameters are detailed in figure 5. The optimization problem consists in modifying the cross section along the stiffeners in order to improve the global stiffness of the structure. A volume constrain is set in order to keep the final volume lower or equal to the initial one V 0 . Due to the symmetry of the problem, only an eighth of the structure is considered.

The optimization procedure is done here in two steps. Height and thickness of the stiffener are optimized oneby-one. Since the thickness is much less sensitive than the height, combining both quantities in a single optimization problem may lead to mostly modify the height of the stiffener. From here on, the shape step denotes the change in height and the sizing step denotes the change in thickness. A first set of design variables h h h updates the height of the stiffener as depicted on figure 10(a). A second set of design variables t t t modifies the thickness as shown in figure 10(b). Obviously, the design parametrization is not unique. Depending on how complex the final shape is allowed to be, a simpler or finer parametrization can be chosen. However, this example shows that using solid-shell NURBS based elements offer an interesting way to optimize through the thickness. Smooth and continuous thickness variations can be imposed. In this work, the optimization problems behind the shape step and the sizing step are different. The shape step consists in minimizing the compliance C with a maximal volume constrain. An additional bound is set in order to limit the height of the stiffener h in the range [25,100]. The sizing step minimizes the volume V under the constrain that the final compliance is lower or equal to the one obtained from the shape step C * . The design space for the thickness t is limited to [0.5, 3.0]. The global optimization process is made of sequences of a shape step followed by a sizing step. The fixed point process iterates until the global convergence is reached; that is, convergence on the design variables, the compliance and the volume. The k th sequence takes the following form:

Minimize C such that 25 ≤ h ≤ 100 V ≤ V 0 h h h k , t t t k h h h k+1 , t t t k C * 1. Shape Step (21) Minimize V such that 0.5 ≤ t ≤ 3.0 C ≤ C * h h h k+1 , t t t k+1 h h h k+1 , t t t k C * 2. Sizing Step (22)
Both optimization problems are complementary. The shape step will use all the available material in order to minimize the compliance. The sizing step rearranges the material distribution through the thickness in order to save material. Then, this saving material can be further use in the shape step in order to continue minimizing the compliance. This sequential optimization strategy is sensitive to the initial starting point and may lead to local optimum. The height and the thickness design variables are decoupled. Thus, to ensure a good quality result, we performed a multi-start algorithm. On top of that, we notice that multidisciplinary optimization strategies could be of interest to properly solve this problem and we refer the interested reader to [START_REF] Balesdent | A survey of multidisciplinary design optimization methods in launch vehicle design[END_REF] Fig. 12 Optimization history of the bi-step approach applied to the stiffened cylinder problem. Three global iterations with successive shape and sizing steps are preformed until convergence. Red zones denote the shape steps where the compliance is minimized and blue zones denote the sizing step where the volume is minimized. [START_REF] Martins | Multidisciplinary Design Optimiza-A Survey of Architectures[END_REF] for futher details regarding this point.

The results of this study are given in figure 11. More precisely, figure 11(a) presents the optimal design when only the height of the stiffener is modified and figure 11(b) depicts the final design when both shape and sizing optimizations are performed. Optimizing through the thickness leads to a better design in the sense of the compliance. It helps to make a better use of the material volume. By reducing the thickness in some locations, the height of the stiffener can be increased in others. The initial compliance C 0 in the case of a straight stiffener is about 1.73e-06. When only shape optimization is employed, the final compliance C 1 is about 1.68e-06. The final compliance C 2 when both shape and sizing optimizations are performed is about 1.64e-06. The compliance gain can be defined as G i = 1 -C i /C 0 . Therefore, the compliance gain G 2 when adding the sizing optimization is almost twice the compliance gain G 1 when only the height of the stiffener is taken into the optimization process. Figure 12 depicts the optimization history. It takes three global iterations to converge. The main design changes are imposed during the first stage (i.e. the first succession of a shape and a sizing step) and the second shape step. For this optimization case, these steps may be enough to get an optimal geometry. Finally, this example shows that the solid-shell NURBS element offers an interesting way to deal with sizing optimization. It extends all the advantages of the NURBS parametrization to sizing optimization problem.

Multi-model optimization

In this last example, we present a multi-model optimization. The idea is to combine the potential of both Kirchhoff-Love and Solid-Shell formulations. In a first stage, the Kirchhoff-Love formulation is used to perform a shape optimization in which the shell structure is varied in the out-of-plane direction. Then, the optimal surface is used to generate a volume model of the structure. In the second stage, thickness optimization is done on the volume model by using the developed solid-shell approach. We notice that this multi-model optimization appears consistent since it corresponds to the design process of structures. The first stage is the prototyping in which the influence of many parameters is studied. It requires a model with low computational cost in or- der to perform many simulations. Once the global shape of the structure is obtained, an high-fidelity volume model is created. Instead of directly using this model to generate the CAD drawings, one can perform a last optimization to apply final adjustments. Thank to the isogeometric approach, the final geometry can directly be transferred to the manufacturing step since no conversion in a CAD format is required.

We show the reliability of the procedure on the following example. Starting from a roof made of a square plate, the objective is to optimize the shape and the thickness of the structure by minimizing the compliance under a maximal volume constraint. This example is based on the work of [START_REF] Kegl | Shape optimization of truss-stiffened shell structures with variable thickness[END_REF] and similar studies can be found for example in [START_REF] Bletzinger | Computational methods for form finding and optimization of shells and membranes[END_REF]; [START_REF] Ikeya | Multi-objective free-form optimization for shape and thickness of shell structures with composite materials[END_REF]. The mechanical parameters of the study are the following: Young modulus E = 210e9, Poisson ratio ν = 0.30, length L = 10, thickness h = 0.10, uniform snow load P = 500. The four corners of the plate are fixed. For the first step, the optimization model is defined by a single quadratic Nurbs patch which discretize the whole structure in 4-by-4 NURBS elements. The surface is parametrized with 32 design variables; the 4 control points located at the corners are kept fixed. A refinement level of 3 is applied in both parametric directions, resulting in an analysis model with 1024 quadratic elements. The volume is constrained to be lower than 110% of the initial volume. Figure 13 shows the final shape. The SLSQP optimizer converges after approximately 40 iterations. The compliance is drastically reduced with a factor of 1.44e-03. The final shape has 4 symmetric planes defined by the origin and the normal vectors X, Y , X+Y and X-Y . This is expected since the problem presents these symmetries. Therefore, some design variables have identical optimal values. The fact that we obtain the symmetries without setting groups of design variables is a meaningful indication to validate the result. Moreover, our result is similar to the one of [START_REF] Kegl | Shape optimization of truss-stiffened shell structures with variable thickness[END_REF] when the control points are only modified in Z-direction (case B in section 5.2 of the cited paper). The optimal positions of the control points are given on the top view and on the side view of the structure on figure 13. This global shape could have been obtained using the solid-shell model as in the previous examples of this paper.

Once the global shape of the structure is obtained, one can build a volume model. Here, we create it by offsetting the previous surface in Z-direction with a distance of h = 0.10. This time, only one quarter of the structure is considered. In order to enrich the design model, knot insertion is performed. The in-plane parametric directions are defined by the same knot vector U = {0 0 0 1 6 2 6 1 2 3 4 1 1 1}. Degree one is chosen in the thickness direction for building the optimization model. Thus, the optimization model counts 25 elements. The design variables modify the distance between the upper and lower surfaces. A total of 36 variables are defined. They parametrize the distance between neighboring control points of both upper and lower surfaces. The design variables are bounded so that the distance between the neighboring control points varies in range t min = 0.025 and t max = 0.20. Since the loading is applied on the upper surface, we keep this surface unchanged. Only the lower control points are moving. The control points located on the symmetric planes are coupled in order to preserve the geometric continuity. Figure 15(a) gives an idea of the initial configuration of the optimization model. A refinement level of two is performed to build the analysis model: it discretizes the structure with 20-by-20 quadratic elements. Figure 15(a) shows the initial displacement field computed with this discretization.

The optimization problem is similar as for the first step. The compliance is minimized under the constraint of keeping the volume lower than the one if the initial configuration. The optimization history is presented in figure 14. The final ratio of the compliance is below 0.55. Therefore, the thickness optimization drastically improves the behaviour of the structure. This gain in stiffness can be observed on the displacement field. Figure 15 depicts the optimal geometry and the displacement field in Z-direction. The overall displacement is reduced in comparison with the initial geometry. Some shape updates are also given in figure 14 where the colormap represents the thickness of the structure. It can be observed that the thickness distribution of the final geometry is quite complex. In case of discrete thicknesses, being able to accurately represent such a distribution may require a high number of elements. Thus, the number of design variables would be higher. Finally, the optimal design (see figure 15) has regions with thin thickness and regions with thick thickness. At the corners, the thickness is equal to 0.20 which is, regarding the slenderness of the structures, relatively important. Using Kirchhoff-Love elements would not be relevant to analyze such a structure.

Conclusion

In the context of isogeometric shape optimization of shells, very little work tackle the influence of the chosen shell formulation. We carry out a range of numerical experiments to compare the results of the proposed solid-shell approach with the commonly adopted technique based on the Kirchhoff-Love formulation. Similar results have been observed for the two shell formulations in terms of final optimized design and in terms of global behavior during the optimization process (convergence speed, intermediary shapes). As a result, this first study accounts for the performance of the developed solid-shell approach for standard shape optimization.

Unlike pure shell elements, we emphasize that the solidshell element used in this work is based on a 3D solid continuum element. With any simple IGA code in hand, one can perform shape optimization of shells without almost any additional developments by imposing identical shape updates on both inner and outer surfaces. This approach leads to identical optimization results as surface shell formulations. In addition, since solid-shell elements are based on a 3D solid continuum element, analytical sensitivities should reasonably be computed even for sophisticated structural response [START_REF] Cho | Isogeometric shape design optimization: exact geometry and enhanced sensitivity[END_REF][START_REF] Qian | Full analytical sensitivities in NURBS based isogeometric shape optimization[END_REF][START_REF] Ha | Generalized isogeometric shape sensitivity analysis in curvilinear coordinate system and shape optimization of shell structures[END_REF][START_REF] Radau | Sensitivity of structural response in context of linear and non-linear buckling analysis with solid shell finite elements[END_REF], even though this point has not been addressed in the present paper.

The major interest of solid-shell elements concerns the geometric representation of the thickness. This is the most important point according to us since it offers new possibilities from the design point of view. Continuous and smooth thickness variations are imposed by modifying the control point coordinates of one outer surface. Thanks to NURBS functions, very few design variables are required to describe quite complex geometries. The direct link with the CAD format provided by the NURBS functions is also highly attractive in this context. The volume representation of the structure makes the optimal shape directly available for CAD drawings, and thus the transfer to the manufacturing step is simplified.

A On the Kirchhoff-Love NURBS shell formulation

A.1 Membrane and bending strain matrices

The membrane and bending matrices introduced in 17 are respectively given by

B m i =     R i,ξ a T 1 R i,η a T 2 R i,η a T 1 + R i,ξ a T 2     (23)
and

B f i =     B i 1 • e 1 B i 1 • e 2 B i 1 • e 3 B i 2 • e 1 B i 2 • e 2 B i 2 • e 3 2 B i 3 • e 1 2 B i 3 • e 2 2 B i 3 • e 3     . (24)
In 23 and 24, (e 1 , e 2 , e 3 ) are the basis vectors of the global Cartesian coordinates system, and the quantities B i 1 , B i 2 and B i 3 are given by

B i 1 = -R i,ξ ξ a 3 + 1 √ a R i,ξ a 1,ξ × a 2 + R i,η a 1 × a 1,ξ + a 3 • a 1,ξ R i,ξ a 2 × a 3 + R i,η a 3 × a 1 , B i 2 = -R i,ηη a 3 + 1 √ a R i,ξ a 2,η × a 2 + R i,η a 1 × a 2,η + a 3 • a 2,η R i,ξ a 2 × a 3 + R i,η a 3 × a 1 , B i 3 = -R i,ξ η a 3 + 1 √ a R i,ξ a 1,η × a 2 + R i,η a 1 × a 1,η + a 3 • a 1,η R i,ξ a 2 × a 3 + R i,η a 3 × a 1 , (25)
with √ a = |a 1 × a 2 |.

A.2 Material tensor

The material tensor within the Voigt formalism is in the form

H = E 1 -ν 2   

Fig. 1

 1 Fig. 1 Isogeometric shape optimization flowchart: overview of the main steps of the process.

Fig. 2

 2 Fig.2Multi-level design approach: optimization and analysis spaces describe the exact same geometry and are initially obtained through different refinement levels of the CAD model (above: NURBS elements, below: control points).

Fig. 3

 3 Fig. 3 Illustration of the design variations of the hemisphere. Middle control points move in radial directions (top) and boundary control points are coupled with their neighbors to ensure G 1 -continuity between patches (bottom).

Fig. 4

 4 Fig. 4 Definition of the design variables in the context of solid-shell isogeometric elements. (a) A first set of design variables updates the master surface and (b) the slave surface follows the shape variation through a second set of variables.

Fig. 5

 5 Fig. 5 Numerical examples to investigate the use of isogeometric solid-shell elements for shape optimization of shell structures: (a) pinched cylinder and (b) pinched hemisphere.

Fig. 6

 6 Fig. 6 Shape optimization of the pinched cylinder: (a) 6 design variables update the optimization model, (b) final shapes obtained by both shell formulations, and (c) comparison of the final shapes depicted as the distance between their mid-surface.

Fig. 7

 7 Fig. 7 Shape optimization of the pinched hemisphere: (a) 25 design variables update the optimization model, (b) final shapes obtained by both shell formulations, and (c) comparison of the final shapes depicted as the distance between their mid-surface.

Fig. 8

 8 Fig. 8 History of compliance and shape during the optimization of the pinched hemisphere when using Kirchhoff-Love NURBS shell elements and NURBS-based Solid-shell elements.

  Fig. 10 Shape and sizing optimization with solid-shell NURBS based elements: the cross-section of the stiffener can be geometrically parametrized and thus optimized. Here, for each section of control points, two variables var1 (a) and var2 (b) impose local design modifications in the height and in the thickness directions respectively

  Fig. 11 Optimization results of the stiffened cylinder: (a) final shape when only the stiffener height is modified, and (b) final shape when both height and thickness are optimized. From left to right are depicted a general 3D view of the final shapes, a top view and several cross sections of the final stiffeners, and the displacement fields for these two optimized designs (annotations on figure 9 can be helpful to understand the detailed views of the stiffener).

Fig. 13

 13 Fig. 13 Optimal shape of a square plate subjected to an uniform snow load: (a) 3D view, (b) top view and (c) side view of the final structure. Red points represent the control points.

Fig. 14

 14 Fig. 14 Variable thickness roof: evolution of the compliance, the volume and the thickness distribution during the optimization.

Fig. 15

 15 Fig. 15 Volume models of the roof: a) initial geometry and initial Z-displacement field, b) optimal geometry and its corresponding Zdisplacement field.

  on the Kirchhoff-Love theory and[START_REF] Benson | Isogeometric shell analysis: The Reissner-Mindlin shell[END_REF] on the Reissner-Mindlin theory, and then, more recently,[START_REF] Echter | A hierarchic family of isogeometric shell finite elements[END_REF] and[START_REF] Dornisch | Isogeometric Reissner-Mindlin shell analysis with exactly calculated director vectors[END_REF]). The structural model approach is based on a discretization of the mean surface alone. It leads to fewer degrees of freedom per element in comparison with solid-shell element. This is particularly true with the Kirchhoff-Love formulation where no rotational degrees of freedom are needed. In return, Kirchhoff-Love shells require C 1 continuity but thanks to NURBS surfaces, this condition is easily satisfied. IGA based on NURBS surfaces offers C p-1 continuity throughout the element. Finally, one should keep in mind that Kirchhof-Love shells neglect transverse shear deformations, and thus they are accurate for analyzing thin structures only.

Table 1

 1 Optimal values of the design variables for the shape optimization of the pinched cylinder.

		x 1	x 2	x 3	r 1	r 2	r 3
	Kirchhoff-Love	60.0	22.65	7.52	0.0	1.74	0.23
	Solid-shell	60.0	23.31	7.58	0.0	1.72	0.51
	Relative Gap [%]	0.0	1.10	0.10	0.0	0.02	0.47

ble 1 gives the optimal values of the design variables for both shell formulations. Results are similar and the relative gap (defined as |s KLs Solid |/60) for each variable does not exceed 1%.

  where: h 12 = νa 11 a 22 + (1 -ν)a 12 a 12 h 33 = 1 2 (1 -ν)a 11 a 22 + (1 + ν)a 12 a 12 (26)with E as the Young's modulus, ν as the Poisson's ratio, and a αβ = a α • a β as the contravariant metric. The corresponding contravariant base vectors a α are defined trough the relation a α • a β = δ αβ where δ αβ is the Kronecker delta. neering 197(33-40):2976-2988 Wang ZP, Turteltaub S (2015) Isogeometric shape optimization for quasi-static processes. International Journal for Numerical Methods in Engineering 104(5):347-371 Wang ZP, Abdalla M, Turteltaub S (2017a) Normalization approaches for the descent search direction in isogeometric shape optimization. CAD Computer Aided Design 82(June):68-78 Wang ZP, Poh LH, Dirrenberger J, Zhu Y, Forest S (2017b) Isogeometric shape optimization of smoothed petal auxetic structures via computational periodic homogenization. Computer Methods in Applied Mechanics and Engineering 323:250-271 Wang ZP, Turteltaub S, Abdalla M (2017c) Shape optimization and optimal control for transient heat conduction problems using an isogeometric approach. Computers and Structures 185:59-74

		a 11 a 11 *	h 12 a 22 a 22	a 11 a 12 a 22 a 12	   
		*	*	h 33


A multi-dimensional percolation approach to characterize sustainable mega-city regions The structure of road networks both translates its past growth dynamics and has a significant impact on the sustainability of territories it irrigates. A method to characterize topologies of these spatial networks is network percolation. Such approaches have been applied to the modeling of urban growth [START_REF] Makse | Modeling urban growth patterns with correlated percolation[END_REF] and to the analysis of street networks for example to extract endogenous urban regions [START_REF] Arcaute | Cities and regions in britain through hierarchical percolation[END_REF] or to characterize the spatial morphology of point patterns [START_REF] Huynh | Characterisation and comparison of spatial patterns in urban systems: A case study of us cities[END_REF]. Existing heuristics however generally focus on a single morphological dimension of networks, and leave out the functional properties of urban systems [START_REF] Burger | Form follows function? linking morphological and functional polycentricity[END_REF].

This communication addresses such a gap by introducing a multi-dimensional percolation heuristic, which is analog to multilayer network percolation [START_REF] Boccaletti | The structure and dynamics of multilayer networks[END_REF]. Given discrete spatial fields, site percolation is operated between two cells given a threshold parameter for each dimension and a distance threshold. We apply the heuristic to urban morphology and road network topology measures in Europe. More precisely, a grid with resolution 50km of population density morphology indicators and road network topology indicators, has been computed on spatial moving windows for all European Union by [START_REF] Raimbault | An urban morphogenesis model capturing interactions between networks and territories[END_REF]. We percolate the population density layer with a network characteristic layer, that we test among number of edges, number of vertices, cyclomatic number and euclidian efficiency, which capture functional properties especially for the two last.

We systematically explore the clusters obtained for 840 parameter configurations. Maps reveal that most configurations resemble the actual distribution of European mega-city regions, which are functionally integrated polycentric urban areas [START_REF] Hall | The polycentric metropolis: learning from mega-city regions in Europe[END_REF]. We use this endogenous definition of regional urban systems the percolation algorithm produced to evaluate their sustainability, in terms of conflicting objectives of economic integration and greenhouse gases emissions.

Applying a gravity model to each region, we estimate transportation flows within each and extrapolate emissions by coupling with the Edgar emission database [START_REF] Janssens-Maenhout | Edgar v4. 3.2 global atlas of the three major greenhouse gas emissions for the period 1970-2012[END_REF] and economic activities with a scaling law of population. We show therein that different population, network and distance thresholds yield different performances in terms of sustainability, exhibiting a Pareto front. This suggests policies in terms of regional integration to increase the sustainability of mega-city regions. Further work will consist in the use of calibration heuristics [START_REF] Reuillon | Openmole, a workflow engine specifically tailored for the distributed exploration of simulation models[END_REF] to find in a more robust way optimal parameter values.

  3611 ISC-PIF et UMR CNRS 8504 Géographie-cités) Keywords : Road network; multilayer percolation; mega-city region