Struct Multidisc Optim DOI 10.1007/s00158-016-1468-4
INDUSTRIAL APPLICATION
Component allocation and supporting frame topology optimization using global search algorithm and morphing mesh Mahsan Bakhtiarinejad 1 & Soobum Lee 2 & James Joo 3
Received: 10 May 2015 / Revised: 14 April 2016 / Accepted: 21 April 2016 # Springer-Verlag Berlin Heidelberg 2016
Abstract This paper proposes a stepwise structural design methodology where the component layout and the supporting frame structure is sequentially found using global search algorithm and topology optimization. In the component layout design step, the genetic algorithm is used to handle system level multiobjective problem where the optimal locations of multiple components are searched. Based on the layout design searched, a new Topology Optimization method based on Morphing Mesh technique (TOMM) is applied to obtain the frame structure topology while adjusting the component locations simultaneously. TOMM is based on the SIMP method with morphable FE mesh, and component relocation and frame design is simultaneously done using two kinds of design variables: topology design variables and morphing design variables. Two examples are studied in this paper. First, TOMM method is applied to a simple cantilever beam problem to validate the proposed design methodology and justify inclusion of morphing design variables. Then the stepwise design methodology is applied to the commercial Boeing 757 aircraft wing design problem for the optimal placement of multiple components (subsystems) and the optimal
* Soobum Lee
[email protected] 1
Department of Mechanical Engineering, University of Maryland Baltimore County, 1000 Hilltop Circle, Baltimore, MD 21250, USA
2
Department of Mechanical Engineering, University of Maryland Baltimore County, 1000 Hilltop Circle, ENG214, Baltimore, MD 21250, USA
3
Design & Analysis Branch, Aerospace Vehicles Division, Air Force Research Laboratory, 2130 8th St., WPAFB, OH 45433, USA
supporting frame structure around them. Additional constraint on the weight balance is included and the corresponding design sensitivity is formulated. The benefit of using the global search algorithm (genetic algorithm) is discussed in terms of finding the global optimum and independency of initial design guess. It has been proved that the proposed stepwise method can provide innovative design insight for complex modern engineering systems with multi-component structures. Keywords Topology optimization . SIMP . Component relocation . FE morphing . Aircraft wing
Abbreviations ρ Material density p Penalization exponent c Compliance F Global force vector U Global displacement vector ue Element displacement vector K Global stiffness matrix Element stiffness matrix ke xt Topology design variables (elemental density) Morphing design variables (component xs displacement) V Material volume V0 Volume fraction CGx x directional center of gravity y directional center of gravity CGy CGLx Lower bound of CGx CGUx Upper bound of CGx CGLy Lower bound of CGy CGUy Upper bound of CGy Wk Weight of kth component Weight of frame structure Wfr
M. Bakhtiarinejad et al.
ak,x ak,y Ae ecxe ecye h Ma Mmax Wa Wf Ww CGa,y CGf,y CGt,y th Fa
x directional centriod point of kth component y directional centriod point of kth component Area of jth non-passive element x directional centroid of eth non-passive element y directional centroid of eth non-passive element Small pertubation of xs Aerodynamics moment Maximum threshold of moment Weight of the aircraft Weight of the fuselage Weight of the wing y directional center of gravity of the aircraft y directional center of gravity of the fuselage y directional center of gravity of the tail Thickness of the elements Aerodynamic force
1 Introduction As modern engineering systems require higher energy efficiency, structural optimization for compact and lightweight structure design becomes increasingly important. A great amount of work has been done in structural topology design community to find the structural layout for the best performance while satisfying various constraints such as limited volume/weight fraction. Over the last few decades researchers provided various applications of the density-based topology optimization to a variety of engineering disciplines. Later topology optimization method has evolved to design structural systems with multiple components; design optimization of both the component locations and the frame structure topology is conducted simultaneously. Since Bendsøe and Kikuchi (1988) introduced the optimization of material distribution through the use of homogenization approach, which can develop the optimal shape as well as the topology of a mechanical element using an artificial composite material with microscopic voids, topology optimization method has been developed as a widely accepted design tool to determine conceptual design of structures. SIMP, a density-based method, has been the most popular methodologies for structural topology optimization since it was developed in the late 1980s. The SIMP method was originally developed by Bendsøe (1989), Zhou and Rozvany (1991), and Rozvany et al. (1992). Bendsøe (1989) introduced a method to remove discrete design variables by introducing an artificial density. The shape of the mechanical component was defined by domains of high density and then homogenization method was utilized to compute the effective intermediate densities by defining periodically dispersed, microscopic void. Rozvany et al. (1992) showed that the SIMP method can efficiently develop optimal topologies for intermediate densities. Zhou and Rozvany (1991) proposed continuum-based optimality
criteria (COC) algorithm for the simultaneous optimization of the topology and geometry of grid-type structures. Topology optimization has been recognized as one of the most effective and influential design methods in aeronautics engineering with various applications, e.g., airframe structures, stiffener ribs for aircraft panels, multi-component layout design for aerospace structural systems, multi-fasteners design for assembled aircraft structures (Zhu et al. 2015). Employing topology optimization in aerospace structures, one can acquire a light-weight, stable, and stress-efficient aircraft component/ system designs (Deaton and Grandhi 2014; Krog et al. 2004). Applications of topology optimization methods to flight vehicle structures can be categorized into two groups: whether aero-elastic coupling is considered or not. The first group obtained classical objective functions such as compliance by aerodynamic load without aero-elastic coupling. For this purpose, pressure distribution over a wing needs to be computed. Balabanov and Haftka (1996) considered a grand structural approach to optimize the internal wing structure of high speed civil transport (HSCT) aircraft. The aircraft internal truss structure was modeled and the overall arrangement of spars and ribs was determined to minimize compliance. Eschenauer and Olhoff (2001) developed the conceptual layout of aircraft wing ribs through topology optimization using two load cases: (i) air loads along the rib length; and (ii) tank pressure. Krog et al. (2004) also proposed several approaches including a global optimization approach and two local optimization approaches to design wing box ribs for a variety of design criteria. A local optimization with free-body loads analysis identified the solution for both the global wing bending loads and local compression rib loads. An alternative local optimization approach used free-body loads analysis to demonstrate the energy in a baseline design (Krog et al. 2008). Wang et al. (2011) proposed a subset simulation based topology optimization method to establish the optimal design of the wing leading-edge ribs with the combination of the advanced sensitivity filtering technique and subset simulation. Choi et al. (2011) carried out Computational Fluid Dynamics (CFD) simulation accompanied by the equivalent static load method (ESL) to obtain both path and dynamic topology optimization of flapping wings. Oktay et al. (2011) also proposed a design approach for aircraft wing structure with density-based topology optimization and CFD analysis. Another group of aircraft applications considered aerodynamic forces on a flexible wing which inherently includes aeroelastic coupling between the aerodynamic and structural analysis. Maute and Allen (2004) developed the geometrical layout of the wing stiffener by the SIMP method and sequential augmented Lagrangian method to solve the large-scale optimization problem. This study was expanded later for the design of mechanisms inside an airfoil in morphing aircraft structures. Design of a quasi three-dimensional section of an adaptive wing was conducted with topology optimization of
Component allocation and supporting frame topology optimization
mechanisms layout and design optimization of coupled fluidstructure problems (Maute and Reich 2006). Stanford and Ifju (2009a) conducted topology optimization based on the SIMP approach to design a membrane (wing skin) and carbon fiber (laminate reinforcement) of micro air vehicle wing with a twomaterial formulation. Topology optimization was performed for variety of objectives including load-augmenting, load-alleviating, pitching moment and efficiency. Optimal loadaugmenting and load-alleviating were obtained to connect stiffness distribution and aerodynamic performance of membrane micro air vehicle wings (Stanford and Ifju 2009b). Stanford and Beran (2010) utilized nonlinear finite element model to obtain a thrust-optimal compliant drive mechanism design for flapping wings via density-based topology optimization. More recently, aeroelastic tailoring was developed as a fundamental procedure in composite material wing design. With maximization of natural frequencies related to the vibration modes, the flutter speed will be increased which is the goal of the structural optimization process for composite laminated wing. De Leon et al. (2012) utilized topology optimization and composite fiber orientation to determine aeroelastic tailoring for design of thin wings with minimum mass under a minimum flutter velocity constraint. Topology optimization has evolved to design structural systems with multiple components. That is, design optimization of the component locations and the frame structure topology is conducted simultaneously. Qian and Ananthasuresh (2004) presented a new embedding technique in topology optimization which maximizes the structural stiffness by connecting the structure to one or more components. The simultaneous optimization of position (and orientation) of components and structure’s topology was conducted through a proposed approach. Location and orientation of embedded objects were assigned as a new set of design variables in 2D design space. Zhu et al. (2008) proposed a coupled shape and topology optimization (CSTO) technique for simultaneous design of components layout and their supporting structures. Topology design optimization of supporting structure and interconnect components are carried out simultaneously. By assigning density variables to predefined density points rather than finite element mesh and embedded meshing techniques, the finite element mesh is allowed to
Fig. 1 Two step approach for component layout design and frame topology design
update with the movement of components (Zhang et al. 1995). Later the position of the system center gravity, mostly in aerospace and aeronautical systems, is considered as an additional constraint and the corresponding effects are analyzed (Zhu et al. 2009). Inertial forces (self-weight) are also introduced into the design optimization for multiple component layout and a special material interpolation model between element stiffness and inertial loads is applied to avoid the localized deformation in the low density region (Zhu et al. 2010). Zhang and Zhang (2009) proposed the finite circle method (FCM) for 2D packaging optimization where a group of circles approximate the exact shape of each component and the nonoverlapping constraints can be transformed into explicit distance constraints between circles. Both genetic and gradient-based algorithms are incorporated in the FCM to solve the problem. Although it has been demonstrated that FCM can detect the overlap between the components and improve the optimization convergence, the large number of constraints affected the optimization efficiency adversely and because of the inherent nature of the NP-hard problem, it is difficult, in certain cases, to obtain a global optimum solution. (Zhang et al. 2011). Recently, Xia et al. (2012a) suggested a new sensitivity analysis approach with respect to location design variables; a material perturbation model (MPM) using Heaviside function is developed with the concept of fixed finite element mesh. In other words, when the locations of components change, only material properties are modified while the finite element mesh itself does not change. Therefore, the sensitivities of components can easily be computed and the computing efficiency is improved. Additionally, Xia et al. (2012b) proposed a new formulation named superelement formulation (SEF) to improve computing time and efficiency of simultaneous design optimization of component layout and frame structure topology. More recently, Zhu et al. (2014) proposed a new multipoint constraints (MPC) based layout and topology optimization method. The components were easily able to float and bound on the surface of the topology design domain by introducing the rivets and bolts connections between the components and their supporting structures. Avoidance of finite elements remeshing and precise geometry of components were obtained with MPC method. A new computational framework
M. Bakhtiarinejad et al.
Fig. 2 Interference check between component i and component j
for topology optimization of structures with moving morphable components was suggested by Guo et al. (2014). Rather than pixel or node point-based topology optimization, the proposed framework applied topology optimization in an explicit and geometrical way, which significantly reduced the computational difficulties accompanied by topology optimization algorithms. Later Gao et al. (2015) introduced a Kreisselmeier-Steinhauser (KS) function based adaptive constraint aggregation approach to deal with the large number of non-overlapping constraints of the integrated layout and topology optimization of multi-component structures. MPC with fixed finite element meshes and analytical sensitivities was also applied to provide the interconnection between moveable components and supporting structures. Zhang et al. (2015) proposed new distance control methods for optimizing the layout of structural systems with embedding components. For this purpose, the minimum and maximum
Fig. 3 Component shift by morphing mesh: a Cantilever beam with void rectangular component; b morphed component in beam domain
distance constraints between the components were considered, which can give a complete control of the layout of embedding component in an explicit way. Although researchers accomplished a variety of designs for topology optimization of multi-component structures, most of the obtained design results were dependent on an approximated FE models such as fixed finite element mesh (Gao et al. 2015; Xia et al. 2012a; Zhu et al. 2014) or density point (Xia et al. 2012b; Zhu et al. 2008, 2009, 2010) and they may cause inaccurate sensitivity. This is a critical issue when a sensitivitybased design algorithm is used. The density point approach reduced the computational cost on remesh, but it still requires local remesh every iteration whenever the component locations change. Additionally, a sensitivity-based approach cannot guarantee global optimality of the solution, especially for component locations which are truly dependent on their initial locations. In order to address the issues described above, this paper proposes a new stepwise structural design methodology where the component layout and the supporting frame structure is sequentially found by global search algorithm and topology optimization. In the component layout design step, the genetic algorithm is used to determine optimal locations of multiple components. Based on the layout design results, a new topology optimization methodology entitled as Topology Optimization based on Morphing Mesh technique (TOMM) is followed. This methodology is basically a density-based topology optimization methodology, but it utilizes a morphing FE mesh for component relocation and for more accurate design sensitivity calculation. The paper is organized as follow. In Section 2, the of stepwise design methodology for component location and frame structure topology is explained, including the design formulation and sensitivity analysis. In Section 3, multiple design cases are studied using the proposed design method—the cantilever beam and the Boeing 757 wing— and the benefit of this approach is discussed in depth. Section 4 concludes the proposed design methodology and its contributions.
Component allocation and supporting frame topology optimization Fig. 4 Flow chart of simultaneous design optimization of component relocation and frame structure topology with Matlab and HyperMesh
Start
1. Define Objective Function Compliance Define Constraint Functions Weight Balance Volume Fraction
Initial Morphing Design Variables
2. Assign Initial Design Variables
Mesh Generation Using HyperMesh Initial Topology Design Variables
HyperMesh Output Files 3. Mesh Generation Update Using Morphing Topology Design Variables
4. Do FE Analysis
Morphing Design Variables
5. Calculate Objective Function and Constraints 7. Update Design Variables
6. Stopping Criteria met?
NO
YES End
2 Design methodology
energy efficiency. The proposed design methodology provides a systematic approach by performing the consecutive two design steps: (i) component allocation, and (ii) structural topology optimization. The scheme of this methodology is shown in Fig. 1. In the first step, the design problem to allocate each component is solved considering system level multiobjective functions including dynamic stability and operation costs. Non-intuitive design results can be
Modern engineering systems are very complex with its various components. Such components in an aircraft system include: control system (flight control, engine control), fuel system, hydraulic system, and structural system. An advanced systematic design approach is needed for effective packaging of multiple components while satisfying load-carrying performance and minimizing weight for
Fig. 5 Design domain of tiploaded cantilever beam with void component and the corresponding FE model H = 20
y
F
x
L = 12
M. Bakhtiarinejad et al. Fig. 6 Topology optimization results for the cantilever beam problem with different locations of void component: a void location 24; b void location 28; c void location 32; d void location 36; e void location 40; f void location 44
(a) c = 40.5513
(b) c = 39.7212
(c) c = 38.4037
(d) c = 40.8060
(e)
(f)
expected by adapting global optimization algorithms such as the genetic algorithm. Component locations are found, not overlapping with each other and with the design boundary. Handling multiobjectives, multiple solutions can be obtained in this step. In the second step, topology optimization for structural robustness is performed based on the solution found in the first step. The design optimization for the determination of frame structures is formulated as compliance minimization problem subject to the constraints on volume fraction and the system level performance found in the first step (e.g., weight balance). The component layout found in the first step provide loading and boundary conditions and no frame structure is placed in them. The detail formulation for each design step is explained in the following subsections. 2.1 Component allocation The design formulation in the component allocation step will determine the location of nc number of components or xs = {x1, y1, … xnc, ync} (2D case) to satisfy system level performance such as the center of gravity (CG) requirement, proximity (inter-distances between components) and interference prevention among the components. The CG requirement indicates that the center of gravity of the whole structure agrees with its target, and it is very important for an aerospace
application. Examples of proximity may include: component 1 (aircraft landing gear) needs to be located close enough to fuselage for landing stability; component 2 (engine) and component 3 (control system) need to be apart from each other in order to prevent overheat. This allocation problem is formulated as a multiobjective design optimization as follows: mini¼1;…D jCGðxs Þ−T ji s:t: gLl ðxs Þ ≤ 0;
gU m ðxs Þ≤ 0;
ð1Þ
gI ðxs Þ≤ 0
where CG(xs) is the center of gravity, T is its target, and i is the dimension index. Depending on the dimension of the problem and the design requirement, D (number of objective functions) can be 1, 2 or 3 (1D, 2D or 3D problem, respectively). l, m are the indexes for the proximity constraints (g Ll ; gU m ) with lower bound and upper bound, respectively. gI is the interference constraint. For applying multiobjective GA algorithm available in the Matlab (gamultiobj), this formulation is rewritten as an unconstrained multiobjective minimization as: n mini;k CGðxs Þ−T i ;
gLl ðxs Þ;
gU m ðxs Þ;
o g I ð xs Þ :
ð2Þ
Among the multiple solutions obtained from this formulation, only the feasible solutions are taken (gLl ðxs Þ ≤ 0; gU m ð xs Þ ≤ 0; gI ðxs Þ ≤ 0 ).
Component allocation and supporting frame topology optimization
(a)
(a)
(b)
(b)
Fig. 7 Topology optimal results of the cantilever beam: a topology without xs (Fig. 6d); b topology with xs
The proximity constraints can be easily handled by the distance function d(xs) such as 0; if −d Ll ðxs Þ þ DLl ≤ 0 g Ll ¼ 1; otherwise ð3Þ U 0; if d U ð x Þ−D ≤ 0 U s m m gm ¼ 1; otherwise for the cases where the lower (DLl ) or upper bound (DU m) of distance is assigned, respectively. They are defined as a discrete function between 0 and 1. The interference constraint is also defined as a discrete function as 0; if no overlap ð4Þ gI ¼ 1; otherwise and it is figured out by checking if there is any intersecting point between any two components or between any component and design boundary. Figure 2 shows an example to check the Table 1
Compliance for topology results of the cantilever beam
Solution
(a) Without xs
(b) With xs
Compliance
40.8060
39.6898
Fig. 8 Iteration history of the cantilever beam: a without xs; b with xs
existence of any intersection point. Components i (ci) and j (cj) are defined as a polygon composed of four (Li,1, … Li,4) and five lines (Lj,1, … Lj,5), respectively. We implemented a double loop code to check if any line belonging to ci intersects any line belonging to cj (Lee and Kwak 2008; Lee et al. 2007) This code also checks interference between components and the design boundary as well. If a point is found, the code returns 1 as indicated in Eq. (4). 2.2 Topology optimization with morphing mesh (TOMM) Based on the solutions found in the component allocation design step, topology optimization is applied to determine optimal frame structure. This section presents a new design methodology, Topology Optimization based on Morphing Mesh technique (TOMM) for simultaneous design optimization of component relocation and frame structure topology based on the general SIMP approach. The variable for component placement (xs) in the first design step is continuously used for refined relocation of the components using morphing mesh. The sensitivity analysis of objective function (compliance) and constraints (volume fraction, weight balance) with respect to both topology and morphing design variables is carried out in Section 2.2.2. The formulation of the compliance minimization of the SIMP density-based topology optimization (Sigmund 2001)
M. Bakhtiarinejad et al. Fig. 9 Von-Mises stress contour for void location 32: a overall von-Mises stress distribution; b local xt distribution at one corner
(a) can be updated for simultaneous design optimization of components layout and structural topology optimization problem. Additional constraints on weight balance needs to be considered in order to maintain the center of gravity found in the first step. When a 2D design space is considered, the design problem can be formulated as
Find
x ¼ fxt ; xs g ¼ xt;1 ; xt;2 ; …xt;N t ; xs;1 ; xs;2 ; …xs;N s to
min
cðxt ; xs Þ ¼ FT Uðxt ; xs Þ
subject to
g1 ¼ CGLx −CGx ðxt ; xs Þ≤ 0 g2 ¼ CGx ðxt ; xs Þ−CGU x ≤ 0 g3 ¼ CGLy −CGy ðxt ; xs Þ≤ 0 g4 ¼ CGy ðxt ; xs Þ−CGU y ≤ 0 g5 ¼ V ðxt ; xs Þ−V 0 ≤0
ð5Þ where x = {xt, xs} is the set of design variables, combination of: topology design variables (elemental density variables) xt ¼ xt;1 ; xt;2 ; …xt;N t ; and morphing design variables xs ¼ xs;1 ; xs;2 ; …xs;N s ¼ fx1 ; y1 ; …xnc ; ync g which are the relative displacements of components as used in the component allocation design step. Nt is the number of topology design variables and Ns is the number of morphing design variables (Ns =2×nc in 2D case). U and F are the global displacement and force vectors, K is the global stiffness matrix. U is the
Fig. 10 Boeing 757 aircraft (Wikipedia 2014)
(b)
function of both xt and xs such that U(xt, xs) = K(xs, xt)− 1F. CGLx (CGLy) is the lower bound of CGx (CGy), and CGUx (CGUy) is the upper bound of CGx (CGy), respectively. The constraints for the weight balance, CGLx ≤ CGx ≤ CGUx and CGLy ≤ CGy ≤ CGUy, are separated into four inequalities for x (g1 and g2) and y (g3 and g4) directions. Note that two additional constraints can be added when 3D design space is considered. g5 is the volume constraint formulated with V(xt, xs) and V0, which are the material volume and its fraction, respectively. The design domain is composed of passive zone and non-passive zone. The components are defined as passive zone (or non-design zone) where the topology design variables are not assigned. Frame structures are designed in non-passive zone (or design zone) where the topology design variables xt are assigned at each finite element. The optimal topology for frame structure is found while each component shifts in random direction by updating morphing design variables. The bounds for xs are assigned such that the Jacobian value, a measure of finite element distortion (Cook 2007), is always higher than the threshold (0.7 chosen in this paper). A Matlab program to control FE Morphing function in HyperMesh (HyperWorks 2014); a pre-processor software for finite element modeling, is implemented in order to realize finite element morphing. Figure 3 shows the use of morphing mesh in a simple cantilever beam with a rectangle void component. In Fig. 3b, the
Component allocation and supporting frame topology optimization Table 2 Descripting on subsystems of Boeing wing
Subsystem index (j)
1
2
3
4
–
Subsystems Weight (kg)
Landing gear 2312
Engine mount 3381
Fuel tank 1 5270
Fuel tank 2 2455
Wing framea 4382
a
To be designed in the second step using TOMM
FE mesh for the beam domain (dark grey) adaptively morphs (distorts) when the component (white) shifts randomly in x direction. Using this mesh one can (i) exactly maintain the boundary of each component and (ii) realize simultaneous optimization on component relocation and topology design.
mx(e) = ecxe × xt,e, my(e) = ecye × xt,e. The material volume V(xt, xs) in g5 (Eq. (5)) is formulated as X xt;e Ae ðxt ; xs Þ V ðxt ; xs Þ ¼ Xe : ð7Þ A ð x ; x Þ e t s e
2.2.1 Constraints formulation The optimization problem in Eq. (5) includes five constraints: four constraints for the center of gravity (weight balance) and one constraint for volume fraction. The center of gravity (CGx, CGy) is a function of xt and xs and they are formulated as: X
X
Ae mxðeÞ W k xs;2k−1 þ ak;x þ W f r : Xe Ax e e t;e CGx ¼ ; X W þ W k f r k ð6Þ X X A my ð e Þ e W x þ ak;y þ W f r : Xe k k s;2k Ax e e t;e C Gy ¼ X W þ W fr k k
k
where Wk and Wfr are the weight of kth component and frame structure, respectively. (ak,x, ak,y) is the centroid point of kth component, and Ae is the area of eth non-passive element. (ecxe, ecye) is the centroid of eth non-passive element and
y
2
4
8.52
1
3
16.8
Fig. 11 Wing design domain and four subsystems (in m)
Design sensitivity analysis is required in a sensitivitybased optimization algorithm. In this section, the sensitivity of objective function and constraints with respect to both topology and morphing design variables are explained. The sensitivity of objective function with respect to topology design variables is written as (Bendsoe and Sigmund 2003; Sigmund 2001) p−1 T ∂c ¼ −p xt;e ue ke ue ∂xt;e
ð8Þ
where ue and ke are the element displacement vector and stiffness matrix, p is the penalization power which is considered to be 3 in this paper. The sensitivity of objective function with respect to morphing design variables can be found numerically by finite difference method as c xs;1 ; …; xs;i þ h; …; xs;N s −c0 ∂c ¼ h ∂xs;i
ð9Þ
where c xs;1 ; …; xs;i þ h; …; xs;N s is the compliance of design variable with small perturbation of xs,i by h, and c0 is the compliance of current design variables. In this paper, h is chosen to be a small number (0.25–1 % of the morphing design variables).
2.33
x
2.2.2 Sensitivity analysis
Table 3 Lower and upper bounds of center of gravity
CGLx
CGUx
CGLy
CGUy
5
7
−6.943
−4.286
M. Bakhtiarinejad et al. Fig. 12 Subsystem layout design results: a Design A; b Design B
(a)
(b)
The sensitivities of weight balance and volume fraction constraints can be written in matrix form, with respect to topology and morphing design variables respectively, such that
The sensitivity of volume fraction constraint with respect to topology design variables (the last four columns in Eq. (10)) is found as
2 ∂g 1 6 ∂x t;1 ∂g 6 ¼6 ⋮ ∂xt 4 ∂g1 ∂xt;N t
Aj ∂g 5 ¼X : ∂xt;e Ae
2 ∂g 1 6 ∂x s;1 ∂g 6 ¼6 ⋮ ∂xs 4 ∂g 1 ∂xs;2N s
⋯ ⋱ ⋯
⋯ ⋱ ⋯
∂g5 ∂xt;1 ⋮ ∂g5 ∂xt;N t
3 7 7 7; 5
∂g5 ∂xs;1 ⋮ ∂g5 ∂xs;N s
ð10Þ
The sensitivities of weight balance constraints with respect to morphing design variables (the first four columns in Eq. (11)) is obtained as
3 7 7 7: 5
ð11Þ
In each matrix, the first four columns indicate the weight balance sensitivity while the fifth column is the volume fraction sensitivity. The sensitivities of weight balance constraints with respect to topology design variables (the first four columns in Eq. (10)) is obtained as h i X X W :A ecx A x − A mx ð e Þ f r e e e t;e e ∂g 1 ∂CGx e e ¼− ¼ − X ; X 2 ∂xt;e ∂xt;e : W þ W A x fr k k e e t;e ∂g 2 ∂CGx ¼ ; ∂xt;e ∂xt;e h i X X W f r :Ae ecye A x − A my ð e Þ e t;e e ∂CGy ∂g 3 e e ¼− ¼ − X ; X 2 ∂xt;e ∂xt;e : W þ W A x k f r e t;e k e ∂CGy ∂g 4 ¼ : ∂xt;e ∂xt;e ð12Þ Table 4
Subsystem layout results for Design A and B
Design
Landing gear Engine mount Fuel tank 1
Design A (0.98, −4.81) (2.62, −1.78) Design B (3.02, −4.82) (0.85, −1.00)
ð13Þ
Fuel tank 2
(4.54, −3.13) (9.22, −5.63) (4.94, −3.32) (10.15, −6.05)
x and y coordinate of upper/left corner as indicated in Fig. 12
∂g 1 ¼ ∂xs;i 2
8 < X −W k ; when i ¼ 2k−1 W þ W k f r : k 0 ; otherwise
¼ 4X
−W k
−W k
3T
0 X 0; … 5 ; W þ W W þ W k fr fr k k 8 k Wk < ; when i ¼ 2k−1 X ∂g 2 ¼ W þ W k f r : k ∂xs;i 0 ; otherwise 2 3T W W k k ¼ 4X 0X 0; … 5 ; W þ W W þ W k f r k f r k 8k −W k < ; when i ¼ 2k X ∂g 3 ¼ W þ W fr : k k ∂xs;i 0 ; otherwise 2 3T −W −W k k ¼ 40 X 0 X ; …5 ; W þ W W þ W k fr fr k k 8 k W k < ; when i ¼ 2k X ∂g 4 ¼ W þ W k f r : k ∂xs;i 0 ; otherwise 2 3T W W k k ¼ 40 X 0X ; …5 : W þ W W þ W k others k f r k k ð14Þ
Component allocation and supporting frame topology optimization
Fig. 13 Finite element model of Boeing wing: a design A; b design B
(a) The sensitivity of volume fraction constraint with respect to morphing design variables (the last four columns in Eq. (11)) is found by finite difference method as g xs;1 ; …; xs;i þ h; …; xs;N s −g 5 xs;1 ; …; xs;i ; …; xs;N s ∂g 5 ¼ 5 ∂xs;i h
ð15Þ 2.2.3 Code implementation In order to solve the optimization problem in Eq. (5), a sensitivity-based algorithm (interior point algorithm) is utilized which is robust and efficient for solving nonlinear convex optimization problems (Byrd et al. 1999). This algorithm follows a barrier method that applies sequential quadratic programming and trust regions to solve the optimization subproblems. The Hessian matrix is approximated using LBFGS method, a limited-memory, large-scale quasi-Newton approximation particularly well suited for optimization problems with a large number of variables (Mathworks 2015). In this paper, the Hessian is calculated using 25 iteration memory (TopOpt Research Group 2015). The detail procedure in the interface between Matlab and HyperMesh for simultaneous design optimization of component relocation and frame structure topology is shown in Fig. 4 After the design optimization problem is formulated (Step 1), the initial guess of the design variables are assigned (Step 2). The components are relocated according to the updated morphing design variables, and the FE mesh is updated using morphing
Fig. 14 Loading and boundary conditions of Boeing wing
(b)
technique (Step 3). After constructing the element stiffness matrix, FE analysis is performed to obtain the global displacement vector (Step 4). The objective function, compliance, and the constraints are calculated (Step 5). The optimality condition and the stopping criteria of optimization algorithm are checked in every iteration (Step 6). If these conditions are not met, the design variables and optimization problem will be updated (Step 7). The optimization procedure is done iteratively until the specified stopping criteria are met.
3 Case studies This section explains two different case studies for the simultaneous design optimization of component location and structural topology and discusses the benefit of the proposed design strategy.
3.1 Cantilever beam The first case study is the cantilever beam with one component (a rectangular hole). This case study is to verify the TOMM method (topology optimization with morphing variable) and to justify inclusion of morphing design variables. A fixed shape of rectangular void is considered as a component, which can be relocated through the design optimization. Figure 5 shows the shape of a tip-loaded cantilever beam and its design domain (grey). The design domain size is 60 × 40 and the void component (a rectangular hole) has a 12 × 20 dimension that is considered as a passive zone. The nonpassive zone is discretized into 2160 finite elements where the topology design variables (xt) are assigned. One design variables xs is additionally assigned for the location of void component along x direction. The component is bounded to move between -5 and 5. The beam is fixed along the left edge and a point load Fy = − 1 is applied in the middle of the free end. 4-node quadrilateral elements with 2 DOF (u, v) at each node are used for finite element modeling.
M. Bakhtiarinejad et al. Fig. 15 Topology optimal results of the Boeing 757 wing
Design B
with
(Case 2)
without
(Case 1)
Design A
The weight balance constraints are not considered in this study. The optimization problem in Eq. (5) is rewritten for this cantilever beam case with 50 % volume fraction as follows: Find x ¼ fxt ; xs g ¼ xt;1 ; xt;2 ; …xt;2160 ; xs;1 to Min
c ðxt ; xs Þ ¼ FT U ðxt ; xs Þ
Subject to
g 1 ¼ V ðxt ; xs Þ − V 0 ¼ V ðxt ; xs Þ − 0:5 ≤ 0 − 5 ≤ xs ≤ 5
ð16Þ The termination tolerance is set as 1×10−6 for the objective function and 1×10−10 for constraints and design variation. The filtering radius is considered to be 1.2. 3.1.1 Results The benefit of including morphing design variable (xs) is explained in this section. If one needs to find the optimum void location without morphing design variable, multiple topology optimization studies with different void locations are required as shown Fig. 6. This figure displays the topology optimization results of cantilever beam with six different locations of void component and the corresponding compliance values. The distances between the void component and the left edge are set as 24, 28, 32, 36, 40 and 44, respectively. Table 5
Compliance for topology results of the Boeing wing
Compliance (J)
Design A
Design B
without xs (Case 1) with xs (Case 2)
9.3063 × 105 9.2267 × 105
1.2787 × 106 1.2363 × 106
As indicated in Fig. 6, the void location of 32 (Fig. 6c) has a minimum of compliance value (c = 38.4037). From this study it is concluded that the optimum void location is around 32, but additional topology optimization studies are required for the refined search on true optimum location. The TOMM method is attempted to find optimum void location and frame structure simultaneously as shown in Fig. 7. Initial guess of the void location is chosen to be 36, the same as the forth fourth solution in Fig. 6d (compliance = 40.8060 when no xs variable is considered). On the other hand, when topology optimization design is performed with xs, the compliance is improved to 39.6898 (Fig. 7b). Both morphing and topology design variables are shown in Fig. 7b. The compliance variables are also displayed in Table 1. The iteration histories of two solutions are shown in Fig. 8; it took 35 and 34 iterations for two solutions to achieve the convergence respectively. Two graphs became different after 17th iteration from where the second solution is converging to smaller compliance. One can understand that both solutions requires almost the same number of iterations but the second solution with xs is converged with smaller value of compliance.
3.1.2 Discussion The results clearly indicate that simultaneous design optimization of component layout and structural topology is more effective than conducting multiple topology design optimizations separately. In other words, rather than employing multiple topology optimization trials with different locations of components, one can carry out single topology optimization while determining the optimal location of components at the same time.
Component allocation and supporting frame topology optimization Table 6 Morphing design variables after TOMM design step
Design
Landing gear
Engine mount
Fuel tank 1
Fuel tank 2
Design A
1.03, −5.01a (+0.05, −0.20)b
2.70, −1.90 (+0.08, −0.12)
4.44, −3.31
9.12, −5.73
(−0.1, −0.18)
(–0.1, −0.1)
4.89, −3.52
10.05, −5.87
(–0.05, −0.20)
(–0.10, +0.18)
Design B
2.92, −4.89 (–0.1, −0.07)
0.75, −0.87(–0.10, +0.13)
a
x and y coordinate of upper/left corner as indicated in Fig. 15
b
x and y coordinate of morphing design variables
As shown in Fig. 7b, using the TOMM method, the void component is shifted to the left (xs = –3.99) and the FE mesh of the frame structure (design domain) is adaptively modified. The result by TOMM method shows more improved compliance 39.69 than the case without the morphing variable (Fig. 7a, compliance = 40.81). The void location found by the TOMM method is 32.01 (=36.0 – 3.99) and the closest void location in Fig. 6 is Fig. 6c (32) where the compliance value is the smallest. This observation verifies that TOMM approach can find the optimal location of void component and structural topology without conducting multiple topology optimization trials. It is noted that the compliance value by TOMM is dependent on the FE mesh quality. The distorted mesh such as Fig. 7b may result in increased compliance value when compared to the regular fixed FE mesh with the same void location. The topology design result shows sharp corners around the four vertexes of rectangular void and one may claim that a significant stress concentration is expected in the corners. However, as shown in Fig. 9a, the stress is evenly distributed by the elements with xt > 0.990. Actually the structural boundary around the corner looks sharp because of intermediate xt values (numbers with grey color in Fig. 9b, 0.1< xt