Meccanica (2010) 45: 693–704

DOI 10.1007/s11012-010-9285-0

S I M U L AT I O N , O P T I M I Z AT I O N & I D E N T I F I C AT I O N

On the Weibull cost estimation of building frames designed

by simulated annealing
Ignacio Paya-Zaforteza · Víctor Yepes ·
Fernando González-Vidosa · Antonio Hospitaler

Received: 4 March 2009 / Accepted: 21 January 2010 / Published online: 16 February 2010
© Springer Science+Business Media B.V. 2010

Abstract This paper proposes a general methodology 1 Introduction

to determine the number of numerical tests required to
provide a solution for a heuristic optimization problem Initial attempts to solve structural optimization prob-
with a user-defined accuracy as compared to a global lems date back to the 1600s, when Leonardo da Vinci
optimal solution. The methodology is based on the ex- and Galileo conducted tests with models and full-scale
treme value theory and is explained through a problem structures [1]. Modern optimization is linked to the
use of computers. The spectacular increase in compu-
of cost minimization for reinforced concrete build-
tational resources available to engineers in the 1980s
ing frames. Specifically, 1000 numerical experiments
greatly improved the structural analysis of complex
were performed for the cost minimization of a two-bay
structures, especially as regards linear elastic analysis
and four-floor frame using the Simulated Annealing of bar structures and 2-D and 3-D finite element mod-
(SA) algorithm. Analysis of the results indicates that eling. Structural analysis software was upgraded with
(a) a three-parameter Weibull distribution function fits the many computer-aided design (CAD) programs de-
the results well, (b) an objective and general procedure veloped in the 1990s. CAD programs allowed users
can be established to determine the number of experi- to tentatively optimize the structure by trial-and-error,
ments necessary to solve an optimization problem with but they did not include any design variables subject
a heuristic which generates independent random solu- to optimization. Current research focuses on structural
tions, and (c) a small number of experiments is enough optimization programs, like the one used in the present
to obtain good results for the structural engineer. study. These programs build on CAD programs and
define the structure in terms of design variables which
the optimization algorithm must modify in the search
Keywords Optimization · Reinforced concrete · for the optimum structure. Despite the potential ca-
Weibull distribution · Extreme value theory pabilities of structural design by optimization, eco-
nomic design of concrete structures is at present highly
conditioned by the experience of the structural engi-
neer. Basically, designing concrete structures is a prob-
I. Paya-Zaforteza () · V. Yepes · F. González-Vidosa · lem of selecting design variables subject to structural
A. Hospitaler constraints. Design variables include material grades,
ICITECH, Departamento de Ingeniería de la Construcción
y Proyectos de Ingeniería Civil, Universidad Politécnica
cross section dimensions and reinforcement. Most cur-
de Valencia, Camino Vera S/N, 46022 Valencia, Spain rent design office procedures still adopt cross-section
e-mail: igpaza@upv.es dimensions and material grades based on sanctioned
694 Meccanica (2010) 45: 693–704

