Discretization of Random Fields Representing Material Properties and Distributed Loads in FORM Analysis

Abstract

The reliability analysis of more complicated structures usually deals with the finite element method (FEM) models. The random fields (material properties and loads) have to be represented by random variables assigned to random field elements. The adequate distribution functions and covariance matrices should be determined for a chosen set of random variables. This procedure is called discretization of a random field. The chapter presents the discretization of random field for material properties with the help of the spatial averaging method of one-dimensional homogeneous random field and midpoint method of discretization of random field. The second part of the chapter deals with the discretization of random fields representing distributed loads. In particular, the discretization of distributed load imposed on a Bernoulli beam is presented in detail. Numerical example demonstrates very good agreement of the reliability indices computed with the help of stochastic finite element method (SFEM) and first-order reliability method (FORM) analyses with the results obtained from analytical formulae.

1. Introduction

In general, the safety of a structure is analyzed in the space ${\mathrm{\Omega }}_{\mathit{X}}=\left\{\mathit{X}\in R^n\right\}$ of basic random variables $\mathit{X}$ . For a given failure mode or serviceability requirement, represented by the limit state surface $g\left(\mathit{X}\right)=0$ , the space ${\mathrm{\Omega }}_{\mathit{X}}$ is divided into the safe subset, ${\mathrm{\Omega }}^S=\left\{\mathit{X}\in R^n;g\left(\mathit{X}\right)>0\right\}$ , and the failure subset, ${\mathrm{\Omega }}^F=\left\{\mathit{X}\in R^n;g\left(\mathit{X}\right)\le 0\right\}$ . If all random variables are continuous with the multivariate joint probability density function $f_{\mathit{X}}\left(\mathit{x}\right)$ , the failure probability is given by the integral

$$P_f=\underset{{\Omega }^F}{\int }{f}_x\left(\mathit{x}\right)d\mathit{x}$$ E1

The integral (Eq. (1)) can be evaluated exactly for a few cases with the most important one: the linear limit state surface and multidimensional normal distribution function of variables $\mathit{X}$ .

Development of reliability methods resulted in variety of powerful algorithms to estimate the probability of failure for complicated mechanical and statistical models of structures. The first-order reliability method (FORM) is the most popular approach applied in practice.

FORM algorithm starts with the nonlinear transformation. In general, non-normal random vector $\mathit{X}$ is transformed into a standard normal (Gaussian) vector $\mathit{Y}$ with zero mean and unit covariance matrix ${\mathit{C}}_{\mathit{YY}}=\mathit{I}$ . The limit state surface $g\left(\mathit{x}\right)=0$ is mapped into a limit state surface $G\left(\mathit{y}\right)=0$ . Next, the design point ${\mathit{y}}^{\ast }$ , that is, the point on the limit state surface with the minimum distance to the origin of the $\mathit{Y}$ space, is determined by solving the nonlinear optimization problem with a nonlinear constraint $G\left(\mathit{y}\right)=0$

$$\beta =\mathit{min}\sqrt{{\mathit{y}}^T\mathit{y}}\mathit{\text{for}}\mathit{y}\mathit{on}G\left(\mathit{y}\right)=0$$ E2

The hyperplane tangential to the limit state surface at the point ${\mathit{y}}^{\ast }$ is given by the formula

$$\beta -{\mathit{a}}^T\mathit{y}=0$$ E3

where $\mathit{a}$ is a unit outward normal vector to the hyperplane and $\mathrm{\beta }$ is the distance between the hyperplane and the origin (Figure 1). Since the random vector $\mathit{Y}=\mathit{Y}\left(\mathit{X}\right)$ has standard normal distribution, the first-order approximation of the failure probability is easily derived as follows

$$P_f\cong P\left[\beta -{\mathit{a}}^T\mathit{Y}\le 0\right]=\mathrm{\Phi }\left(-\beta \right)$$ E4

The nonlinear constrained optimization problem (Eq. (2)) can be solved with many standard procedures as well as algorithms developed especially for this purpose, for example, algorithm for the case of independent, non-normal random variables [1], algorithm for problems with incomplete probability information [2].

All such solvers are iterative: for the assumed value of design point ${\mathit{x}}_{\left(\mathrm{k}\right)}^{\ast }$ , the values of limit state function $g\left({\mathit{x}}_{\left(\mathrm{k}\right)}^{\ast }\right)$ and its gradient $\nabla g\left({\mathit{x}}_{\left(\mathrm{k}\right)}^{\ast }\right)$ are determined. Next, a new position of design point ${\mathit{x}}_{\left(\mathrm{k}+1\right)}^{\ast }$ is derived and the process continues until the convergence criteria are fulfilled. If the safety of mechanical problem is described by the limit state function with analytical form, then the gradient can be evaluated easily and one of the algorithms solving the optimization problem (Eq. (2)) can be applied directly. However, if the stochastic variability of material properties and loads is to be taken into account, SFEM approach must be applied.

Figure 1.

Definition of the FORM reliability index $\beta$ .

In general, the limit state function $g\left(\mathit{X}\right)=g\left(\mathit{R}\left(\mathit{X}\right),\mathit{S}\left(\mathit{X}\right)\right)$ can be represented in terms of two vectors: resistance variables $\mathit{R}$ and load effects $\mathit{S}$ . The elements of resistance variables vector $\mathit{R}$ (e.g., yield stress, allowable strain or allowable displacement), are prescribed to finite elements or nodes and can be treated as deterministic or random variables. In the latter case, the vector $\mathit{R}$ corresponds to the part of the vector of basic random variables $\mathit{X}$ . The vector of load effects $\mathit{S}$ (e.g., stresses, displacements and deformations) contains functions of basic random variables $\mathit{X}$ such as material properties, geometrical quantities or loads. The relation $\mathit{S}=\mathit{S}\left(\mathit{X}\right)$ is called the mechanical transformation. In most practical cases, the load effects $\mathit{S}$ have to be evaluated by using numerical algorithms, for example, FEM.

Two main problems are to be solved in order to apply FEM in FORM analysis:

• discretization of random fields of material properties and random fields of loads

• determination of the gradient of the limit state function $\nabla g\left({\mathit{x}}^{\ast }\right)$ , when the load effect is defined by means of the implicit mechanical transformation $\mathit{S}=\mathit{S}\left(\mathit{X}\right)$

The solution of the second problem is presented in many papers and books [3].

2. Probabilistic description of random fields

2.1. Basics definitions

The spatial probabilistic variability of physical quantities such as Young's modulus, thickness of a plate and intensity of a distributed load can be described by means of random fields, $w\left(\mathit{z}\right)$ , where $\mathit{z}$ is the vector of space coordinates. One-dimensional random fields can be defined for beams, bars and columns, two-dimensional random fields for plates or shells, and three-dimensional random fields for bodies.

For any specific location $\mathit{z}$ , random field $w\left(\mathit{z}\right)$ is a random variable with the cumulative distribution function

