Sub-space approximations for MDO problems with disparate disciplinary variable dependence

An approach to solving multidisciplinary design optimisation problems using approximations built in sub-spaces of the design variable space is proposed. Each approximation is built in the sub-space significant to the corresponding discipline while the optimisation problem is solved in the full design variable space. Since the approximations are built in a space of reduced dimensionality, the computational budget associated with building them can be reduced without compromising their quality. The method requires the designer to make assumptions on which design variables are significant to each discipline. If such assumptions are deficient, the resulting approximations suffer from errors that are not possible to reduce by additional sampling. Therefore a recovery mechanism is proposed that updates the values of the insignificant variables at the end of each iteration to align with the current best point. The method is implemented within a trust region based optimisation framework and demonstrated on a multidisciplinary optimisation of a thin-walled beam section subject to stiffness and impact requirements.


Introduction
In multidisciplinary design optimisation (MDO) problems the individual disciplines are often evaluated using different computational models, making it not only MDO but also multi-model optimisation. At the same pace as computational performance is increasing, more complexity is built into the computational models, making them more computationally expensive. Thus when dealing with MDO problems one may face a very computationally expensive problem. For reviews on MDO architectures, see for instance (Cramer et al. 1993), Sobieszczanski-Sobieski and Haftka (1997), and most recently (Martins and Lambe 2013) This work is focused on solving MDO problems using the multidisciplinary design feasible (MDF) architecture as formulated by Cramer et al. (1993) assisted by approximations. The reader is referred to Barthelemy and Haftka (1993), Wang and Shan (2007), Forrester and Keane (2009), and Viana et al. (2014) and Viana et al. (2010) for reviews on metamodel-based optimisation. The main challenge of using approximations is ensuring their quality. The number of training points required to obtain approximations of desired quality increases superlinearly with the number of design variables unless the responses are linear in the design variable space. For MDO problems the total number of design variables can become very large leading to the need for a large number of computationally expensive simulations.
Given that there are variables in the MDO problem that are not shared between all disciplines, i.e. have their own relevant subsets of design variables, approximations can be built on only the sub-set of the design variables that is significant to each discipline whilst the optimisation problem is solved in the full design variable space. There are several examples of the use of this approach in the automotive industry, e.g. by Sobieszczanski-Sobieski et al. (2001), Kodiyalam et al. (2004), Ollar et al. (2014) and Ollar et al. (2015) and Ryberg et al. (2015). The benefit is that sampling and approximation building can be carried out in a space of reduced dimensionality. This can greatly reduced the required computational budget for obtaining an approximation of sufficient quality. However, deficient assumptions, made when identifying significant variables, lead to approximation errors that can not be reduced by additional sampling.
In this work, an adaptive discipline-related sampling and approximation building technique is introduced within the trust region framework known as the Mid-Range Approximation Method (MAM). In contrast with the previous work, it is part of a trust region strategy that accounts for erroneous assumptions by means of a recovery strategy that updates the values of the variables classed as insignificant but remain present in the full variable space during the optimisation stage of the trust region strategy.