common practice. Once the structure is defined, it fol- times the algorithm should be run to achieve sufficient
lows the CAD structural analysis and computation of accuracy: 10, 50, 100, 200 or 500 times. Most pub-
passive and active reinforcement. Should the dimen- lished papers on structural optimization do not even
sions or material grades be insufficient, the structure is report the number of runs performed by their algo-
redefined on a trial-and-error basis. Such process leads rithms, and in most cases researchers simply report
to safe designs, yet the cost of concrete structures is, best predictions and compare algorithms on the basis
therefore, very much linked to the experience of the of these best predictions. A second crucial question is
structural designer. how close the results are to the unknown global opti-
There are two main methods for structural opti- mum. In practice, it is necessary to achieve a balance
mization. Exact methods are generally based on math- between the quality of the results and the amount of
ematic programming (Hernández and Fontan [2]), computing time required to obtain them. Further, it is
while heuristics are artificial intelligence strategies necessary to assess how close the heuristic method so-
usually emulating natural processes. Heuristic meth- lutions are to the global optimum, given that the global
ods include a large number of algorithms, such as ge- optimum solution is usually unknown for the problems
netic algorithms (GA henceforth), tabu search, sim- that justify the use of a heuristic method. Any solution
ulated annealing (SA), ant colony optimization and given by a heuristic stochastic procedure is generally
neuronal networks, among others [3–6]. Heuristics a good quality solution regarding the overall solution
usually provide near-optimal solutions. Among the
space. Each independent solution of a stochastic al-
first studies on heuristic optimization applied to struc-
gorithm can be considered a rare event, and we can
tures, the contributions of Goldberg and Samtani [7]
construct a rare event distribution by multiple runs of
and Coello et al. [8] are worthy of mention. These
the algorithm. The occurrence of rare events is studied
pioneering authors applied GA to the weight opti-
in this paper by the extreme value theory (EVT hence-
mization of a steel truss structure and to the cost of
forth). The use of the EVT to estimate the global opti-
a simply supported reinforced concrete (RC) beam,
mum solution to a heuristic method has been described
respectively. Recent research using GA strategies can
previously in McRoberts [25], Golden and Alt [26],
be found in Panigrahi et al. [9], Bassir et al. [10] and
Sid et al. [11]. Regarding concrete structures, the use and more recently in Bettinger et al. [27] for large
of non-evolutionary methods as an alternative to GA combinatorial problems.
is rare, but applications can be found in Balling and In this line of study, this paper describes a method
Yao [12] and Ceranic et al. [13]. Recently, our research to define the minimum number of computer runs an
group has also used non-evolutionary algorithms such optimization algorithm must be performed to ensure
as SA, threshold acceptance and ant colonies to op- that the best result differs to a certain degree from the
timize earth retaining walls, frame bridges, building global optimum estimation. One thousand numerical
frames, road vaults and bridge piers [14–18]. Further, runs are carried out so that the cost of a RC building
the multiobjective optimization of building frames frame was optimized by SA. The EVT is applied to
has been examined by the authors [19]. Kicinger the results obtained. The paper is organized as follows.
et al. [20] and Adeli and Sarma [6] have recently pub- Section 2 specifies the optimization problem, that is,
lished surveys on structural evolutionary optimization structural system analyzed, objective function, vari-
of steel and RC structures. Other recent applications ables, parameters, and optimization algorithm. Sec-
of structural optimization to civil, naval, mechani- tion 3 describes the three-parameter Weibull statistical
cal and aerospace engineering are discussed in Nieto distribution function. This distribution fits the results
et al. [21], Marannano and Mariotti [22], Callegari and of the optimization problem very well as described in
Palpacelli [23], and Guadagni [24]. Sect. 4. A methodology based on the estimation of the
Heuristic methods are not deterministic since they location parameter of the three-parameter Weibull dis-
provide a different result every time they are run in tribution is proposed as well in Sect. 4 to determine the
the computer. This is so because they include a large minimum number of computer runs necessary to solve
number of iterations governed by random decisions. the optimization problem with a required accuracy. Fi-
In these terms, they provide as many different results nally, the main conclusions of the work are drawn in
as algorithms run. The key question then is how many Sect. 5.
Meccanica (2010) 45: 693–704 695

2 Optimization problem units (concrete, steel, formwork, etc.) to calculate the

cost of a RC frame. Each unit has a price pi and a mea-
2.1 Structures studied surement mi . Thus, the cost of the frame is the sum of
the unit prices multiplied by the measurement. Con-
This example is related to RC frames commonly used struction units considered in this study and their prices
in building construction. Building frames (see Fig. 1) are detailed in Table 1. Constraints in expression (2)
typically have horizontal beams measuring 5.00 to are all the service and ultimate limit states that the
10.00 m in horizontal span. These beams sustain the structure must satisfy as well as the geometrical and
vertical loads of the floors and transfer them to vertical constructibility constraints of the problem. The com-
columns measuring 3.00 to 5.00 m in height. Moder- ponents of the X  vector are the design variables of the
ate horizontal loads are usually included in the design, optimization problem.
but high levels of horizontal loading are usually trans-
ferred to adjacent shear walls. Building frames are cal- 2.3 Design variables, parameters and constraints
culated to sustain the loads prescribed by the codes and
must satisfy all the limit states required as a RC struc- Figure 2 illustrates the typical symmetrical frame of 2
ture. bays and 4 floors which will be used for the analy-
ses in Sect. 4. This frame has a total of 77 design
2.2 Objective function variables which have been previously described in de-
tail in [19, 28]. All variables in this analysis are dis-
The structural concrete design problem established in crete. These variables include (1) the type of reinforc-
this study involves an economic optimization, which ing steel which has a characteristic yield stress of 400
aims to minimize the objective function F of expres- or 500 MPa, (2) the type of concrete for the columns
sion (1), satisfying as well the constraints of expres- and beams of each floor, which may vary from a min-
sions (2): imum characteristic compressive strength of 25 MPa
 to a maximum of 50 MPa in steps of 5 MPa, (3) the
 =