$$F_w\left(w\left(\mathit{z}\right)\right)=P\left[w\left(\mathit{z}\right)\le w\right]$$ E5

which is called the first-order distribution of the random field $w\left(\mathit{z}\right).$

The $m-\mathit{th}$ order distribution, that is, the joint cumulative distribution function of the random vector $\mathit{w}={\left[w\left({\mathit{z}}_1\right),\dots ,w\left({\mathit{z}}_m\right)\right]}^T$ , is defined as follows:

$$F_{w\dots w}\left(\mathit{w}\right)=P\left[w\left({\mathit{z}}_1\right)\le w_1,\dots ,w\left({\mathit{z}}_m\right)\le w_m\right]$$ E6

The first- and second-order probability density functions of the random field $w\left(\mathit{z}\right)$ are defined accordingly

$$\begin{array}{c}{f}_w\left(w\left(\mathit{z}\right)\right)=\frac{d}{\mathit{dw}}{F}_w\left(w\left(\mathit{z}\right)\right)\\ f_{\mathit{ww}}\left(w\left({\mathit{z}}_1\right),w\left({\mathit{z}}_2\right)\right)=\frac{{\partial }^2}{\partial w\left({\mathit{z}}_1\right)\partial w\left({\mathit{z}}_2\right)}{F}_{\mathit{ww}}\left(w\left({\mathit{z}}_1\right),w\left({\mathit{z}}_2\right)\right)\end{array}$$ E7

Following the well-known definition [4] with the help of the first-order probability density function (Eq. (6)) and the second-order probability density function (Eq. (7)), the second-order representation of the random field is defined by using the following functions: the mean value function ${\mu }_w\left(\mathit{z}\right)$ , the variance function ${\sigma }_w^2\left(\mathit{z}\right)$ , the covariance function, $C_w\left({\mathit{z}}_1,{\mathit{z}}_2\right)$ and the correlation function ${\rho }_w\left({\mathit{z}}_1,{\mathit{z}}_2\right)$ .

A random field $w\left(\mathit{z}\right)$ is called strict-sense homogeneous, if its statistics are invariant to the translation of the origin and in particular, the $n-\mathit{th}$ order density function has the property

$$f_{w\dots w}\left(w_1\left({\mathit{z}}_1\right),\dots ,w_n\left({\mathit{z}}_n\right)\right)=f_{w\dots w}\left(w_1\left({\mathit{z}}_1+\mathit{z}\right),\dots ,w_n\left({\mathit{z}}_{\mathit{n}}+\mathit{z}\right)\right)$$ E8

for any separation vector $\mathit{z}$ .

A random field $w\left(\mathit{z}\right)$ is called wide-sense homogeneous or second-order homogeneous if its mean value and variance are constant,

$${\mu }_w\left(\mathit{z}\right)={\mu }_w{\sigma }_w^2\left(\mathit{z}\right)={\sigma }_w^2$$ E9

and its covariance function as well as correlation function depends only on the separation vector $\mathit{z}$ ,

$$C_w\left(\mathit{z},,,0\right)=C_w\left({\mathit{z}}_1,{\mathit{z}}_1+\mathit{z}\right)=C_w\left(\mathit{z}\right)\mathit{for}\mathit{any}\mathit{z}$$ E10

A random field that is homogeneous in time is referred to as stationary process.

2.2. Ensemble average versus spatial average of random field

In order to estimate the statistical parameters of a random field, the sample (realization) must be collected in separate experiments. If the sample size is sufficiently large, the estimators of statistical parameters can be computed at each point of the random field domain. For example, at the location ${\mathit{z}}_1$ , the estimator of mean value and the estimator of variance are equal to ensemble averages

$${\widehat{\mu }}_w\left({\mathit{z}}_1\right)=\frac{1}{k}\sum _{i=1}^k{w}_{\left(i\right)}\left({\mathit{z}}_1\right)$$ E11

$${\widehat{\sigma }}_w^2\left({\mathit{z}}_1\right)=\frac{1}{k-1}\sum _{i=1}^k{\left[w_{\left(i\right)}\left({\mathit{z}}_1\right)-{\widehat{\mu }}_w\left({\mathit{z}}_1\right)\right]}^2$$ E12

where $k$ is the number of realizations (also called as the sample size) and $w_{\left(i\right)}\left(\mathit{z}\right)$ is the $i-\mathit{th}$ measurement of a random field.

Ensemble averages usually depend on the location vector. However, if the limit of the ensemble averages are invariant with respect to location

$${\widehat{\mu }}_w\left({\mathit{z}}_1\right)={\mathit{lim}}_{k\to \infty }\frac{1}{k}\sum _{i=1}^k{w}_{\left(i\right)}\left({\mathit{z}}_1\right)={\mathit{lim}}_{k\to \infty }\frac{1}{k}\sum _{i=1}^k{w}_{\left(i\right)}\left({\mathit{z}}_2\right)={\mu }_w$$ E13

then the random field can be considered as homogeneous, in strict- or wide-sense, which depends on the order of probability function invariant to location vector.

On the other hand, the spatial averages over the domain can be computed for every realization (measurement) of random field. For example, the average taken along with any single realization of a one-dimensional random field is equal to

$${\mu }_w\left(i,L\right)=\frac{1}{L}\underset{0}{\overset{L}{\int }}{w}_{\left(i\right)}\left(z\right)\mathit{dz}$$ E14

and it usually depends on the character of field and the length of averaging interval $L$ .

A homogeneous random field is called ergodic, if all statistical information can be obtained from one realization of the random field. This means that ensemble averages are invariant with respect to the location vector and the spatial averages are equal to the ensemble averages. Thus, in case of a homogeneous one-dimensional random field, the ensemble and spatial averages are equal in the limit

$$\underset{k\to \infty }{\mathit{lim}}\frac{1}{k}\sum _{i=1}^k{w}_{\left(i\right)}\left(\mathit{z}\right)=\underset{L\to \infty }{\mathit{lim}}\frac{1}{L}{\int }_0^L{w}_{\left(i\right)}\left(\mathit{z}\right)d\mathit{z}$$ E15

In general, it is usually difficult to prove that a random field is homogeneous, and it is even more difficult to prove that a random field is ergodic. Great number of samples over a sufficiently large domain should be collected. These conditions are rarely fulfilled. Thus, the homogeneity and ergodicity is usually assumed. Most of the concepts and methods developed in the reliability analysis are based on these assumptions.

2.3. One-dimensional homogeneous random field

A one-dimensional homogeneous random field is often used in the reliability analysis of linear elements such as beams, bars and frames. All the above-mentioned definitions are valid in this case.

The variance reduction function ${\gamma }_w\left(L\right)$ , which has been presented in detail in [5, 7], describes the correlation of the moving average $w_L\left(z\right)$ of one-dimensional homogeneous random field $w\left(z\right)$

$$w_L\left(z\right)=\frac{1}{L}{\int }_{z-L/2}^{z+L/2}w\left(l\right)\mathit{dl}$$ E16

