AP08_5.vp


1 Introduction
The inherent variability of a rock mass is difficult to model

and for this reason engineers very often have to ask the ques-
tion “What value should be used in analyses?” The answer to
such question requires a probabilistic approach in the evaluat-
ing the uncertainty in the input parameters in geologic me-
dia. In the recent times, specialists have encountered prob-
lems in the input parameters derived from uncertainty mod-
elling based on the Fuzzy set theory, Monte Carlo Simulation,
Latin Hypercube Sampling etc.

2 Fuzzy methods
In the design of underground structures, it is very difficult

to take into account the inherent variability of rock mass using
the current rock mass classifications. One of the main means
of improving rock mass classification is by accounting for
variations in the individual parameters using fuzzy mathe-
matics [1]. The fuzzy set was first introduced in 1965 by Lofti
A. Zadeh [2] as a mathematical way to represent linguistic
vagueness. In a classical set, an element either belongs to or
does not belong to a set. The concepts and definitions of
the fuzzy set theory are described in many publications.
Dubois and Pride [3] provide the following definition: “Fuzzy
set is a generalization of ordinary or classical set theory. It
consists of mathematical tools developed to model and pro-
cess incomplete and/or gradual information, ranging from
interval – valued numerical data to symbolic and linguistic
expressions”. Contrary to crisp (or ordinary) sets, fuzzy sets
have no sharp or precise boundaries. In crisp sets, element x
belongs to or does not belong to a set A and the membership
function (or degree of membership) �A(x) is unique. Fuzzy
sets assign the membership function �A(x) for each element x
as a range over the interval (0 to 1). This type of membership
function is characterised by a smooth transition from ”be-
longing to a set” (1) to ”not belonging to a set” (0) and gives
fuzzy sets flexibility in modelling based on linguistic expres-

sions of engineering practice (such as ”fairly rough surface”).
Membership functions can also be represented in analytical
methods.

Two types of variables are used in the fuzzy model: fuzzy
variables and fuzzy numbers. Fuzzy variables are defined di-
rectly as fuzzy sets based on a group of reference linguistic
terms. Generally, a fuzzy variable xj can acquire any value be-
tween the reference terms. For the linguistic terms the fuzzy
variables are fuzzy singletons. Fuzzy variables need not to
be associated with any numerical universe. This can be be-
cause their values are qualitative in nature, or because they
are treated as qualitative for the sake of convenience. Fuzzy
numbers are defined over a continuous domain – triangular
or trapezoidal membership functions (or numbers). Trape-
zoidal number is defined by the height and four distinct
elements of the fuzzy set interval. Fuzzy Logic provides a
number of functions for performing fuzzy arithmetic. The
fuzzy arithmetic functions are a little different from the rest of
Fuzzy Logic’s functions in that they operate on lists of fuzzy
numbers.

The application of fuzzy arithmetic in the rock mass classi-
fication is direct and generates a fuzzy number representing
the classification value. Fuzzy mathematics even introduces
the uncertainty in the evaluating of parameters in the rock
mass classification. Let us take, for example, the Index Q rock
mass classification. This classification was established in 1974
in Norway (by Barton, Lien, Lunde and Loset) and is based
on six parameters [4]:

Q
RQD

J
J
J

J
SRFn

r

a

w
� � � , (1)

where RQD � rock quality designation, Jn � joint set num-
ber, Ja � joint alteration number, Jw � joint water reduction
number, SRF � stress reduction factor.

By applying fuzzy logic to the equations for the Index Q
(Equation 1), we obtain results in the fuzzy classification value
with non-linear distributions. The convex nature of the flanks

©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 3

Acta Polytechnica Vol. 48  No. 5/2008

Statistical Analysis of Input Parameters
Impact on the Modelling of

Underground Structures
M. Hilar, J. Pruška