F (X) 
pi · mi (X), (1) cross section dimensions of each column of the frame,
i=1,r (4) the cross-section dimensions of the beams, and
 ≤ 0.
gj (X) (2) (5) the reinforcement steel variables which follow a
standard reinforcement setup. Figure 3 shows a typi-
The objective function in expression (1) is the cost cal longitudinal reinforcement setup of the beams of
of the structure. The structure is split in r construction the structure. It includes basic top and bottom bars in

Fig. 1 Structural systems analyzed in this study: 3D view (left) and real structure during construction (right)
696 Meccanica (2010) 45: 693–704

Table 1 Unit prices considered for reinforced concrete framed structures in this study

Unit Descriptiona Cost (€)

kg Steel B400S 1.27

kg Steel B500S 1.30
m3 Concrete HA-25 78.40
m3 Concrete HA-30 82.79
m3 Concrete HA-35 98.47
m3 Concrete HA-40 105.93
m3 Concrete HA-45 112.13
m3 Concrete HA-50 118.60
m2 Formwork in beams 25.05
m2 Formwork in columns 22.75
m2 Scaffolding for beams 38.89
a Notation for steel and concrete follows the Spanish Code for Structural Concrete [29] and is explained in Sect. 2.2

factors. Structural constraints followed standard pro-

visions for the design of this type of structure [29, 30],
which include checks of the service and ultimate limit
states of flexure, shear and instability for the stress
envelopes due to the vertical loads and the horizon-
tal wind loads. Vertical loads amount to a total uni-
form distributed load of 35 kN/m (7.00 kN/m2 of self
weight plus live load and 5.00 m of spacing between
parallel frames). The total uniform distributed load of
wind loads is 4.5 kN/m. Stress resultants and reactions
were calculated by an internal matrix method program
using a 2-D mesh. Deflections were limited to 1/250
of the horizontal span for the total load and to 1/400
for the active deflection, which is the part of the de-
Fig. 2 Typical RC building frame with 2 bays and 4 floors flection measured after constructing the elements that
can be damaged by vertical displacements.
addition to positive and negative extra reinforcements
of standard length. Variables for beam stirrups include 2.4 Simulated annealing optimization algorithm
3 zones in the left, central and right positions of trans-
verse reinforcement. The longitudinal reinforcement The search method used in this research is a version of
of columns includes 330 possible values which vary the single objective SA proposed originally by Kirk-
from a minimum of 4φ12 to a maximum of 34φ25, patrick et al. [31] and independently by Cerny [32].
whereas the transverse reinforcement of columns in- These researchers showed how a model for simulat-
cludes 21 possible values. It is worth noting the com- ing the annealing of crystals could be used for opti-
plexity of optimizing the cost of this structure given mization problems, in which the objective function to
the large number of variables (77) and the size of the be minimized corresponds to the crystal energy states.
solution space of the problem (10116 ) as detailed in The SA algorithm is based on the analogy of crystal
Payá-Zaforteza [28]. formation from masses melted at high temperatures
The main parameters are the horizontal spans of the and left to cool slowly. At high temperatures, con-
bays, the vertical height of the columns, the vertical figurations of greater energy than previous ones may
and horizontal loads considered and the partial safety form randomly, but as the mass cools, the probability
Meccanica (2010) 45: 693–704 697

Fig. 3 Typical longitudinal reinforcement bars of the beams in RC building frames

of higher energy configurations forming diminishes. 40%, the initial temperature is halved; and if they are
The process is governed by the Boltzmann expression less than 20%, the initial temperature is doubled.
exp(−E/T ), where E is the increment in the en- The SA algorithm used in this research was coded
ergy of the new configuration and T is the tempera- in Compaq Visual Fortran Professional 6.6.0. The run-
ture. The algorithm starts with a feasible solution ran- ning time of one execution of the algorithm on a
domly generated and a high initial temperature. The PC Pentium IV 3.20 GHz computer was, on average,
initial working solution is changed by a small random 25 minutes. The calibration of the algorithm recom-
move of the values of the variables. The new current mended Markov chains of 70000 iterations, a cooling
solution is evaluated in terms of cost. Greater cost so- coefficient of 0.80 and two Markov chains without im-
lutions are accepted when a 0 to 1 random number is provement in the current solution as stop criterion. The
smaller than the expression exp(−E/T ), where E most efficient move identified was a random variation
is the cost increment and T is the current temperature. of 3 or up to 3 variables of the 77 in the problem.
The current solution is then checked against structural
constraints and if it is feasible, it is adopted as the new
working solution. The initial temperature is decreased 3 The Weibull distribution function
geometrically (T = kT ) by means of a cooling coeffi-
cient k. A number of iterations called Markov chains The Weibull cumulative distribution function (cdf
is allowed for each step in temperature. The algorithm henceforth) is given by:
stops when there are no improvements in the current
solution after a certain number of Markov chains. The FX (x0 ) = Prob{X ≤ x0 }