where $L$ denotes the length of the averaging segment.

The mean value function and the variance function of the random field $w_L\left(z\right)$ are easy to determine

$$E\left[w_L\left(z\right)\right]=E\left[w\left(z\right)\right]={\mu }_w$$ E17

$$\mathit{Var}\left[w_L\left(z\right)\right]={\sigma }_L^2={\gamma }_w\left(L\right){\sigma }_w^2$$ E18

where ${\gamma }_w\left(L\right)$ is the variance reduction function, and ${\mu }_w$ and ${\sigma }_w^2$ are mean value and variance value of the one-dimensional homogeneous random field $w\left(z\right).$ The variance reduction function ${\gamma }_w\left(L\right)$ demonstrates how fast the point variance ${\sigma }_w^2$ is reduced under local averaging. This dimensionless function has the following properties:

$$\begin{array}{ccc}{\gamma }_w\left(L\right)={\gamma }_w\left(-L\right)\ge 0& & {\gamma }_w\left(0\right)=1\end{array}$$ E19

and is related to the correlation function ${\rho }_w\left(z\right)$ of the one-dimensional homogeneous random field $w\left(z\right)$ by the integral

$${\gamma }_w\left(L\right)=\frac{1}{L^2}\underset{0}{\overset{L}{\int }}\underset{0}{\overset{L}{\int }}{\rho }_w\left(l_1-l_2\right)d{l}_{1}d{l}_2=\frac{2}{L}\underset{0}{\overset{L}{\int }}\left[1-\frac{l}{L}\right]{\rho }_w\left(l\right)\mathit{dl}$$ E20

Another useful scalar measure of the correlation is the scale of fluctuation ${\theta }_w$ defined by the limit value of the variance reduction function

$${\theta }_w={\mathit{lim}}_{L\to \infty }L{\gamma }_w\left(L\right)$$ E21

It can be proved that the scale of fluctuation ${\theta }_w$ is related to the correlation function ${\rho }_w\left(z\right)$

$${\theta }_w=2\underset{0}{\overset{\infty }{\int }}{\rho }_w\left(l\right)\mathit{dl}$$ E22

The variance reduction function and the scale of fluctuation are especially useful in the discretization procedure of the homogeneous random field.

Table 1 presents four correlation models. It should be noticed that the rectangular and triangular models are not proper correlation functions for the homogeneous random field, since they do not fulfill the basic condition of weak-homogeneity. However, they are quite often assumed, mostly as visualization tools. Triangular model demonstrates the meaning of the scale of fluctuation in a simple way, that is, correlation between values of random field at points separated by greater distance than the scale of fluctuation is equal to zero. The rectangular model constitutes the upper limit for variance reduction functions. The simple form of the exponential correlation function makes analytical computation of many integrals possible. On the other hand, similarity between the squared exponential model and the triangular model allows the simple physical interpretation of the scale of fluctuation, that is, the correlation functions are equal to zero for a separation interval greater than the scale of fluctuation.