The behaviour of a geomechanical model and its final results are strongly affected by the input parameters. As the inherent variability of rock
mass is difficult to model, engineers are frequently forced to face the question “Which input values should be used for analyses?” The correct
answer to such a question requires a probabilistic approach, considering the uncertainty of site investigations and variation in the ground.
This paper describes the statistical analysis of input parameters for FEM calculations of traffic tunnels in the city of Prague. At the
beginning of the paper, the inaccuracy in the geotechnical modelling is discussed. In the following part the Fuzzy techniques are
summarized, including information about an application of the Fuzzy arithmetic on the shotcrete parameters. The next part of the paper is
focused on the stochastic simulation – Monte Carlo Simulation is briefly described, Latin Hypercubes method is described more in details. At
the end several practical examples are described: statistical analysis of the input parameters on the numerical modelling of the completed
Mrázovka tunnel (profile West Tunnel Tube km 5.160) and modelling of the constructed tunnel Špejchar – Pelc Tyrolka.

Keywords: Finite element method, input parameters, rock mass, Fuzzy, LHS sampling, Monte Carlo method, numerical model,
underground structures.



has the effect of increasing the possibility that the conditions
will be worse than a single computed index Q.

3 Monte Carlo simulation
Monte Carlo simulation is a well-known tool that is used to

analyze random phenomena. In the Monte Carlo simulation,
a random problem is transformed into several deterministic
problems that are much easier to solve – sample inputs are
used to generate sample outputs with statistical or probabilis-
tic information about the random output quantity. Monte
Carlo simulation is simple to use and therefore has found
much favour in geomechanics, particularly in stability analysis
of rock slopes [5].

The simplest sampling scheme of a Monte Carlo simula-
tion approach is to use a pseudorandom number generator to
select random numbers between 0 and 1 and use them to
generate values for each variable which is an input to the
calculation. However, this simple (and best-known) random
sampling scheme requires many samples for good accuracy
and repeatability – in practice, generating a probability distri-
bution of the safety factor of a rock slopes requires a minimum
of 200 up to 2000 selections (depending on the desired
accuracy).

The simulation output (random variable which depends
upon random input variables, fields and processes) may be
presented in several ways.

One way is to define the probability that a safety factor F is
less than a prescribed value F0:

P F F
n
N

( )� �0 , (2)

where n � number of trials in which F<F0 , N � number of
selections.

Another approach [6] is to plot a cumulative distribution
of F, which can be used to determine the probability of F
being less than a given value.

4 Latin hypercube sampling (LHS)
Other sampling methods have been developed to reduce

the number of samples required for good accuracy in the
Monte Carlo simulation. One of the best is Latin hypercube
sampling [7]. Latin hypercube sampling preserves marginal
probability distributions for each simulated variable (Fig. 1).
To fulfilled this aim, Latin hypercube sampling constructs a
highly depend joint probability density function for the ran-
dom variables in the problem, which allows good accuracy in
the response parameter even when using only a small number
of samples.

In Latin hypercube sampling, the bounded interval 0 1,
is uniformly divided into N non-overlaying intervals – Fig. 2
(with the same probability of accuracy). One value per interval
is generated afterwards. This can be performed by initially
generating N random numbers within range 0 1, . These
values are linearly transformed to the random numbers in the
non-overlaying intervals for each random variable:

4 ©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/

Acta Polytechnica Vol. 48  No. 5/2008

0.5�1 �0.5 0 1

0

F
t(
)

Range of the input parameters ´qi

0.5

0.6

0.8

1

0.4

0.2

�1

Fig. 2: Dividing the distribution function into intervals

0.45

0

0.4

0.35

0.25

0.15

0.05

0.3

0.2

0.1

Standard variable t

f
t(
)

�3 �2 �1 0 1 2 3

Fig. 1: Probability density function – Gaussian distribution



u
u
N

i
Ni � �
�( )1 , (3)

where i n�1 2, , ,� (n number of values), u � a random
number (u � 0 1, ), ui � random numbers in the ith interval,
N � number of the non-overlaying intervals.

From the above equation, it follows that we generate only
one value for each random variable (this value is randomly
selected within each of N intervals):

( )i
N

