2nd EUROPEAN MADYMO USERS' CONFERENCE, Stuttgart 1999
|
Astrid Watermann
|
Ralf Reuter
|
|
|
BMW AG
|
EASi Engineering GmbH
|
|
|
Abteilung: EG-221
|
Siemensstr. 12
|
|
|
D-80788 München
|
D-63755 Alzenau
|
|
|
|
|
|
Hier können Sie den Vortrag als PDF-Dokument downloaden:
|
A major challenge in Computer Aided Engineering is to produce reliable simulation results, which can be used for the prediction of realistic system behavior. As conventional simulation methods are based on a deterministic approach, they do not reflect uncertainties in the design parameters, i.e. scatter of input variables.
An easy way to make models more realistic is to use stochastic simulations, based on Monte-Carlo methods. The core of this method is to assume input parameters to be stochastic variables with a given realistic distribution. Then a number of simulation runs are conducted where these parameters are varied according to their distributions. The corresponding clouds of results show in which range the real system's response has to be expected, when physical uncertainties are considered. Furthermore, measures of variation can be computed, which enable the engineer to evaluate the distribution and to assess the reliability of the results.
Of course this valuable information is not free of charge. MC-methods require a significant amount of additional computational effort, as a considerable number of simulations is required to obtain such results. MADYMO simulations are particularly well suited for this type of simulations, as in most cases CPU time for a single model is fairly low.
This paper describes how to integrate stochastic simulations into the CAE process to manage uncertainties. It will be shown how stochastic simulations can be used for system design and optimization.
During the last years simulation models have grown considerably in terms of size. Models have become very detailed and simplifying modeling techniques (e.g. flat or scaled FE airbag models with approx. 1,000 elements) were replaced by more sophisticated ones (e.g. folded finite element airbag models with more than 10,000 elements). Certainly, simulation models have become more realistic as a result of this development towards a higher model resolution, which was mainly enabled by the availability of affordable high performance computers (Watermann/Altes 1996).
Today we seem to have reached a point, where an increase in model resolution will hardly improve the predictiveness of simulations. The main reason for this phenomenon is that the current models ignore certain parts of the underlying physics. Therefore, improvements rather result from improving, for instance, the thermodynamic airbag model, than from meshing an airbag with even more elements.
Another physical phenomenon that is not represented in typical simulation models, is the presence of scatter or uncertainty. Today's deterministic models simply ignore the fact that scatter is a natural part of any real-world physical system. Such models merely analyze one out of an infinite number of possible configurations of a system. In fact they typically analyze an ideal (or nominal) state of a system. The probability of the real system being in this particular state is extremely low. As Daniel S. Goldin, head of the NASA, put it: "Our models are limited, so we oversimplify the real world, relying on separate test programs to establish worst-case operating and failure conditions. To account for uncertainty and to quantify the risk level we must move from traditional deterministic methods to nondeterministic methods, ..." (Goldin/Venneri/ Noor (1998))
In some cases the analyst tries to model the worst case. For highly complex systems, such as a car, however, it is nearly impossible to determine in advance which actually is the worst case. The worst case for each component, which sometimes might be relatively easy to determine, is not necessarily the worst case for the whole system. Most engineers are aware of this fact. The consequence is to introduce high safety margins to make sure that the system will not fail in the actual worst case (which has not been analyzed!).
In the light of the above statements, it seems highly appropriate to include scatter in simulation models. Only by doing so, models actually simulate the real world physics and not an ideal state which is very unlikely to occur. Only this type of models will give the engineer enough confidence to use simulation results as a basis for major decisions. Of course, simulation experts are also aware of computational scatter that is added to the real world scatter (Marczyk 1999, p. 11). Considering this it is even more important to have confidence in numerical models before using them for decision making.
A straightforward approach to evaluate the consequences of scatter is the use of Monte-Carlo methods. The basic concept of this method is rather simple:
A sample of N sets of the scattered parameters of the model is created. Within the sample each parameter follows a predefined statistical distribution (e.g. a Gaussian distribution). In order to reach an agreeable degree of fit with the theoretical distribution, a certain number of sets is required. Typical sample sizes range between N=50 and N=150.
Each set of parameters is used to run a deterministic simulation of the model. These N jobs can be run in arbitrary order, as each is independent of all other jobs. This characteristic, which is called "trivial parallelism" (Marczyk, 1997, p. 5), allows to efficiently distribute the sample over a large number of CPUs, thus being able to run even larger samples within a reasonable time frame.
The result of this method are N solutions of the (stochastic) model. So, instead of a single deterministic result, one receives a population of solutions, which provides a deeper insight into the model. The results can be analyzed with statistical methods. Histograms show the distribution of results. Rather than a single point, a range of possible results is analyzed. The following figure depicts a typical case.