Model	Correlation function	Variance function
Rectangular	$\rho \left(z\right)=\left\{\begin{array}{cc}1& \left|z\right|\le \frac{\theta }{2}\\ 0& \left|z\right|>\frac{\theta }{2}\end{array}\right.$	$\gamma \left(L\right)=\left\{\begin{array}{cc}1& L\le \frac{\theta }{2}\\ \left(\frac{\theta }{L}\right)\left[1-\frac{\theta }{4L}\right]& L>\frac{\theta }{2}\end{array}\right.$
Triangular	$\rho \left(z\right)=\left\{\begin{array}{cc}1-\frac{\left|z\right|}{\theta }& \left|z\right|\le \theta \\ 0& \left|z\right|>\theta \end{array}\right.$	$\gamma \left(L\right)=\left\{\begin{array}{cc}1-\frac{L}{3\theta }& L\le \theta \\ \left(\frac{\theta }{L}\right)\left[1-\frac{\theta }{3L}\right]& L>\theta \end{array}\right.$
Exponential	$\rho \left(z\right)=\mathit{exp}\left(-2\frac{\left|z\right|}{\theta }\right)$	$\gamma \left(L\right)=\frac{1}{2}{\left(\frac{\theta }{L}\right)}^2\left[\frac{2L}{\theta }-1+\mathit{exp}\left(-2\frac{L}{\theta }\right)\right]$
Squared exponential	$\rho \left(z\right)=\mathit{exp}\left(-\pi {\left(\frac{z}{\theta }\right)}^2\right)$	$\gamma \left(L\right)=\frac{2\theta }{L}\left[\mathrm{\Phi }\left(\frac{L\sqrt{2\pi }}{\theta }\right)-0,5\right]+A$; $A=\frac{{\theta }^2}{\pi L^2}\left[\mathit{exp}\left(-\pi {\left(\frac{L}{\theta }\right)}^2\right)-1\right]$

Table 1.

Description of four correlation models.

Note: $\theta$ is the scale of fluctuation; $L$ is the separation (averaging) distance; and $\mathrm{\Phi }\left(z\right)$ is the Laplace function.

A special case of random field is the Gaussian random field, in which the random variables $w\left({\mathit{z}}_1\right),\dots ,w\left({\mathit{z}}_n\right)$ for any points ${\mathit{z}}_1,\dots ,{\mathit{z}}_n$ are jointly normal distributed. This random field is completely determined by two functions such as the mean value function and the covariance function. The n-th order probability density function has the joint normal density.

Figure 2 presents correlation functions and Figure 3 presents variance reduction functions for the correlation models described in Table 1.

Figure 2.

Four correlation functions for correlation models defined in Table 1.

Figure 3.

Four variance reduction functions for correlation models defined in Table 1.

3. Discretization of random fields representing material properties

A vast amount of papers deal with the problem how to develop the accurate and numerically efficient discretization methods for random fields of material properties.

The variability of random field is usually more accurately represented, if the number of random field elements or the number of series components is increased. However, greater number of random variables leads to longer computation time for realistic problems. Therefore, it has been an important issue to find out the optimal size of random field elements with respect to the scale of fluctuation, that is the scalar correlation measure. The accuracy of different methods is discussed by Zeldin and Spanos [9].

The reliability analysis of more complicated structures usually deals with FEM models. The random fields (material properties and loads) have to be represented by random variables assigned to random field elements. The adequate distribution functions and covariance matrices should be determined for a chosen set of random variables. This procedure is called discretization of a random field.

Two groups of methods for discretization of material random fields can be distinguished:

1. Random field elements

The value of a material property for any finite element is represented by a single random variable, constant within a random field element. Mean value, standard deviation and covariance as well as distribution function can be assigned to those random variables according to different procedures:

• the spatial averaging method [10]

• the midpoint method [8, 9, 11]

• the interpolation method [13]

2. Random series

The random field is described in terms of series of deterministic functions and random coefficients. Two examples of this approach are as follows:

• series composed of deterministic shape functions and random variables [16]

• the Karhunen-Loeve orthogonal expansion [15, 17, 20]

Two discretization methods, namely the spatial averaging method and the midpoint method, are presented in detail.

3.1. The spatial averaging method of one-dimensional homogeneous random field

The spatial averaging method has been developed by Vanmarcke [5]. We consider a one-dimensional homogeneous random field $w\left(z\right)$ that represents the spatial random variability of a material property, for example, modulus of elasticity along beam. In general, the domain of the random field can be divided into finite elements of lengths $L_i$ . The material property within the $i-\mathit{th}$ element is represented by a random variable which is assumed to be equal to the spatial average over the $i-\mathit{th}$ finite element

$$W_i=\frac{1}{L_i}\underset{0}{\overset{L_i}{\int }}w\left(z\right)\mathit{dz}$$ E23

The mean value of a random variable $W_i$ is equal to the mean value of the random field $w\left(z\right)$

$$E\left[W_i\right]=\frac{1}{L_i}\underset{0}{\overset{L_i}{\int }}E\left[w\left(z\right)\right]\mathit{dz}={\mu }_w$$ E24

and the variance of random variable $W_i$ is expressed in terms of the variance function ${\gamma }_w\left(L\right)$ of the random field $w\left(z\right)$

$$\mathit{Var}\left[W_i\right]={\gamma }_w\left(L_i\right){\sigma }_w^2$$ E25

The formula for the covariance between two random variables $W_i$ and $W_j$ related to the i-th and j-th random elements is more complicated

$$\mathit{Cov}\left[W_i,W_j\right]=\frac{{\sigma }_w^2}{2L_i{L}_j}\sum _{k=0}^3{\left(-1\right)}^k{L}_k^2\gamma \left(L_k\right)$$ E26

where the distances $L_k$ are defined in Figure 4.

Figure 4.

Definition of intervals used in the calculation of covariance between the spatial averages related to two random field elements.

Eqs. (25 and 26) can be generalized for a random field defined in two- or three-dimensional spaces [5]. Eq. (26) depends on the variance reduction function ${\gamma }_w\left(L\right)$ , which expresses a relation between the variance of the spatial average and the size of the averaging interval $L$ .

Since full information about the variability of the random field is seldom available, Vanmarcke [5–6] suggested using in the practical analysis the approximation of the variance reduction function by its asymptotic form

$${\gamma }_w\left(L\right)=\left\{\begin{array}{cc}1& L\le {\theta }_w\\ \frac{{\theta }_w}{L}& L>{\theta }_w\end{array}\right.$$ E27

where ${\theta }_w$ is the scale of fluctuation.

FORM analysis demands the knowledge about distribution functions of basic random variables. The spatial averaging method results in the normal random variables for the Gaussian random field $w\left(z\right)$ , since the integration is a linear operation. However, for non-Gaussian random field, it is difficult to derive the distribution function of a random variable $W_i$ defined by Eq. (23).

Der Kiureghian [18] has suggested a heuristic model for the distribution of random variable $W_i$ , which is based on the concept of the weighting the random field by the shape function, which results in the weighted variance reduction function. Figure 5 shows that the weighted variance function has much smaller values than the original variance reduction function. If the averaging interval $L_i=L_{\mathit{\text{long}}}$ is assumed many times longer than the scale of fluctuation ${\theta }_w$ , the distribution of random variable $W_i$ tends to have the normal distribution according to the central limit theorem. Then, the variance reduction function can be approximated as follows:

$$\begin{array}{cc}{\gamma }_w\left(L_{\mathit{\text{long}}}\right)\approx \frac{{\theta }_w}{L_{\mathit{\text{long}}}}\approx \frac{1}{n}& \begin{array}{cc}\mathit{\text{for}}& L_{\mathit{\text{long}}}\gg {\theta }_w\end{array}\end{array}$$ E28

where the parameter $n\gg 1$ should be determined by calculations or judgment.

According to Figure 3, which shows the variance functions of four models, as well as taking into account Figure 5, $n\approx 10$ could be assumed, which means that for the element length which is 10 times longer than the scale of fluctuation, the normal distribution can be assigned to the random variable $W_i$ . For shorter elements, non-normal distribution should be considered. On the other hand, for very short element length, the distribution of random variable $W_i$ is close to the first-order distribution of random field. Taking into account both limits the approximate density distribution function of $W_i$ has been proposed by Der Kiureghian [18]

$$f_{W_i}\left(w_i\right)=\alpha \frac{\sqrt{n}}{{\sigma }_w}\phi \left(\frac{w_i-{\mu }_w}{{\sigma }_w/\sqrt{n}}\right)+\left(1-\alpha \right)f_{w\left(z\right)}\left(w_i\right)$$ E29

where $\phi \left(u\right)$ is the standard normal probability density function and $f_{w\left(z\right)}\left(w\right)$ is the first-order probability density function of random field $w\left(z\right)$ .

Figure 5.

The original and weighted variance reduction functions for the random field with exponential correlation function.

The weight parameter $0\le \alpha \le 1$ can be determined for the current length $L_i$ of random field element by requiring that the variance of random variable $W_i$ with the approximate density function ((Eq. (29)) is equal to the variance of the spatial average of random field over the length $L_i$ , that is, ${\gamma }_w\left(L_i\right){\sigma }_w^2$ . It can be shown that for the element of length $L_i$

$$\alpha =\frac{n}{n-1}\left(1-{\gamma }_w\left(L_i\right)\right)$$ E30

3.2. Midpoint method of discretization of random field

The midpoint method [3] corresponds to the interpolation method [13] with constant interpolation function and is suitable for discretization of non-Gaussian random fields. In general, a nonhomogeneous random field $w\left(z\right)$ is discretized with the help of a set of random variables defined as follows:

$$W_i=w\left(z_c^{\left(i\right)}\right)$$ E31

where $z_c^{\left(i\right)}$ determines the position of the centroid of the i-th element.

The mean value of random variable $W_i$ and the covariance between random variables $W_i\text{and}{W}_j$ are given below, correspondingly

$$\begin{array}{c}E\left[W_i\right]=E\left[w\left(z_c^{\left(i\right)}\right)\right]\mathit{Cov}\left[W_i,W_j\right]\end{array}=C_w\left(z_c^{\left(i\right)},z_c^{\left(j\right)}\right)$$ E32

where $C_w\left(z_1,z_2\right)$ is the covariance function of nonhomogeneous random field.

In the midpoint discretization method, the probability distribution of the random variable $W_i$ is equivalent to the first-order distribution of the random field $w\left(z\right)$

$$F_{W_i}\left(w_i\right)=P\left[W_i\le w_i\right]=P\left[w\left(z_c^{\left(i\right)}\right)\le w_i\right]=F_w\left(w_i\left(z_c^{\left(i\right)}\right)\right)$$ E33

Thus, in case of a homogeneous random field, the probability distribution function does not depend on the location of the centroid of the i-th element, $F_{W_i}\left(w_i\right)=F_w\left(w_i\right)$ .

3.3. Selection of the optimal size for random field mesh representing variability of a material property

Both discretization procedures described in the previous section are based on assumption that property (e.g., modulus of elasticity) is constant within a finite element, see [19].

In deterministic FEM, the variability of material properties is modeled by means of sufficient number of finite elements. The structural finite element size is chosen with respect to the gradient of the stress field.

In the same way, the variability of the random field is usually more accurately represented if many random field elements are used in the analysis. However, finer random field mesh increases the number of random variables, which leads to longer computation times for the reliability analysis. Therefore, it has been an important issue to find out the optimal size of the random field elements.

The scale of fluctuation of the random field has been shown to be a very important measure of the correlation, since it governs the optimal size of a random field mesh. Der Kiureghian and Ke [12] have shown that a sufficiently accurate value of the reliability index is obtained if the random field element size is between one-half and one-quarter of the scale of fluctuation for the midpoint method with the exponential correlation function. Hisada and Nakagiri [11] have presented similar results.

If random field elements shorter than one-quarter of the scale of fluctuation are chosen, then a singular correlation matrix can be obtained, indicating linear dependency of the random variables. Then, the nonhomogeneous linear transformation to the set of uncorrelated, normalized basic random variables must be proceeded by an extra transformation, which decreases the dimension of the random variable space. However, this extra transformation is not unique. Improper choice of the transformation can lead to numerical difficulties in the iteration procedure for determining the reliability index. Therefore, too small random field elements should be avoided. Liu and Liu [19] have derived a simple rule of thumb regarding the selection of an appropriate random field mesh: a coarse mesh should be assumed in an area where the gradient of the limit state function with respect to the random variable representing random field is small and a finer random field mesh in an area with large gradient.

However, in many case a random element size equal to the scale of fluctuation may be considered as adequate with respect to the reliability analysis accuracy.

4. Discretization of random fields representing distributed loads

Discretization of random field loads has not been a subject of many studies. The finite element modeling introduces the well-known procedure for representing distributed forces $\mathit{p}\left(\mathit{s}\right)$ by a set of equivalent nodal forces. The random field of distributed loads has to be discretized according to the structural finite element mesh. Thus, if the distributed loads are random, all equivalent nodal forces become random variables. This approach seems obvious and it can be applied directly to study the response variability of stochastic engineering problems [14]. However, if discretized random field load is a part of FORM calculations, then both the discretization procedure for random field loads as well as FORM/FEM analysis should be modified in order to get accurate results.

The general approach for the discretization procedure of the distributed body forces $\mathit{p}\left(\mathit{s}\right)$ is presented below [20]. A similar algorithm can be applied for other types of distributed loads (surface forces and initial stresses). For the purpose of FORM analysis, it is convenient to assume that the nodal forces for i-th element are random variables

$${\mathit{Q}}_{\mathit{p}}^{\left(i\right)}=\underset{V}{\int }{{\mathit{N}}^{\left(i\right)}}^T\mathit{p}\left(\mathit{s}\right)\mathit{dV}$$ E34

where ${\mathit{N}}^{\left(\mathrm{i}\right)}$ is the displacement interpolation matrix for the i-th element in the local coordinate system. Thus, the mean value vector of the load vector ${\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{i}\right)}$ is equal to

$$E\left[{\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{i}\right)}\right]=\underset{\mathrm{V}}{\int }{{\mathit{N}}^{\left(\mathrm{i}\right)}}^{\mathrm{T}}E\left[\mathit{p}\left(\mathit{s}\right)\right]\mathit{dV}$$ E35

The covariance matrix of the vectors of equivalent forces ${\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{i}\right)}\mathit{\text{and}}{\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{j}\right)}$ corresponding to the i-th and j-th finite elements has the form

$$\mathit{Cov}\left[{\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{i}\right)},{{\mathit{Q}}_{\mathit{p}}^{\left(\mathrm{j}\right)}}^{\mathrm{T}}\right]=\underset{{\mathrm{V}}^{\left(\mathrm{i}\right)}}{\int }\underset{{\mathrm{V}}^{\left(\mathrm{j}\right)}}{\int }{{\mathit{N}}^{\left(\mathrm{i}\right)}}^{\mathrm{T}}{\mathbf{C}}_{\mathit{pp}}{\mathit{N}}^{\left(\mathrm{j}\right)}d{V}^{\left(\mathrm{i}\right)}d{V}^{\left(\mathrm{j}\right)}$$ E36

where the matrix ${\mathit{C}}_{\mathit{pp}}$ contains the cross-covariance functions between different components of the vector random field.

In general, a load effect can be an internal force, stress or strain at any point of structure which does not coincide with a nodal point (or Gaussian integration point). Thus, the load effect $S^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ at a point $\mathrm{s}$ of i-th element can be represented as a sum of a general solution, $S_{\mathrm{u}}^{\left(\mathrm{i}\right)}$ (s) and a particular solution, $S_p^{\left(\mathrm{i}\right)}$ (s). The general solution $S_{\mathrm{u}}^{\left(\mathrm{i}\right)}$ (s) is the load effect as a function of geometry, material properties and equivalent nodal forces applied at all nodal points. Whereas the particular solution ${\mathrm{S}}_{\mathrm{f}}^{\left(\mathrm{i}\right)}$ (s) is the load effect at the point $\mathrm{s}$ of the i-th element due to the distributed body forces ${\mathrm{f}}^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ and reactions $Q_p^{\left(\mathrm{i}\right)}$ at the i-th element.

In the FORM analysis, the vector of basic random variables $\mathit{X}$ also contains nodal equivalent forces $Q_p^{\left(\mathrm{i}\right)}$ (where the parameter ${}^{"}{\mathrm{i}}^{"}$ runs over all finite elements). Thus, in the search for the most likely failure point, the current values of equivalent nodal forces $Q_p^{\left(\mathrm{i}\right)}$ have to be determined at each iteration step of FORM algorithm. Those equivalent nodal forces, valid at a specific step of FORM iteration, correspond to unknown functions of distributed body forces $p\left(\mathrm{s}\right)$ . In order to determine the particular solution $S_p^{\left(\mathrm{i}\right)}$ (s), the function of the distributed body forces $p^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ within the $i-\mathit{th}$ finite element must be known. One way to solve this problem is to assume that the distributed body forces can be approximated with the help of shape functions.

$${\mathit{p}}^{\left(i\right)}\left(\mathrm{s}\right)={{\mathit{N}}^{\left(i\right)}}^T\mathit{b}$$ E37

The matrix $\mathit{b}$ can be determined from the condition that the equivalent nodal force $Q_p^{\left(\mathrm{i}\right)}$ at the i-th finite element due to body forces $p^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ should be equal to the calculated equivalent nodal forces in the FORM algorithm. The vector ${\mathit{b}}_k$ , which is the k-th column of the matrix $\mathit{b}$ and corresponds to the k-th component of the vector $p^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ , is determined by solving the system of linear equations

$$\left(\underset{V}{\int }{\mathit{N}}_k^{{\left(i\right)}^T}{\mathit{N}}_{\mathit{k}}^{\left(i\right)}\mathit{dV}\right){\mathit{b}}_{\mathit{k}}={\mathit{Q}}_{\mathit{p}}^{\left(i\right)}$$ E38

where $N_{\mathrm{k}}^{\left(\mathrm{i}\right)}$ is the k-th row of the matrix $N^{\left(\mathrm{i}\right).}$

In this way, the function of the distributed body forces $p^{\left(\mathrm{i}\right)}\left(\mathrm{s}\right)$ as well as the particular solution $S_p^{\left(\mathrm{i}\right)}$ (s) is determined as functions of nodal equivalent forces $Q_p^{\left(\mathrm{i}\right)}$ .

For the distributed loads represented by the Gaussian random field, the components of vector $Q_p^{\left(\mathrm{i}\right)}$ , which are determined by means of a linear transformation (Eq. (34)), have the multidimensional normal distribution. For the non-Gaussian random field, the probability distribution function $F_{{\mathrm{Q}}_{\mathrm{p}}^{\left(\mathrm{i}\right)}}\left(\mathrm{r}\right)$ of the vector of nodal equivalent forces $Q_p^{\left(\mathrm{i}\right)}$ cannot be determined easily. The first possible choice is to assume the normal distribution on the basis of the central limit theorem. Another approximate solution has been developed [20] on the basis of the approach presented by Der Kiureghian [18].

5. Discretization of transverse distributed load for a Bernoulli-Euler beam

In case of Bernoulli-Euler beam, four shape functions are applicable

$$\begin{array}{cc}{\mathrm{N}}_1\left(\mathrm{s}\right)=1-3{\left(\frac{\mathrm{s}}{\mathrm{L}}\right)}^2+2{\left(\frac{\mathrm{s}}{\mathrm{L}}\right)}^3& {\mathrm{N}}_2\left(\mathrm{s}\right)=\mathrm{s}{\left(1-\frac{\mathrm{s}}{\mathrm{L}}\right)}^2\\ {\mathrm{N}}_3\left(\mathrm{s}\right)={\left(\frac{\mathrm{s}}{\mathrm{L}}\right)}^2\left(3-2\frac{\mathrm{s}}{\mathrm{L}}\right)& {\mathrm{N}}_4\left(\mathrm{s}\right)=\frac{{\mathrm{s}}^2}{\mathrm{L}}\left(\frac{\mathrm{s}}{\mathrm{L}}-1\right)\end{array}$$ E39

The nodal forces equivalent to the transverse distributed load $\mathrm{q}\left(\mathrm{s}\right)$ are defined by the integrals

$$\begin{array}{cc}{Q}_j=\underset{0}{\overset{\mathrm{L}}{\int }}{\mathrm{N}}_j\left(\mathrm{s}\right)\mathrm{q}\left(\mathrm{s}\right)\mathrm{ds}& \mathrm{j}=1,\dots ,4\end{array}$$ E40

The distributed load $\mathrm{q}\left(\mathrm{s}\right)$ is assumed to be a homogeneous random field, with constant mean ${\mathrm{\mu }}_{\mathrm{q}}$ value, constant variance ${\mathrm{\sigma }}_{\mathrm{q}}^2$ and the covariance function for the correlation function ${\mathrm{\rho }}_{\mathrm{q}}$

$$\mathrm{Cov}\left[\mathrm{q}\left(\mathrm{s}\right),\mathrm{q}\left(\mathrm{t}\right)\right]={\mathrm{\sigma }}_{\mathrm{q}}^2{\mathrm{\rho }}_{\mathrm{q}}\left(\left|\mathrm{t}-\mathrm{s}\right|\right)$$ E41

Thus, the mean value of the nodal force ${\mathrm{Q}}_{\mathrm{i}}$ is just equal to

$$\mathrm{E}\left[Q_i\right]={\mathrm{\mu }}_{\mathrm{q}}\underset{0}{\overset{\mathrm{L}}{\int }}{\mathrm{N}}_i\left(\mathrm{s}\right)\mathrm{ds}$$ E42

and the covariance between nodal force ${\mathrm{Q}}_{\mathrm{i}}^{\left(\mathrm{m}\right)}$ at the $m-\mathit{th}$ finite element and the nodal force ${\mathrm{Q}}_{\mathrm{i}}^{\left(\mathrm{n}\right)}$ at the $n-\mathit{th}$ finite element is defined by the double integral

$$\mathrm{Cov}\left[Q_i^{\left(\mathrm{m}\right)},Q_j^{\left(\mathrm{n}\right)}\right]={\mathrm{\sigma }}_{\mathrm{q}}^2\underset{{\mathrm{s}}_{1\mathrm{m}}}{\overset{{\mathrm{s}}_{2\mathrm{m}}}{\int }}\underset{{\mathrm{s}}_{1\mathrm{n}}}{\overset{{\mathrm{s}}_{2\mathrm{n}}}{\int }}{\mathrm{N}}_{\mathrm{i}}\left(\mathrm{s}\right){\mathrm{N}}_{\mathrm{j}}\left(\mathrm{t}\right){\mathrm{\rho }}_{\mathrm{q}}\left(\left|\mathrm{t}-\mathrm{s}\right|\right)\text{dsdt}$$ E43

where ${\mathrm{s}}_{1\mathrm{m}},{\mathrm{s}}_{2\mathrm{m}}$ and ${\mathrm{s}}_{1\mathrm{n}},{\mathrm{s}}_{2\mathrm{n}}$ are the coordinates of the two ends of the finite elements, defined in a common coordinate system.

If a beam is divided into N finite elements, then the random distributed load $\mathrm{q}\left(\mathrm{s}\right)$ is modeled by $\mathit{4}\mathit{N}$ random variables.

A typical limit state function for a beam can be defined at the $\mathit{m}-\mathit{th}$ finite element

$$g^{\left(m\right)}\left(\mathrm{s}\right)={\mathrm{R}}_{\mathrm{b}}^{\left(m\right)}\left(\mathrm{s}\right)-{\mathrm{M}}^{\left(m\right)}\left(\mathrm{s}\right)$$ E44

where ${\mathrm{R}}_{\mathrm{b}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ is the bending resistance at cross-section $\mathrm{s}$ (usually assumed as a basic random variable) and ${\mathrm{M}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ is the bending moment due to external loads at cross-section s, which is a random function depending on the other basic random variables, for example, random nodal equivalent forces ${\mathrm{Q}}_1^{\left(1\right)},{\mathrm{Q}}_2^{\left(1\right)},{\mathrm{Q}}_3^{\left(1\right)},{\mathrm{Q}}_4^{\left(1\right)},\dots {\mathrm{Q}}_1^{\left(\mathrm{N}\right)},{\mathrm{Q}}_2^{\left(\mathrm{N}\right)},{\mathrm{Q}}_3^{\left(\mathrm{N}\right)},{\mathrm{Q}}_4^{\left(\mathrm{N}\right)}$ .

The bending moment ${\mathrm{M}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ can be represented as a sum

$${\mathrm{M}}^{\left(m\right)}\left(\mathrm{s}\right)={\mathrm{M}}_{\mathrm{u}}^{\left(m\right)}\left(\mathrm{s}\right)+{\mathrm{M}}_{\mathrm{q}}^{\left(m\right)}\left(\mathrm{s}\right)$$ E45

The general solution ${\mathrm{M}}_{\mathrm{u}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ depends on the nodal displacements, which are the functions of all nodal forces ${\mathrm{Q}}_1^{\left(1\right)},{\mathrm{Q}}_2^{\left(1\right)},{\mathrm{Q}}_3^{\left(1\right)},{\mathrm{Q}}_4^{\left(1\right)},\dots {\mathrm{Q}}_1^{\left(\mathrm{N}\right)},{\mathrm{Q}}_2^{\left(\mathrm{N}\right)},{\mathrm{Q}}_3^{\left(\mathrm{N}\right)},{\mathrm{Q}}_4^{\left(\mathrm{N}\right)}$ imposed to FEM model. The particular solution ${\mathrm{M}}_{\mathrm{q}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ is the bending moment within the $\mathit{m}-\mathit{th}$ element due to the distributed load $\mathrm{q}\left(\mathrm{s}\right)$ and reactions $-{\mathrm{Q}}_1^{\left(\mathrm{m}\right)},-{\mathrm{Q}}_2^{\left(\mathrm{m}\right)},-{\mathrm{Q}}_3^{\left(\mathrm{m}\right)},-{\mathrm{Q}}_4^{\left(\mathrm{m}\right)}$ at the $\mathit{m}-\mathit{th}$ element,

$${\mathrm{M}}_{\mathrm{q}}^{\left(m\right)}\left(\mathrm{s}\right)=Q_2^{\left(m\right)}-Q_1^{\left(m\right)}\mathrm{s}-\underset{{\mathrm{s}}_{1\mathrm{m}}}{\overset{\mathrm{s}}{\int }}\mathrm{q}\left(\mathrm{t}\right)\left(\mathrm{s}-\mathrm{t}\right)\mathrm{dt}$$ E46

At the $\mathit{k}-\mathit{th}$ iteration step of FORM algorithm, the function of the distributed load is unknown for the corresponding nodal forces. Therefore, the function of distributed load imposed on the $\mathit{m}-\mathit{th}$ finite element is assumed as a linear combination of the shape functions

$$\mathrm{q}\left(\mathrm{s}\right)\cong {\mathrm{q}}^{\left(m\right)}\left(\mathrm{s}\right)=\sum _{\mathrm{i}=1}^4{b}_i{\mathrm{N}}_i\left(\mathrm{s}\right)$$ E47

and the unknown coefficients ${\mathrm{b}}_{\mathrm{i}}$ , which have to be determined by requiring that the equivalent nodal forces (Eq. (40)) for the function ${\mathrm{q}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ defined by relation (Eq. (47)) should be equal to the current equivalent nodal forces. The coefficients ${\mathrm{b}}_{\mathrm{i}}$ as well as the distributed load ${\mathrm{q}}^{\left(\mathrm{m}\right)}$ can be obtained as functions of the current equivalent nodal forces ${\mathrm{Q}}_1^{\left(\mathrm{m}\right)},{\mathrm{Q}}_2^{\left(\mathrm{m}\right)},{\mathrm{Q}}_3^{\left(\mathrm{m}\right)},{\mathrm{Q}}_4^{\left(\mathrm{m}\right)}$ at the $\mathit{m}-\mathit{th}$ finite element, by solving the system of linear equations

$$\begin{array}{cc}\sum _{i=1}^4\left[\underset{s_{1m}}{\overset{s_{2m}}{\int }}{\mathrm{N}}_i\left(s\right){\mathrm{N}}_j\left(s\right)\mathit{ds}\right]\mathrm{b}=Q_j^{\left(m\right)}& j=1,\dots ,4\end{array}$$ E48

Finally, the bending moment ${\mathrm{M}}_{\mathrm{q}}^{\left(\mathrm{m}\right)}\left(\mathrm{s}\right)$ is determined as a function of equivalent nodal forces

$${\mathrm{M}}_{\mathrm{q}}^{\left(m\right)}\left(\mathrm{s}\right)={\mathrm{m}}_1\left(\mathrm{s}\right)Q_1^{\left(m\right)}+{\mathrm{m}}_2\left(\mathrm{s}\right)Q_2^{\left(m\right)}+{\mathrm{m}}_3\left(\mathrm{s}\right)Q_3^{\left(m\right)}+{\mathrm{m}}_4\left(\mathrm{s}\right)Q_4^{\left(m\right)}$$ E49

where

$$\begin{array}{cc}{m}_1\left(s\right)& =-s\left(1-8\frac{s}{L}+20{\left(\frac{s}{L}\right)}^2-20{\left(\frac{s}{L}\right)}^3+7{\left(\frac{s}{L}\right)}^4\right)\\ m_2\left(s\right)& =1-60{\left(\frac{s}{L}\right)}^2+200{\left(\frac{s}{L}\right)}^3-225{\left(\frac{s}{L}\right)}^4+84{\left(\frac{s}{L}\right)}^5\\ m_3\left(s\right)& =-s\left(2\frac{s}{L}-10{\left(\frac{s}{L}\right)}^2+15{\left(\frac{s}{L}\right)}^3-7{\left(\frac{s}{L}\right)}^4\right)\\ m_4\left(s\right)& =-30{\left(\frac{s}{L}\right)}^2+140{\left(\frac{s}{L}\right)}^3-195{\left(\frac{s}{L}\right)}^4+84{\left(\frac{s}{L}\right)}^5\end{array}$$ E50

In order to determine the reliability index ${\mathrm{\beta }}_{\text{FORM}}$ , the gradient of the limit state function (Eq. (44)) has to be calculated. The partial derivatives of function (Eq. (44)) with respect to the nodal forces are given below

$$\begin{array}{cc}\frac{\partial g^{\left(m\right)}}{\partial Q_i^{\left(n\right)}}=-\sum _{j=1}^4\frac{\partial M_u^{\left(m\right)}}{{\partial u}_j^{\left(m\right)}}\frac{\partial u_j^{\left(m\right)}}{\partial Q_i^{\left(n\right)}}& \begin{array}{cc}\mathit{\text{if}}& \begin{array}{cc}m\ne n& i=1,\dots ,4\end{array}\end{array}\end{array}$$ E51

$$\begin{array}{cc}\frac{\partial g^{\left(m\right)}}{\partial Q_i^{\left(m\right)}}=-\sum _{j=1}^4\frac{\partial M_u^{\left(m\right)}}{{\partial u}_j^{\left(m\right)}}\frac{\partial u_j^{\left(m\right)}}{\partial Q_i^{\left(m\right)}}-m_i\left(s\right)& \begin{array}{cc}& \begin{array}{cc}& i=1,\dots ,4\end{array}\end{array}\end{array}$$ E52

where ${\mathrm{u}}_{\mathrm{j}}^{\left(\mathrm{m}\right)}$ is the nodal displacements in the local coordinate system of the $\mathit{m}-\mathit{th}$ element and the derivatives $\frac{\partial {\mathrm{u}}_{\mathrm{j}}^{\left(\mathrm{m}\right)}}{\partial {\mathrm{Q}}_{\mathrm{i}}^{\left(\mathrm{n}\right)}}$ are computed with the help of the SFEM algorithm (Liu Der Kiureghian, 1991).

5.1. Example of discretization of random distributed load for a simply supported beam

The deterministic simply supported beam of length L is subjected to the transverse homogeneous random load $\mathrm{q}\left(\mathrm{s}\right)$ with mean value ${\mathrm{\mu }}_{\mathrm{q}}$ and variance ${\mathrm{\sigma }}_{\mathrm{q}}^2$ . Assuming the exponential correlation function with the scale of fluctuation ${\mathrm{\theta }}_{\mathrm{q}}$ ,

$$\rho \left(\tau \right)=\text{exp}\left(-2\frac{\left|\tau \right|}{{\theta }_q}\right)$$ E53

the mean value and the variance of the bending moment function can be derived analytically:

$$\mathrm{E}\begin{array}{c}\left[\mathrm{M}\left(s\right)\right]={\mu }_q\frac{s\left(L-s\right)}{2}\mathrm{Var}\left[\mathrm{M}\left(s\right)\right]={\sigma }_q^2\left[A\frac{{\theta }_q^4}{8}+B\frac{L{\theta }_q^3}{4}+C\frac{L^3{\theta }_q}{3}\right]\end{array}$$ E54

where

$$A=\left(1-\frac{s}{L}\right)\left[\frac{s}{L}\mathit{exp}\left(-2\frac{L}{{\theta }_q}\right)-\mathit{exp}\left(-2\frac{s}{{\theta }_q}\right)\right]+\frac{s}{L}\left[\frac{s}{L}-\mathit{exp}\left(-2\frac{L-s}{{\theta }_q}\right)-1\right]+1$$

$$\begin{array}{cc}\begin{array}{cc}B=-\left(1-\frac{s}{L}\right)\frac{s}{L}& \end{array}& \begin{array}{cc}& C={\left(1-\frac{s}{L}\right)}^2{\left(\frac{s}{L}\right)}^2\end{array}\end{array}$$

We consider the linear limit state function $\mathrm{g}\left(\mathrm{s}\right)={\mathrm{R}}_{\mathrm{b}}-\mathrm{M}\left(\mathrm{s}\right)$ where ${\mathrm{R}}_{\mathrm{b}}$ is the deterministic bending moment capacity, constant along the beam.

If the random field $\mathrm{q}\left(\mathrm{s}\right)$ is Gaussian, then the FORM reliability index is equivalent to the Cornell reliability index

$${\mathrm{\beta }}_{\mathit{\text{FORM}}}\left(\mathrm{s}\right)={\beta }_C\left(s\right)=\frac{R_b-\mathrm{E}\left[M\left(s\right)\right]}{\sqrt{\mathrm{Var}\left[M\left(s\right)\right]}}$$ E55

On the other hand, the FORM reliability indices have been computed for the finite element model of a simply supported beam. The distributed load random field $\mathrm{q}\left(\mathrm{s}\right)$ has been discretized according to the procedure described earlier.

The calculations have been carried out for the following data:

${\mathrm{\mu }}_{\mathrm{q}}=1000\mathrm{N}/\mathrm{m}$ , ${\mathrm{\sigma }}_{\mathrm{q}}=200\mathrm{N}/\mathrm{m}$ , ${\mathrm{R}}_{\mathrm{b}}=5000\mathrm{Nm}$ , $\mathrm{L}=6\mathrm{m}$ and 5 cross-sections: $\mathrm{s}=\left[1.02;1.5;2.04;2.52;3\right]\left(\mathrm{m}\right)$ . Three finite element sizes have been considered: ${\mathrm{L}}_{\mathrm{e}}=0.5;1;3\left(\mathrm{m}\right)$ . The scale of fluctuation has been assumed as: ${\mathrm{\theta }}_{\mathrm{q}}=0.5;1;4;100\left(\mathrm{m}\right)$

The results of the FORM analysis presented in Table 2 are in very good agreement with the reliability indices computed according to analytical formulae (Eq. (53–55)). Moreover, the reliability indices computed with the help of the described discretization procedure are insensitive to the scale of fluctuation and the finite element size.

Finite element length	Position of cross-section (m)
(m)	1.02	1.50	2.04	2.52	3.00
Scale of fluctuation = 0.5
Analytical	14.7350	7.3059	3.6040	2.1244	1.6821
0.5	14.7347	7.3055	3.6038	2.1243	1.6820
1.0	14.7345	7.3055	3.6038	2.1243	1.6820
3.0	14.7511	7.3166	3.6047	2.1247	1.6819
Scale of fluctuation = 1.0
Analytical	10.7242	5.3037	2.6112	1.5378	1.2172
0.5	10.7240	5.3036	2.6112	1.5377	1.2172
1.0	10.7239	5.3036	2.6112	1.5377	1.2172
3.0	10.7314	5.3092	2.6116	1.5378	1.2172
Scale of fluctuation = 4.0
Analytical	6.5774	3.2516	1.5997	0.9415	0.7451
0.5	6.5773	3.2517	1.5997	0.9415	0.7451
1.0	6.5773	3.2516	1.5997	0.9415	0.7451
3.0	6.5781	3.2523	1.5997	0.9415	0.7451
Scale of fluctuation = 100.0
Analytical	4.9153	2.4423	1.2062	0.7114	0.5633
0.5	4.9153	2.4423	1.2062	0.7114	0.5633
1.0	4.9153	2.4423	1.2062	0.7114	0.5633
3.0	4.9153	2.4423	1.2062	0.7114	0.5633

Table 2.

Reliability indices computed according to analytical formula and FORM algorithm.

6. Conclusions

The reliability analysis of more complicated structures usually deals with the FEM models. The random fields (material properties and loads) have to be represented by random variables assigned to random field elements. The adequate distribution functions and covariance matrices should be determined for a chosen set of random variables. This procedure is called discretization of a random field.

The chapter presents the discretization of random field for material properties with the help of the spatial averaging method of one-dimensional homogeneous random field, and midpoint method of discretization of random field.

The second part of the chapter deals with the discretization of random fields representing distributed loads. The discretization of distributed load imposed on a Bernoulli beam is presented in detail. An example shows that the presented procedure for discretizing random fields representing distributed loads is very efficient, that is, the reliability indices computed with the help of SFEM and FORM analysis are in very good agreement with the results of analytical calculations.