u
Ni

�
� �

1 1 , (4)

where i n�1 2, , ,� (n number of values), ui � a random num-
bers in the ith interval, N � number of the non-overlaying.

The n values obtained for the first random variable x1 are
paired in a random manner with the values of the second
random variable. These n pairs of (x1, x2) are combined in a
random manner with n values of the third random variable.
This combination results in the n triads (x1, x2, x3) which are
combined with n values of the 4th random variable and so on
until (K�2) triads are formed (K is a number of random vari-
ables). Thus we can assemble an N K� matrix. It can be seen
that certain statistical correlations between columns of the
matrix may have a significant influence on the results of simu-
lation [8].

5 Construction of tunnels in Prague
The north-western sector of the City Circular Road in

Prague contains three major road tunnels of large cross
section area. While the Strahov and Mrázovka tunnels are
already open for traffic, the Špejchar – Pelc Tyrolka tunnel
is still in the planning phase. The very difficult geological
conditions in the Prague Ordovician make it necessary to
use NATM mehods. Stabilization of the deformations was
achieved only after closure of the whole primary lining, but
the following interim supporting technical measures were
performed: anchoring, widening of the top heading legs,
micropile support under the top heading legs, reinforcing
grouting performed in advance from the exploratory gallery,

and closing the top heading by a temporary invert. Due to
increased deformations occurring at the excavation with a
horizontally divided face, a sequence with a vertical pattern
was used in the further course of excavating the Mrázovka
tunnel.

6 Modelling of the Mrázovka tunnel
The outputs of the Mrázovka tunnel modelling were veri-

fied by the Latin Hypercubes method [9]. The numerical
analysis was carried out by means of the PLAXIS program
system. The 3D behaviour of the excavation face area, and
correct description of the influence on the deformations and
the state of the massif were simulated by the common proce-
dure of loading the excavation and lining using the so-called
�-method (Fig. 3). The modelled area of the profile was
approximately 200 m wide and 110 m high and was divided
into eight basic sub-areas according to the types of rock en-
countered (Fig. 4). The rock mass behaviour was approxi-
mated by means of the Mohr-Coulomb model. The input
geotechnical parameters of the rock mass (Edef, �, c, �, �) were
determined on the basis of the engineering-geological inves-
tigation results – Tables 1, 2 and 3 [10]. A comparison of
theoretically determined deformations with the values ob-
tained by monitoring [11] was used to verify the applicability
of the mathematical model. The results of the statistical study
of the West Tunnel Tube (WTT) in profile km 5.160 are shown
in Table 4. They show that the probability of final settlements
being between 71 mm and 213 mm is 95 %. The range of
settlement without including deformations caused by the
excavation of the pilot adit is 65 mm to 198 mm. The
predicted range was consistent with the measured value of
194 mm (value does not include effect of pilot adit).

7 Modelling of the Špejchar –
Pelc-Tyrolka tunnel
The Špejchar – Pelc-Tyrolka project is 4.320 m long with

the length of the tunnelled section being 3.438 m [12]. As well

©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 5

Acta Polytechnica Vol. 48  No. 5/2008

Fig. 3: Numerical modelling stages



as the tunnel, the project includes underground garages in
Letná, four underground Technology Centres and the Trója
Bridge. As designed, the excavation will be carried out using
the New Austrian Tunnelling Method (NATM). Due to the
predicted conditions, a vertical excavation sequence will be
adopted in the three-lane tunnel excavation. In the two-lane
tunnels, a horizontal excavation sequence is expected.

The position of the exploratory drift has been selected to
coincide with the complicated sections of the tunnel to pro-
vide information to predicted whether the vertical sequence

should be applied to the top heading only, or to the entire
cross section. Mechanical rock breaking is expected, but in
combination with the drill-and-blast when passing through
the Quartzites. A transition zone between the Quartzites
exists at the foot of the slope falling from the Letná Plain. The
gradient parameters mean that the tunnels run are in the
vicinity of fully saturated Quaternary sediments. A multi-
criteria analysis resulted in the selection of pre-excavation