Figure 1: Histogram
Assume that for a certain vehicle a limit for the HIC value is set to 660 (indicated by the dashed line in fig. 1). The cross in Figure 1 shows the result of the deterministic model (HIC=555). A Monte-Carlo study for this model is conducted. The HIC distribution (histogram) is shown in Figure 1. While the deterministic result would indicate that the system performs well, as far as the given limit is concerned, we learn from the stochastic model that (with the given scatter of input parameters) there is a significant chance that the HIC is above the limit. Before accepting this system one would certainly try to reduce the probability to violate the limit.
Therefore stochastic simulations, which deliberately take scatter into account, are a means for the engineer to get confidence in his models and his decisions. Only when the failure probabilities are known, serious decisions are possible.
In addition to establishing result distributions and failure probabilities, the Monte-Carlo method is a useful tool for various engineering tasks, such as parameter studies, sensitivity analyses, optimization etc. The subsequent chapters of this paper show some applications of the method to MADYMO occupant simulations. MADYMO simulations are very well suited for Monte-Carlo studies, as typical MADYMO jobs are relatively small (<2hrs.). So, it is possible to run studies with up to a few hundred shots in reasonable times. An industrial tool to perform this studies is the ST-ORM (STochastic Optimization and Robustness Management) package. It can be used in conjunction with any solver that uses ASCII input files and result files. It allows the user to assign distribution functions to any variable in the input file, it automatically creates and submits the jobs for the Monte-Carlo study and evaluates the results statistically.
One major question in engineering problems is to identify those parameters which mainly influence the performance of the system under investigation. Such parameters may either be design parameters or uncertainties (scatter) in the system. Frequently, only the design parameters are taken into account. However, in some cases it is possible that the uncertainties are even more important for the performance. This means that even if one finds a good solution for the design parameter, the performance may still be poor because of the scatter in the system. In such cases it is important to work on reduction of the scatter rather then on improvement of the design parameters. Only if the scatter is under control (e.g. by decreasing the tolerances), the design parameters will deploy their potential.
In other cases, one may find that a particular parameter does not influence the system's performance at all. In such cases the tolerances may be relaxed, which might have significant economic effects.
The first example shows a sensitivity analysis. A relatively simple MADYMO model was used for this purpose (Figure 2). For analysis of the scatter the following parameters have been varied:
Normal distributions with realistic standard deviations were applied to the parameters.

Figure 2: MADYMO Model
The results were evaluated both graphically and analytically. Graphical analysis can be done by means of scatterplots (ant-hill plots). Such plots show the relationship between 2 (or 3) parameters. Figure 3 depicts an example.

Figure 3: Scatterplot - HIC(Hole Size)
The scatterplot shows a dependency of the HIC on the size of the airbag hole. With increasing hole size the airbag becomes softer and the HIC decreases. This effect is independent from small variations of the other parameters, meaning there is a significant influence of small changes of the airbags hole size.
Correlation coefficients measure the degree (and direction) of a relation. The most frequently used correlation coefficient is the Bravais-Pearson coefficient. This coefficient, however, only measures linear relations. A coefficient that is capable of measuring non-linear (but strictly increasing or decreasing) relations is Spearman's rank correlation coefficient rs.

The rank of xi equals the number of values xn for which xn < xi, plus 1. The following table illustrates the calculation of ranks and mean ranks.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Table 2 shows the correlation matrix with the rank correlation coefficients for all the independent variables.
|
|
|
|
|
|
H-point (x)
|
|
|
|
|
H-point (z)
|
|
|
|
|
Vent size
|
|
|
|
|
Friction
|
|
|
|
|
Firing time
|
|
|
|

Figure 4: Parameters influencing HIC
Figure 4 shows the influence of the examined parameters on the HIC. It can be seen that in this case the HIC is mainly influenced by the size of the airbag hole as well as the H-point position (x). However, one has to recognize that in this simplified example not all uncertain parameters have been investigated, and there might be more parameters influencing the result.

Figure 5: NCAP Rating
The scatter of the mentioned variables influences the NCAP rating of the car as depicted in Figure 5. Even though the nominal model achieved a 5 star rating, the stochastic analysis shows that the result cloud crosses the border to a 4 star rating. This shows that for the regarded example there is a significant chance that only 4 stars will be reached in the actual test. It is therefore important to reduce the result scatter to ensure a 5 star rating.
The results described above help identify those parameters which mainly cause this scatter. These parameters need to be better controlled. For example, to reduce the HIC scatter it is necessary to minimize the scatter in the airbag hole-size and the scatter in the H-point position (x) (Figure 4). While the first parameter requires improvement of the manufacturing process or of the airbag fabric, the second parameter is a test condition, which, if testing is done externally (as with the NCAP test), is out of the manufacturers' control. Therefore, in this case the goal would be
From this example one can summarize the following
Another task which can be approached with stochastic simulation is the multi-parameter optimization of systems. In the example presented here, the task was to optimize the characteristic of a belt load limiter for a rear passenger. This case is relatively simple, but it gives a good idea of the method that was used.
Figure 6 shows the simplified characteristic that was optimized. The first part of the function has a parabolic shape which remained unchanged. The load level at which the belt yields is the first parameter (P1) that was investigated. Two further points (P2,P3) describe the load limiter characteristic.

Figure 6: Load Limiter Characteristic
The objective was to find a characteristic that performs well for both 50th%ile and 95th%ile occupants. Therefore, a model was built by merging two existing models into one single model, containing both dummies independently (Figure 7). This enabled the simultaneous evaluation of both cases.