SA method is capable of surpassing local optima at 1 − exp{−( x0 η−γ )β }, x0 > γ ,
high-medium temperatures, and it gradually converges = (3)
as the temperature reduces to zero. The SA method 0, x0 ≤ γ ,
requires the calibration of the initial temperature, the
where
length of the Markov chains, the number of Markov
chains without improvement to stop the algorithm, and η, β > 0 (4)
the cooling coefficient. The initial temperature is ad-
justed following the method proposed by Medina [33], where γ is the location parameter, η is the scale para-
which consists in choosing an initial value and check- meter, and β is the shape parameter.
ing whether the acceptances of higher energy solutions The distribution was developed by the Swedish
fall between 20 and 40 percent. If they are greater than physicist Waloddi Weibull [34] to describe the strength
698 Meccanica (2010) 45: 693–704

behavior of materials. It is well-known that the Weibull €3503.48, with a confidence interval of ±€2.654 for
distribution belongs to the extreme value distribu- a 0.05 level of significance; the standard deviation of
tions. It represents the distribution of the smallest the sample is €42.82; the median is €3,495.67; the
or largest values in random samples of increasing kurtosis coefficient is €90.64, and the skewness co-
size. Suppose a number of independent samples of efficient is €6.67. Therefore, the sample has a posi-
size m is taken from a parent continuous population tive asymmetric distribution and a leptokurtic distrib-
which is bounded from below by γ . Fisher and Tip- ution (a high degree of values concentrated around the
pett [35] proved that, as m grows, the distribution mean). There is no reason to rule out the hypothesis
of min(X1, X2 , . . . , Xm ) approaches a three-parameter that the histogram corresponds to the probability den-
Weibull distribution with γ as the location parame- sity function of the Weibull distribution.
ter. This is applicable to our study since the result of
each computer run has the lowest cost from a large 4.3 Three-parameter Weibull distribution fitting
sample of solutions considered throughout the opti-
mization search. Regarding continuity and although This section checks the hypothesis that the 1000 SA
the building frame studied has an extremely high but results fit a three-parameter Weibull distribution func-
finite number of potential solutions, it is important to tion. Three conditions must be fulfilled to confirm this
note that we assume that the population space approx- hypothesis. Firstly, it must be verified that there is no
imates a continuous space well since the number of reason to reject the Weibull hypothesis. Secondly, it
variables is large and it is represented in an almost must be shown that the 1000 solutions of minimal cost
continuous manner. found by the SA algorithm are independent (see Fisher
and Tippett [35]). Finally, the correlation coefficient of
the Weibull distribution that best fits the 1000 numeri-
4 Results and discussion cal results must be high enough.
Kolmogorov-Smirnov and Chi-squared statistics
4.1 Overview (see, e.g., Conover [36]) were computed assuming in-
dependence to verify that the SA solutions are Weibull
If the statistical distribution of the local optima given distributed. These tests are non-parametric tests used
by SA fits a three-parameter Weibull distribution, then to compare a sample with a reference probability dis-
the estimated location parameter γ can be used as an tribution. These statistics fell far below the critical
estimation of the global optimum. In addition, a num- value at the 0.05 level of significance in both cases,
ber of computer runs based on the difference between indicating that there is no reason to reject the Weibull
γ and the best solution found by the algorithm can be hypothesis.
proposed. Section 4.2 describes the statistical proper- One of the main assumptions underlying the use of
ties of the sample of solutions obtained with 1000 runs EVT is that each SA solution can be considered an in-
of the SA algorithm. Section 4.3 proves that the three- dependent sample from the whole population of 1000
parameter Weibull distribution fits the 1000 results samples, because each search process starts with a dif-
well. Finally, in Sect. 4.4 we propose a method to de- ferent, randomly-defined solution. A Wald-Wolfowitz
fine the number of runs required for a given accuracy. runs test was applied to the 1000 SA solutions taken in
order of occurrence so as to confirm that SA best so-
4.2 Statistical description of 1000 SA results lutions are independent. This is a non-parametric test
that may be used as a test for randomness in a se-
To verify the applicability of the EVT to structural op- quence, and it is based on the total number of runs and
timization, the cost of the frame in Fig. 2 was opti- the number of cases on the same side of a cut point
mized 1000 times with the SA algorithm described in (see, e.g., Conover [36]). This test dichotomizes the
Sect. 2.4. Figure 4 shows the histogram of the sam- sequence of observations by classifying each observa-
ple of 1000 minimal cost solutions found by SA. The tion as either above or below the sample median to
statistical description of this sample is the follow- ensure that there is no order to the resulting sequence.
ing: the maximum and minimum values are €4224.19 The 2-tailed significance value of this randomness test
and €3460.00, respectively; the sample mean value is is the probability of obtaining a Z statistic as extreme
Meccanica (2010) 45: 693–704 699