6 ©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/

Acta Polytechnica Vol. 48  No. 5/2008

Number of calculation a b In �1 In � 2 In � 3 In � 4 In � 5

Value qij’ �1 1 �0.64 �0.26 0.00 0.26 0.64

Ground 1 25 27 25.4 25.6 26.0 26.4 26.6

Ground 2 27 29 27.4 27.6 28.0 28.4 28.6

Ground 3 22 25 22.5 22.9 23.5 24.1 24.5

Ground 4 25 30 25.9 26.5 27.5 28.6 29.1

Ground 5 30 35 30.9 31.5 32.5 33.6 34.1

Ground 6 35 38 35.5 35.9 36.5 37.1 37.5

Ground 7 35 38 35.5 35.9 36.5 37.1 37.5

Ground 8 25 28 25.5 25.9 26.5 27.1 27.5

Table 1: Input � values for various geotechnical layers based on random intervals

Calculation q E1 � q2 � � q3 � � q c4 � q5 � �

1 4 2 2 5 3

2 1 3 5 2 1

3 3 5 3 4 5

4 5 1 4 3 4

5 2 4 1 1 2

Table 2: Random permutations of the input parameters qi

Parameter qij Edef [MPa] � [-] � [kN/m3] c [kPa] � [°]

Interval number In 4 2 2 5 3

Ground 1 4.8 0.30 17.3 11.3 26.0

Ground 2 13.6 0.30 18.8 5.3 28.0

Ground 3 24.2 0.35 20.2 14.1 23.5

Ground 4 79.7 0.31 22.6 22.3 27.5

Ground 5 171.0 0.29 24.3 45.5 32.5

Ground 6 271.0 0.26 25.1 91.0 36.5

Ground 7 421.0 0.20 25.8 91.0 36.5

Ground 8 24.5 0.33 23.1 13.2 26.5

Table 3: Input parameters for the first run

Fig. 4: Geometry of generated mesh in WTT km 5.160



grouting. The outputs from the exploratory drift near the
transition zone (km 5.900 STT) were verified by Latin Hyper-
cube Sampling method [13]. The 2D and 3D numerical
analyses were carried out by means of the CESAR – LCPC
program system. The rock mass behaviour was approximated
by means of the Mohr-Coulomb model. The intervals of the
input parameters of the rock mass were determined on the
basis of the engineering-geological investigation results. The
results of statistical study of South Tunnel Tube (STT) at pro-
file km 5.900 showed that final settlements of the tunnel
lining will be between 13 mm and 347 mm, with a probability
of 95 %. The range of predicted surface settlements above
excavated tunnels is from 7 mm to 26 mm.

8 Conclusion
The results of analyses show the effect of the variation of

input parameters describing rock mass behaviour on the re-
sults of FEM and demonstrate that the differences from this
influence can be significant. The differences in the results flow
from variations in the geological conditiond. However, the
importance for determination of the final structure behaviour
of at least a basic study of the variation of the input para-
meters can be clearly seen. The Latin Hypercube Sampling
method is a procedure with advantages for the qualified
statistical evaluation of FE calculation. These methods make
significant time saving possible in common statistical
methods (Monte Carlo method, estimations of probability
moments etc.).

9 Acknowledgment
The authors would like to thank to Ministry of Education

of the Czech Republic for research project VZ 03 CEZ MSM
6840770003 “Developments of the Algorithms of Computa-

tional Simulations and their Application in Engineering” for
assistance with a research.

References
[1] Hudson, J. A.: Rock Engineering Systems: Theory and Prac-

tice. Chichester: Ellis Horwood, 1998.
[2] Zadeh, L. A.: Fuzzy Sets, Inform Control, Vol. 8 (1965),

p. 338–353.
[3] Dubois, D., Prade H.: Fundamentals of Fuzzy Sets, Norwell,

MA, Kluwer, 2000.
[4] Bieniawski, Z. T.: Enngineering Rock Mass Classification.