Figure 7: Merged Models
For each occupant size an objective function was established:

The two objective functions were combined:
![]()
Several constraints were also established. Among these was the condition that a minimum clearance between the head and the back-rest of the front seat was required
( dist(Head - frt. Seat)>cmin ).
The optimization was carried out in several steps. Initially, large ranges for each parameter were investigated. For better visualization, the results can approximated by surfaces. Figure 8 shows the results (OFcomb) plotted over the load values at points P1 and P2.
One can clearly identify regions with good performance in Figure 8. Parts of these regions were not acceptable because some constraints were violated. In order to establish adequate ranges for a second optimization step, the ten best solutions which violated no constraint were evaluated. Thus, the parameter ranges could be narrowed down for a finer optimization. The dark fields in Figure 8 indicate the reduced search region for this optimization step. This fine tuning delivered an optimized set of parameters.

Figure 8: Result Surface
While normally an optimization would stop at this point, we were clearly aware of the fact that this optimal solution was not necessarily a robust one. It is apparent from testing that a real load limiter does not follow exactly the ideal load curve that we found. Therefore it was necessary to find a design point at which the system still performs well under deviations from the ideal state. The robustness of the optimized model was evaluated by introducing scatter. A random deviation of ±0.1 kN (uniformly distributed) at each point was assumed. With this - relatively small - amount of scatter, it was found that the probability of failure (=violation of one or more constraints) was approximately 18%. The reason for failure in all cases was a violation of the minimum clearance constraint. As this was not acceptable, the load values at each point were slightly increased. A second evaluation of robustness with this characteristic showed a low probability of failure (no failure in 50 runs with a random deviation of ±0.1 kN).
Compared to the baseline model, the mean result (OFcomb) of the robust optimized model was reduced by 35%.
This example shows the main advantages of the method: confidence in and insight into the model. It does not only find the optimal value, it also shows the environment at this point. This case showed that even small deviations could cause violation of constraints. Being aware of this, and being aware that there will always be some scatter in a real world system, it is clear that the "theoretical" optimal point is not a good operating point. Considering this, it makes more sense to define a somewhat sub-optimal operating point, including tolerances which ensure a low probability of failure.
Typically at this point a trade-off has to be made between system performance and tolerances. Obviously, by reducing the tolerances the operating point comes closer to the optimal point. Reducing tolerances, however, normally increases manufacturing costs and in some cases there are even technical limitations to the reduction of tolerances.
In this situation the method described here helps identifying those parameters that determine the system performance and therefore should be kept within narrow tolerance bands. At the same time those parameters that have little or no influence on the performance are also identified. The tolerances for these parameters can then be relaxed, which may lead to significant cost savings. In our example we found that the load at point P3 had little influence on the results. While the results were clearly correlated to the loads at P1/P2 (correlation coefficients 0.73/0.66), the correlation to the load at P3 was low (0.4). Therefore we tripled the tolerance at P3, finding that this did not worsen the results.
Figure 9 depicts a detail of the optimized load curve. It shows the optimal solution as well as the operating point. The shaded area indicates the tolerance at the operating point.

Figure 9: Optimized Characteristic
Compared to other optimization methods, the stochastic approach provided additional information that was highly valuable for translating simulation results into practical solutions.
Stochastic simulations, based on Monte-Carlo methods as outlined in this paper, are an important supplement to conventional deterministic simulation methods. Due to limited hardware resources it is not yet feasible to entirely replace deterministic analysis.
But considering the valuable additional information that the method yields, it is worthwhile to spend some additional computation time. Especially in the (nearly) test-free development processes that are targeted throughout the automotive industry it is becoming increasingly important to manage uncertainty in computer simulation models. Managing uncertainty means to actively address the fact that uncertainty exists in all physical systems, instead of ignoring it or avoiding its consequences by designing systems with large safety margins.
As the second example showed, the method has the potential to save manufacturing costs by relaxing tolerances or by reducing overly high safety margins. In the case of mass production, as in the automotive industry, even small reductions in manufacturing costs for a single product will result in large savings which easily outweigh additional costs for doing a few hundred simulation runs.
Goldin, D.S./Venneri, S.L./Noor, A.K.: Beyond Incremental Change. In. Computer, Vol. 31, No. 10, October 1998, pp. 31-39.
Marczyk, J.: Statistical Mechanical Designs - Uncertainty Management in CAE via Monte Carlo Simulation. In: BENCHmark 1/1999, pp. 11-15.
Marczyk, J.: Meta-Computing and Computational Stochastic Mechanics. In: Marczyk, J. (Ed.): Computational Stochastic Mechanics in a Meta-Computing Perspective, Barcelona 1997, pp. 1-18.
Watermann, A./Altes, J.: Finite Element Calculations on the Parallel Computers Intel/Paragon, IBM SP-2, and a Cluster of DEC-ALPHA Workstations. In: Proceedings of the 3rd International Conference on High Performance Computing, Trivandrum, India, IEEE Computer Press, 1996, pp. 88-94.
Last updated: 8. Juli 1999
© EASi Engineering GmbH 1999