Mid-range approximation method
The mid-range approximation method (MAM), also known as the multi-point approximation method, was originally reported by Toropov (1989Toropov ( , 1992 and Toropov et al. (1993). The MAM solves a typical constrained optimisation problem in the form where f 0 (x) is the objective function, f j (x) is the j -th constraint, x is the vector of design variables and A and B are the upper and lower bounds respectively on the design variables. The optimisation problem (1) is replaced by a sequence of approximate sub-problems defined as where k denotes the current iteration number,f 0 k (x) is an approximation of the objective function andf j k (x) is an approximation of the j -th constraint function, both considered valid only in the current trust region. A k i and B k i are the bounds of the current trust region where the sub-problem (2) is solved for the current iteration. The solution procedure for each sub-problem consists of sampling, creating approximations, solving the approximate optimisation problem and determining a new location and size of the trust region for the next iteration. The trust region will move and change size after each iteration until the termination criterion is reached. Figure 1 illustrates the history of trust regions through the sequence of sub problems in two dimensions. The trust region strategy has gone through several developments to account for the presence of numerical noise in the response function values Toropov et al. 1996) and occasional simulation failures (Toropov et al. 1999).
In the beginning of each iteration sampling, according to a small-scale design of experiments (DOE), is carried out within the trust region. Available sampling techniques include a technique based on extensible lattice sequences proposed by Hickernell et al. (2000) and a non-collapsible randomised technique by Korolev et al. (2015). The sampling process benefits from considering existing points, evaluated in previous iterations, however, points that are far away from the trust region may spoil the approximation and are not considered. The region in which points are considered is referred to as the recycle region and is defined as an enlargement of the trust region according to Fig. 2.
Approximations are created for each response using the evaluated sampling points. Available approximation  techniques include approximation assemblies by Polynkin and Toropov (2012), where intrinsically linear and rational functions are assembled into a single approximation using linear regression, the moving least squares method (Lancaster and Salkauskas 1981) and kriging based on a computationally efficient hyper-parameter optimisation method by Mortished et al. (2016).
The approximate sub-problem (2) is solved using the sequential quadratic programming (SQP) method developed by Madsen et al. (2002), based on the work of Powell (1985), starting from several points in order to increase the chance of finding a good solution for multi-modal problems. Candidate solutions of the approximated sub-problem (2) are evaluated in parallel with additional sampling points as introduced by Korolev et al. (2015) in order to fully utilise parallel high performance computing hardware. The best candidate point is compared to the previously obtained best point and, if superior, stored as the current best point.
Finally, termination criteria are checked. The current best point may not be located at the bound of the trust region, the approximations must be of sufficient quality, and the trust region must be sufficiently small compared to the initial design region. If the criteria are met, the process is terminated. Otherwise the position of the trust region is moved so that its centre point coincides with the current best point. The trust region strategy changes the size of the trust region depending on the same factors as the termination criteria and the process is repeated.

Sub-space approximations
Suppose that the responses belonging to a discipline in the MDO problem only depend on a subset of the full set of design variables, i.e. a set of the variables has none or very little influence on the responses of the particular discipline. An example of this from the automotive industry can be seen in Fig. 3. The figure shows a simulation of a vehicle subjected to a front crash load case. Each element is coloured according to its level of internal energy. It can be concluded that, as can be expected, the energy absorption is concentrated in the front of the vehicle and is not much affected by the rest of the structure. A conceptual partitioning based on design variable dependence for four common automotive disciplines, Front Crash, Side Crash,  Fig. 3c. With this partitioning discipline-related approximations, hereafter denoted sub-space approximations, can be built, depending only on the set of significant variables for each discipline.

Formulation of sub-space approximations
Here a mathematical formulation for introducing sub-space approximations in approximation assisted MDO is given. Unlike in previous work, all response functions are assumed to be defined in the full variable space of the optimisation problem in order to control insignificant variables. It will be shown that when sub-space approximations are used within a trust region framework, this becomes necessary to account for possible deficient assumptions on partitioning.
Consider solving the optimisation problem (1) using approximations. The optimisation problem becomes wheref 0 (x) is an approximation of the objective function andf j (x) is an approximation of the j -th constraint function.
Given that the design variables in the optimisation problem can be categorised either as significant or insignificant for each related response. A projection can then be defined, for each response j , from the design variable space onto the space of the significant variables for that particular response. This is denoted as where n is the number of design variables in the optimisation problem and s j is the number of significant variables for the response j . A projection onto the space of the insignificant variables is defined in the same manner as From here on the projections are described according to the following convention noting that the components of ξ j and ψ j are present in x in arbitrary order. The responses in the optimisation problem can then be described as where the values of ψ j can be chosen arbitrarily since they are deemed to be insignificant to the response. The approximations of the responses, therefore, may now be defined in the space of only the significant variables which allows a re-writing of the approximate optimisation problem according to where each approximated response is defined only in the space of variables that are significant to the response. The optimisation problem, however, is defined in the full design variable space. This has the benefit that as each approximation is defined only in the space of the significant variables, the sampling of training points only needs to be carried out in that space while the values of the insignificant variables are kept constant as demonstrated by Fig. 4. Hence if the number of significant variables is small compared to the number of design variables, the density of the training points will increase leading to a better quality approximation as compared to what would have been achieved otherwise. Note that even though there is one projection per response, practicalities may require groups of responses to use the same projection, e.g. due to several responses being evaluated from the same simulation, or even responses belonging to different responses but connected through multi-physics coupling such as fluid structure interaction.

Sub-space approximations in trust regions
In this section an approach to building sub-space approximations within the MAM is proposed. A recovery mechanism for erroneous identification of significant variable for the sub-space partitioning is also suggested. Sub-space approximations can be introduced in the MAM framework by re-writing the sequence of optimisation problems (2) as: noting that the mid-range approximations created in each iteration are now functions of the significant variables only. The significant variables for each discipline are identified by the designer. Such judgement may be based on, for instance, engineering experience or design variable ranking studies. In the simplest form of sub-space approximations, used in the work by (Sobieszczanski-Sobieski et al. 2001), Kodiyalam et al. (2004), Ollar et al. (2014) and Ryberg et al. (2015), deficiencies in sub-space partitioning, i.e. by failing to identify significant variables, can result in approximation errors that cannot be resolved by additional sampling. This is due to changes in response values as a consequence of changes in the insignificant variables. Such changes can not be included in the approximations as they are not defined in the space of the insignificant variables. Regardless of how carefully the partitioning of variables is made, there is always a risk that significant variables will be incorrectly identified as insignificant. Therefore a recovery mechanism for such errors is needed. This can be implemented in the trust region strategy by making sure that the vector of insignificant variables for the individual response is updated at the end of each iteration depending on the current best point as proposed by Ollar et al. (2015). Let x * k−1 denote the solution vector to the previous iteration (k − 1) for the optimisation problem (9). Recalling (4) and (5), the projection of the current best point onto the space of the significant variables can be written as and onto the space of the insignificant variables as Any change in the response values as a consequence of changes in the insignificant variables ψ * j will not be accounted for by the approximations and will mean that even if the approximations are of excellent quality in the space in which they are defined, there will be a discrepancy between the result they deliver and the one returned by the verification of the current best point. In this work this discrepancy is recovered in the following iteration by setting the constant value of the insignificant variables during sampling according to the current best solution according to where ψ k denotes the values of insignificant variables during sampling in the current iteration and ψ * (k−1) the current best solution in the previous iteration. Subscript j denoting the response number has been omitted for brevity. Figure 5 demonstrates how the value of ψ n changes from the previous iteration, k − 1, to the current iteration, k, for the two dimensional case. As the approximations in the current iteration are built with the values of the insignificant variables updated according to the current best solution, any changes in response values due to changes in insignificant variables from the previous iteration will be taken into account by the approximations in the current iteration.
With the possibility of significant variables being identified as insignificant, there is a risk that existing points located within the recycle region but far away from the current value of the insignificant variables may spoil the resulting approximation. The recycle-region size is therefore multiplied by a reduction factor b s , for the insignificant variables, as shown in Fig. 6.

Demonstration on a finite element model
This section describes a simple benchmark MDO study of a thin-walled beam structure subject to static and impact requirements. The impact requirements have disparate but overlapping variable dependence which makes it suitable for demonstration of the proposed sub-space approximation approach.

Beam model
The beam is 1.2 m long and is made of two thin walled hat sections welded together along two flanges as shown in Fig. 7. It is modelled using shell elements with an average mesh size of 10 mm. The welds are modelled using single Fig. 7 The thin-walled beam consisting of two hat sections spotwelded together along the flanges. Top part of structure has been detached at two locations hexahedron elements connected to the shell elements using a fixed contact formulation. The structure is divided into 26 components, shown in Fig. 8, with individual thickness values which are to be determined. The starting thickness is 3 mm for each component and the allowable thickness range is 1-5 mm.

Load cases
The beam is subjected to one static load case and three dynamic load cases as outlined below and shown in Fig. 9.
(a) Torsional stiffness A free node on the lefthand side of the beam is constrained in all translational degrees of freedom as well as from rotation around the longitudinal axis of the beam. The node is connected to the edges of component c25 (see Fig. 8) using one dimensional elements (RBE3) that distributes forces from the free node to the structure. At the right hand side of the beam a force is applied to a lever arm, connected to the beam with the same type of connection, resulting in a moment around the longitudinal axis of the beam. The response is the resulting rotation at the point where the moment is applied. The mass of the structure is also used as a response and is extracted from this load case. The load case is analysed using a static implicit solution procedure in Altair OptiStruct (Altair Engineering, Inc., 2014a). (b) Cylinder impact left (CL) A heavy cylinder with an imposed velocity impacts the left section of the beam in vertical direction. The beam is supported in vertical direction by a rigid plane. At impact the beam is compressed by the heavy cylinder. The response is measured as the vertical deformation of the top part of the beam, caused by the cylinder impact. The load case, which is a dynamic explicit load case, is analysed using Altair RADIOSS (Altair Engineering, Inc., 2014b). (c) Cylinder impact center (CM) The load case definition is the same as for (b) except for the impact position of the cylinder which is now at the centre of the beam.

Optimisation problem set up
The objective of this MDO problem is to minimise the mass of the beam subject to sufficient torsional stiffness and not exceeding maximum intrusion from the cylinder impact.
The target values are summarised in Table 1 together with initial values, which are normalised to target, for each response. The mass and torsion response have available gradients which are used during the optimisation to build gradient-enhanced approximations. As the mass response is linear in the space of the design variables a simple linear regression using the least squares method (LSM) is used to approximate this response. The remaining responses are approximated using the moving least squares method (MLSM). Four types of optimisations were carried out as outlined below. All optimisations were carried out using the MAM with 10 candidate points in each iteration obtained as the best out of 20 SQP optimisations starting from different points, sampled within the trust-region. The number of significant variables identified for each discipline and the corresponding number of points per iteration is presented in Table 2 for each optimisation. The convergence criteria was set such that the approximation quality must be below 5 % and the trust region size must be smaller than 10 % of the design region. The DOE method used for sampling of training points and start points for the SQP is based on a random number generator and will hence produce a different set of points depending on the initially chosen seed. In order to account for this uncertainty 50 optimisations with varying seed were carried out for each of the optimisation types described below.

Optimisation 1 -Full space approximations
The first optimisation is set up as a point of comparison using the MAM without sub-space approximations. The number of points evaluated per iteration is chosen as np = 1.5 × n = 39 per load case, leading to a total of 156 points per iteration.

Optimisation 2 -Same budget
In the second optimisation sub-space approximations with the proposed recovery mechanism is used within the MAM.
Partitioning of the beam model is shown in Fig. 10, where insignificant variables are shown transparent. The optimisation is carried out with the same computational budget per iteration as the first iteration in an attempt to increase the quality of the approximations and hence a greater chance of the trust region moving in the right direction early on in the optimisation, leading to fewer iterations required for convergence.

Optimisation 3 -Decreased budget
In the second optimisation sub-space approximations with the proposed recovery mechanism is used within the MAM. The partitioning of variables is the same as in the second optimisation. Here the number of points are determined individually for each discipline as ns = 1.5 × ns, where ns is the number of significant variables for the discipline. The aim of the optimisation is to keep the number of required evaluations per iteration to a minimum.

Optimisation 4 -Erroneous partitioning
In the final optimisation sub-space approximations with the proposed recovery mechanism is used within the MAM. The partitioning of variables is the same as in the second and third optimisation with one exception. Variable c7, shown in Fig. 11, which is significant to discipline CM, is deliberately identified as insignificant in order to test the proposed recovery mechanism outlined in Section 3.2.

Results
The results of the study is presented in Fig. 12. For each optimisation the median, upper and lower quartiles and the minimum and maximum value is shown for the number of  iterations, number of evaluations, the objective function and the maximum constraint violation. The first optimisation which was used as a comparison for the other three optimisations finished on average in 10.5 iterations, having spent 1407 evaluations, and with an average reduction in objective function of 13.2 %. The second optimisation, which was using sub-space approximations with the same number of evaluations per iteration as the first optimisation, finished on average in 8.5 iterations, 2 iterations less than the first, having spent 1242 evaluations, 165 less than the first optimisation, and with an average reduction in the objective function of 13.5 %, 0.3 % more than the first optimisation. The third optimisation, which was carried out using sub-space approximations with individual allocation of the number of points per iteration for each discipline, finished on average in 11 iterations, slightly more than the first optimisation but having spent 986 evaluations, 421 less than the first optimisation and 256 less than the second. The average reduction in objective function was 13.2 %, just like in the first optimisation. The fourth and final optimisation, where a significant variable was deliberately identified as insignificant, finished on average after 10.5 iterations, having spent 1052 evaluations, with an objective reduction of 12.6 %, 1 % less than the first optimisation.
It can be concluded from the results that, by using sub-space approximations for the presented example, it is possible to reduce the number of required iterations and/or the number of evaluations for carrying out the optimisation without compromising the results. By maintaining the number of points that are needed for full space optimisation when using sub-space approximations, as in optimisation 2, the number of iterations can be reduced. If instead, the number of evaluations per iteration is individually allocated, as in optimisation 3, for each discipline, a reduction in the number of evaluations can be reduced. It can also be included from the fourth optimisation that the proposed recovery mechanism makes sure that erroneous identification of significant variables does not lead to constraint violations, but instead to a slight penalty in the objective function.

Conclusions
An approach to carrying out multidisciplinary design optimisation using approximations built in individual subspaces for each discipline while carrying out the optimisation in the full variable space has been proposed. The sub-spaces in which the approximations are built are defined by the sets of significant variables for the individual disciplines. The main benefit of the technique is the dimensionality reduction of the approximations and sampling. This enables reducing the computational budget required to obtain approximations of sufficient quality.
The method relies on the designer to make assumptions on which variables are significant for each response. If such assumptions are deficient, approximation errors can occur that cannot be reduced by additional sampling. Therefore a technique was proposed that can recover from such errors within a trust region based optimisation framework. By updating the values of the variables identified as insignificant, but remain present in full space, according to the best current solution in each iteration, the technique can recover from such errors.
The approach was demonstrated on a finite element example of a thin-walled beam. A reduction in the computational budget was shown for optimisations carried out using the sub-space approximation approach compared to conventional optimisations. In addition, a test where a significant variable was deliberately identified as insignificant was carried out and showed the validity of the developed recovery mechanism.