Fig. 4 Histogram showing 1000 cost minimization results

as or more extreme (in absolute value) than the ob- parameter estimation and the rank regression on Y
tained value, if the order of ratings above or below the estimation (according to the least squares principle,
median is purely random. In this case, the total num- which minimizes the vertical distance between the
ber of cases is 1000 and the 2-tailed significance value data points and the probability density function) gave
is 0.849. This p-value does not allow us to reject the a γ = €3435.40 value for the location parameter. This
null hypothesis that the ratings are randomly ordered.
value is the estimation that the EVT provides for the
Hence, there is no indication that the total number of
global optimum of the problem using the results from
runs above and below the median is less than that re-
sulting from a random sequence. Therefore, the runs the SA algorithm. Figure 5 shows a representation of
test does not give any reason to reject the hypothesis the cdf obtained from the results of the numerical ex-
of randomness. periments and the Weibull cdf estimated by means of
Finally, the parameters of the Weibull distribution rank regression on Y (γ = 3435.40, η = 75.6579, and
that best fit the 1000 numerical tests must be obtained, β = 2.5764). The Weibull fit has a correlation coef-
and the suitability of the fit must be quantified to deter- ficient of ρ = 0.9343, which is quite high for numeri-
mine how well it satisfies the given data. Several pa- cal results. The difference between the minimum value
rameter estimation methods are available, such as the obtained with the 1000 SA runs and the extreme value
method of moments, maximum likelihood, minimum
estimated using the Weibull distribution is €60.27,
chi-square, least squares, etc. (see Dannenbring [37],
hardly 0.72% of the theoretical minimum value. From
Golden and Alt [26], and Vasko and Wilson [38]). In
this work ReliaSoft’s Weibull++7 software [39] was the point of view of the structural engineer, this indi-
used to estimate the three parameters of the Weibull cates that the difference is small enough to make the
distribution function. Both the maximum likelihood solution provided by the proposed SA acceptable.
700 Meccanica (2010) 45: 693–704

Fig. 5 Theoretical 3 parameters Weibull c.d.f. versus experimental values

4.4 Number of runs criterion γ (γ1 , . . . , γ9 ) were calculated. Additionally, each se-
ries of these nine γ values was ranked, starting from
The number of times the SA algorithm is run should a minimum value (γmin ) and ending with a maximum
be large enough to ensure that the difference between value (γmax ). Figure 7 and Table 2 show the values of
the minimum value obtained from all runs and the ex- γmin , γmax and Cmin . It can be inferred that the vari-
treme value estimated using the theoretical distribu- ability of the γ parameter is given by the difference
tion is less than a given threshold. Then, a new prob- between the estimations of its value provided by γmin
lem arises: the estimation of the γ parameter involves and γmax . This difference is 3.66% of γmin for 5 SA
variability, since its value depends on the number of runs, 1.53% for 9 SA runs, and only 0.22% for 15 SA
run solutions used to estimate it. Eight series of nine runs. The reduction in the difference is practically neg-
samples with replacement were extracted from the ini- ligible from 15 runs on. Secondly, the percentage dif-
tial population of 1000 solutions to study the relation- ference between γmin and Cmin is minimal for 15 SA
ship between the estimated value of γ and the number runs (0.69%) and for 25 SA runs; it is not possible to
of data used for the Weibull fit. Each one of these se- reduce the difference in question (0.72%).
ries corresponded with the results of 5, 9, 15, 25, 50, Therefore, the number of algorithm runs can be
100, 500 or 1000 SA algorithm runs. The first series based on the fulfillment of two conditions. Firstly, the
contained nine samples and each sample was com- difference between the best solution found to date and
prised of the results of five SA runs (see Fig. 6); the the theoretical global optimum estimated by adjusting
second series had nine samples each one containing a three-parameter Weibull cdf to the results obtained,
the results of nine SA runs and so on. Therefore, 72 which must be less than a given level. Thus, the op-
(8 × 9) samples were analyzed. First of all, the min- timum found is guaranteed to be sufficiently close to
imum cost (Cmin ) solution for each one of the eight the theoretical global optimum. Secondly, the value of
series and the γ values corresponding to the Weibull the difference γmax − γmin is less than another prefixed
cdf which best fit each one of the 72 samples were value, so the error in the estimation of γ is limited.
obtained. Therefore, for each series of 5 to 1000 al- Considering that when performing the algorithm, only
gorithm runs, one value of Cmin and nine values of those solutions obtained from previous performances,
Meccanica (2010) 45: 693–704 701