New York: John Wiley and Sons, 1989,
ISBN 0-47160172-1

[5] McMahon, B. K.: A Statistical Method for the Design of
Rock Slopes. Proceedings 1st Australia – New Zealand Con-
ference on Geomechanics, Vol. 1, Melbourne, Australia,
1971, p. 314–321.

[6] Mahtab, M. A., Grasso, O.: Geomechanics Principles in
the Design of Tunnels and Caverns in Rock. Amsterdam:
Elsevier, 1992, ISBN 0-4440883088.

[7] Bažant, Z. P., Křístek, V.: The Effect of Shear Lag on
Creep Deformation. Journal of Structural Engineering,
Vol. 113 (1987), No. 3, p. 557–569.

[8] Et al.: Cesta k pravděpodobnostnímu posudku bezpeč-
nosti, provozuschopnosti a trvanlivosti konstrukcí. In:
Spolehlivost konstrukcí 2003 (in Czech), Ostrava. Dům
techniky, 2003, ISBN 80-02-01410-3 .

[9] Barták, J., Hilar, M., Pruška, J.: Statistical Analysis of In-
put Parameters Influence to the Tunnel Deformations
Modelling. In: Fourth International Conference on Compu-
tational Stochastic Mechanics [CD-ROM]. Houston: Rice
University, 2002, p. 280–292.

©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/ 7

Acta Polytechnica Vol. 48  No. 5/2008

Total settlement (mm) Pilot adit settlement (mm) Settlement without the pilot adit
effect (mm)

Calculation Surface Tunnel Surface Adit Surface Tunnel

1 79 107 15 8 64 99

2 153 197 29 14 124 183

3 92 125 18 10 74 115

4 85 111 15 8 70 103

5 134 171 25 12 109 159

X (average) 108.60 142.20 20.40 10.40 88.20 131.80

s (standard deviation) 29.41 35.61 5.64 2.33 23.80 33.31

X s p� �2 95 45( . %) 167.42 213.42 31.69 15.06 135.81 198.43

X s p� �2 95 45( . %) 49.78 70.98 9.11 5.74 40.59 65.17

X s p� �( . %)6827 138.01 177.81 26.04 12.73 112.00 165.11

X s p� �( . %)6827 79.19 106.59 14.76 8.07 64.40 98.49

VAR X (variance) 5.69 6.22 2.65 1.75 5.15 6.02

Table 4: Results of the Mrázovka tunnel calculations



[10] Hudek, J.: The Mrazovka tunnels – Measurement and Moni-
toring at Construction of the WTT – Presiometric Checking on
the Success of Saving Grouting – Profiles 003 to 009. Praha:
PÚDIS, 1999.

[11] Kolečkář, M., Zemánek I.: Monitoring of the Mrazovka
Tunnel. In: Volume of papers of the Int. Conf. Underground
Construction 2000, Praha, ITA/AITES, 2000, p. 427–433.

[12] Butovič, A.: Exploratory drift for the Blanka twin tube.
Tunel, Vol. 1 (2004), CTuK Praha, 2004, p. 13–18.

[13] Louženský, T.: Mathematical Modeling of Primary Lin-
ing Influence on Rock Mass Stress, Prague, CTU in
Prague, FCE, Diploma Thesis, 2006 (in Czech).

Doc. Ing. Matouš Hilar, M.SC., Ph.D.
phone/Fax: +420 241 443 411

D2 Consult Prague s.r.o.
Zelený pruh 95/97
147 00 Praha 4, Czech Republic

Doc. Dr. Ing. Jan Pruška
phone: +420 224 344 547
Fax: +420 224 344 556
e-mail: Pruska@fsv.cvut.cz

Department of Geotechnics

Czech Technical University in Prague
Faculty of Civil Engineering
Thákurova 7
166 29 Praha 6, Czech Republic

8 ©  Czech Technical University Publishing House http://ctn.cvut.cz/ap/

Acta Polytechnica Vol. 48  No. 5/2008