Fig. 6 Sample with replacement procedure for estimating minimum cost (Cmin ) and location parameters (γmin , γmax ) for nine samples
of five SA algorithm runs

Fig. 7 Value of the minimum costs (Cmin ) and the location parameter (γ ) for different number of SA algorithm runs

rather than the set of all the solutions, are available, the nique has been successfully applied to solve problems
use of the bootstrap technique (Efron [40]) is proposed that would be too complicated for traditional statisti-
for estimating the γmax and γmin values. The boot- cal techniques or to situations in which the classical
strap is a powerful technique for assessing the accu- techniques are not valid (Zoubir and Boashash [41]).
racy of a parameter estimator. If a representative sam- The bootstrap technique involves processing the sam-
ple is available, the bootstrap randomly reassigns the ple data as though they constituted the data of the en-
observations and recalculates the estimate. This tech- tire population. In other words, they are used as the
702 Meccanica (2010) 45: 693–704

Table 2 γmin , γmax , Cmin and the percentage differences for different numbers of SA runs

Number of tests γmin γmax (γmax − γmin )/γmin (%) Cmin (Cmin − γmin )/γmin (%)

5 3330.12 3452.11 3.66 3462.15 3.96

9 3403.17 3455.39 1.53 3462.15 1.73
15 3437.52 3445.22 0.22 3461.33 0.69
25 3436.72 3445.58 0.26 3461.33 0.72
50 3435.40 3440.16 0.14 3460.00 0.72
100 3435.40 3438.18 0.08 3460.00 0.72
500 3435.40 3436.86 0.03 3460.00 0.72
1000 3435.40 3435.85 0.01 3460.00 0.72

universe from which samples with replacements will

be extracted. For the example given, and together with
the original sample, eight additional random samples
with replacement were extracted and used to estimate
each one of the γ parameters.
For this research, the following methodology is
proposed. A set of prior experiments is carried out
using a certain heuristic, thereby permitting an ini-
tial sample of results to be obtained. From this sam-
ple, the theoretical minimum corresponding to a three-
parameter Weibull distribution is estimated. However,
to increase the significance of the estimator in ques-
tion, new samples are extracted with the bootstrap
Fig. 8 Bootstrap technique to estimate minimum cost (Cmin )
technique (see Fig. 8). These samples are then used and location parameters (γmin , γmax ) for nine samples of five
to calculate the values of γmin , γmax and Cmin . If SA algorithm runs
γmax − γmin and Cmin − γmin are lower than two pre-
viously set limits (for example, less than 0.5% of γmin 5 Conclusions
and 1% of γmin , respectively), additional runs of the
algorithm are not necessary. Otherwise, a new numer- The following conclusions can be drawn from the re-
ical experiment is carried out with the heuristic and the search described herein:
procedure is repeated. Figure 9 and Table 3 show the
results of the solutions found for the building frame • The set of solutions found using SA for the cost op-
timization of a 2-bay and 4-floor frame is a sam-
studied.
ple of independent values that the three-parameter
Employing the aforementioned criterion and mak-
Weibull distribution fits well.
ing a linear interpolation based on the data featured in
• The extreme value theory can be used to estimate
Tables 2 and 3, it is possible to deduce that the min-
the global optimum achievable via a heuristic pro-
imum number of performances is: (a) 14, if 1000 re- cedure such as SA when minimizing the cost of a
sults are used (something that cannot be done in prac- concrete structure.
tice) or (b) 18, if the bootstrap technique is used. In the • An objective procedure can be established to deter-
latter cases, the number of SA runs would be 16 if only mine the number of experiments necessary with a
the condition referring to γmax − γmin is considered. In heuristic such as the SA, which generates indepen-
this particular sample, it is worth noting both criteria dent random solutions. Fulfillment of the follow-
lead to a similar final decision about the algorithm stop ing two criteria is sufficient: (1) the difference be-
criterion. tween the cost of the best solution found to date and
Meccanica (2010) 45: 693–704 703

Fig. 9 Value of the minimum cost (Cmin ) and the location parameter (γ ) for different number of SA runs. Obtained using the bootstrap
technique

Table 3 γmin , γmax , Cmin and percentage differences for different numbers of SA runs. Obtained using the bootstrap technique

Number of tests γmin γmax (γmax − γmin )/γmin (%) Cmin (Cmin − γmin )/γmin (%)

9 3412.92 3443.85 0.91 3468.53 1.63

15 3430.28 3449.10 0.55 3468.53 1.12
25 3437.64 3440.05 0.07 3464.70 0.79
50 3437.64 3440.05 0.07 3462.28 0.72
100 3436.72 3437.64 0.03 3461.33 0.72
500 3435.40 3435.85 0.01 3460.00 0.72
1000 3435.40 3435.40 0.00 3460.00 0.72

the theoretical value estimated using the bootstrap Acknowledgements This study was funded by the Spanish
technique is less than a certain level—for exam- Ministry of Education (research project BIA2006-01444). The
authors are grateful to Debra Westall for her thorough revision
ple 1%—, and (2) the difference between the val-
of the manuscript.
ues of γmin and γmax obtained using the bootstrap
technique must be less than another preestablished
level, for example 0.5%.
References
• The methodology described is general and can be
easily adapted to other optimization problems, if the 1. Grierson DE (1994) Practical optimization of structural
results of the heuristic algorithm fit a Weibull distri- steel frameworks. In: Adeli H (ed) Advances in design op-
bution. Nevertheless, the threshold levels that define timization. Taylor & Francis, London
the stop criteria are problem dependent and should 2. Hernández S, Fontan A (2002) Practical applications of de-
sign optimization. WIT Press, Southampton
be chosen by the individual who solves the opti- 3. Dreo J, Petrowsky A, Siarry P, Taillard E, Chatterjee A
mization problem in terms of the accuracy required (2006) Metaheuristics for hard optimization. Methods and
and the time needed to perform one algorithm run. case studies. Springer, Berlin
704 Meccanica (2010) 45: 693–704

4. van Laarhoven PJM, Aarts EHL (1987) Simulated anneal- 20. Kicinger R, Arciszewski T, De Jong K (2005) Evolutionary
ing: theory and applications. Kluwer Academic, Dordrecht computation and structural design: a survey of the state of
5. Goldberg DE (1989) Genetic algorithms in search, op- the art. Comput Struct 83:1943–1978
timization, and machine learning. Addison-Wesley, New 21. Nieto F, Hernández S, Jurado JA (2009) Optimum de-
York sign of long-span suspension bridges considering aeroelas-
6. Adeli H, Sarma KC (2006) Cost optimization of structures. tic and kinematic constraints. Struct Multidiscip Optim
Fuzzy logic, genetic algorithms and parallel computing. 39(2):133–151
Wiley, Chichester 22. Marannano G, Mariotti GV (2008) Structural optimization
7. Goldberg DE, Samtani MP (1986) Engineering optimiza- and experimental analysis of composite material panels for
tion via genetic algorithms. In: ASCE proceedings of the naval use. Meccanica 43(2):251–262
ninth conference on electronic computation, New York, 23. Callegari M, Palpacelli MC (2008) Prototype design of a
pp 471–482 translating parallel robot. Meccanica 43(2):133–151
8. Coello CA, Christiansen AD, Santos F (1997) A simple ge- 24. Guadagni L (2008) Development of a collaborative opti-
netic algorithm for the design of reinforced concrete beams. mization tool for the sizing design of aerospace structures.
Eng Comput 13(4):185–196 Int J Simul Multidiscip Des Optim 2(3):187–192
9. Panigrahi SK, Chakraverty S, Mishra BK (2009) Vi- 25. McRoberts K (1971) A search model for evaluating combi-
bration based damage detection in a uniform strength natorially explosive problems. Oper Res 19:1331–1349
beam using genetic algorithm. Meccanica. doi:10.1007/ 26. Golden BL, Alt FB (1979) Interval estimation of a global
s11012-009-9207-1 optimum for large combinatorial problems. Nav Res Lo-
10. Bassir DH, Zapico JL, González MP, Alonso R (2007) gist Q 26(1):69–77
Identification of a spatial linear model based on earthquake- 27. Bettinger P, Boston K, Kim YH, Zhu J (2007) Landscape-
induced data and genetic algorithm with parallel selection. level optimization using tabu search and stand density-
Int J Simul Multidiscip Des Optim 1(1):39–48 related forest management prescriptions. Eur J Oper Res
11. Sid B, Domaszewski M, Peyraut F (2007) An adjacency 176:1265–1282
representation for structural topology optimization using 28. Payá-Zaforteza I (2007) Optimización heurística de pór-
genetic algorithm. Int J Simul Multidiscip Des Optim ticos de edificación de hormigón armado (Heuristic opti-
1(1):49–54 mization of reinforced concrete building frames). PhD the-
12. Balling RJ, Yao X (1997) Optimization of reinforced con- sis, Departamento de Ingeniería de la Construcción, Uni-
crete frames. ASCE J Struct Eng 123(2):193–202 versidad Politécnica de Valencia, Valencia (in Spanish)
13. Ceranic B, Fryer C, Bines RW (2001) An applica- 29. Ministerio de Fomento (1998) EHE code of structural con-
tion of simulated annealing to the optimum design of crete. Ministerio de Fomento, Madrid (in Spanish)
reinforced concrete retaining structures. Comput Struct 30. Ministerio de Fomento (1988) NBE AE-88. Code about
79(17):1569–1581 the actions to be considered in buildings. Ministerio de Fo-
14. González-Vidosa F, Yepes V, Alcalá J, Carrera M, mento, Madrid (in Spanish)
Perea C, Payá-Zaforteza I (2008) Optimization of re- 31. Kirkpatrick S, Gelatt CD, Vecchi MP (1983) Optimization
inforced concrete structures by simulated annealing. by simulated annealing. Science 220(4598):671–680
In: Tan CM (ed) Simulated annealing. I-Tech Edu- 32. Cerny V (1985) Thermodynamical approach to the travel-
cation and Publishing, Vienna, pp 307–320. Available ing salesman problem: an efficient simulation algorithm.
in http://intechweb.org/book.php?id=37. Accessed 6 Feb J Opt Theory Appl 45(1):41–51
2009 33. Medina JR (2001) Estimation of incident and reflected
15. Paya-Zaforteza I, Yepes V, Hospitaler A, González-Vidosa waves using simulated annealing. ASCE J Waterw Port
F (2009) CO2-optimization of reinforced concrete frames Coast Ocean Eng 127(4):213–221
by simulated annealing. Eng Struct 31(7):1501–1508 34. Weibull W (1951) A statistical distribution function of wide
16. Yepes V, Alcalá J, Perea C, González-Vidosa F (2008) applicability. ASME J Appl Mech Trans 18(3):293–297
A parametric study of optimum earth-retaining walls by 35. Fisher R, Tippett L (1928) Limiting forms of the frequency
simulated annealing. Eng Struct 30(3):821–830 distribution of the largest or smallest member of a sample.
17. Perea C, Alcalá J, Yepes V, González-Vidosa F, Hospitaler Proc Camb Philos Soc 24:180–190
A (2008) Design of reinforced concrete bridge frames by 36. Conover WJ (1971) Practical nonparametric statistics. Wi-
heuristic optimization. Adv Eng Softw 39(3):676–688 ley, New York
18. Martínez F, Yepes V, Hospitaler A, González-Vidosa F 37. Dannenbring DG (1977) Procedures for estimating optimal
(2007) Ant colony optimization of reinforced concrete solution values for large combinatorial problems. Manag
bridge piers of rectangular hollow section. In: Proceedings Sci 23(12):1273–1283
of the ninth international conference on the application of 38. Vasko FJ, Wilson JR (1984) An efficient heuristic for large
artificial intelligence to civil, structural and environmental set covering problems. Nav Res Logist Q 31:163–171
engineering (AICIVIL-COMP2007), St. Julians, Malta, pa- 39. ReliaSoft (2007) Weibull++7 user’s guide. Reliasoft, Tuc-
per 38 son
19. Payá I, Yepes V, González-Vidosa F, Hospitaler A (2008) 40. Efron B (1979) Bootstrap methods. Ann Stat 7:1–26
Multiobjective optimization of reinforced concrete build- 41. Zoubir AM, Boashash B (1998) The bootstrap and its ap-
ing frames by simulated annealing. Comput-Aided Civ In- plication in signal processing. IEEE Signal Process Mag
frastruct Eng 23:596–610 15(1):56–67