2. Safety assessment of the NPP risk

The International Atomic Energy Agency (IAEA) set up a programme [8] to give guidance to its member states on the many aspects of the safety of nuclear power reactors. The IAEA standards [8–10] and the US Nuclear Regulatory Commission (NRC) [4, 11] define the principal steps for the calculation of the risk of the NPP performance by simulations using Latin hypercube sampling (LHS) probabilistic method as follows:


(a magnitude 7.1), a 15-m tsunami disabled the power supply and cooling of three Fukushima

In view of the analysis results of the Fukushima accident [5, 6], the owner and operator of nuclear power plants (NPPs) in Slovakia re-evaluated the safety and reliability of structures and technology of all objects, according to the recommendations in "Stress tests" on all units in operation or under construction [7]. The main plant technological equipment (except for main circulation pumps) for all units was manufactured either in Czech Republic or in Slovakia [3]. In the process of design, construction and operation of plants, significant improvements in safety were made compared to the original project, including enhanced resistance to external hazards [1–3, 8]. Lately, in addition to previous safety improvements, many proposals have been implemented in the reconstruction process of NPP structures with reactor VVER 440/V213 in Slovakia for mitigation of severe accidents [3]. Fifty-one per cent of the overall

The NPP buildings with the reactor VVER 440/213 consist of the turbine hall, middle building, reactor building and bubble condenser [1]. The NPP containment is limited by the area of reactor building and bubble tower (Figure 1). In the case of the loss-of-coolant accident (LOCA), the

Source 2007 % 2008 % 2009 % Nuclear power plants 15,335 51.80 16,704 56.00 14,081 51.42 Thermal power plants 5421 18.30 5647 18.93 4768 17.41 Hydroelectric power plants 4485 15.10 4284 14.36 4662 17.02 Other power plants 2666 9.00 2674 8.96 2563 9.36 Import saldo 1725 5.80 521 1.75 1312 4.79 Consumption SR 29,632 100.00 29,830 100.00 27,386 100.00

Production SR 27,907 29,309 26,074

Table 1. Source ratio on annual electricity consumption in SR.

Figure 1. Section plane of the containment with reactor VVER440/213.

production of electricity in SR comes from the nuclear power plants (Table 1).

Daiichi reactors, hence causing a nuclear accident.

70 Recent Improvements of Power Plants Management and Technology

The final stage of the probabilistic safety analysis (PSA) is the compilation of the outputs of the first four steps into an expression of risk. The risk integration is shown in matrix formulation in Figure 2. The approximate numbers of plant damage states (PDSs), accident progression bins (APBs), and source-term groups (STGs), and the number of consequences are presented in the different reports [1–3, 12] and standards [4, 10–13].

## 3. Safety of the NPP structure resistance

In the case of the loss-of-coolant accident (LOCA) [3, 17, 18], the steam pressure expands from the reactor hall to the bubble condenser. Hence, the reactor hall and the bubble condenser are the critical structures of the NPP hermetic zone.

In the past, the calculation of the structural reliability for containment of the type of VVER 440/213 was carried out determining the probability density function for the ultimate pressure. A basis for the calculations consisted in results of linear and non-linear analysis, depending on whether the modelling considered the update of the properties.

The present work analyses the impact of combination of pressure load with thermal load that can arise in extreme situations related to severe accident progression. To achieve proper results, a detailed finite-element analysis of the concrete structure was carried out using ANSYS software and the programs CRACK [3, 14–18] were employed to solve this task. The basis for the probabilistic evaluation is reviewed, considering the uncertainties connected with loads and material properties. In the "Result" section, a comparison with linear evaluation is also mentioned.

## 3.1. Scenario of the hard accident

Most of the countries rate the safety of the NPP hermetic structures through the double-ended break guillotine test of the largest pipe in the reactor coolant system. During the design process of NPP structures with reactors VVER-440/230 type, the guillotine break of 2 500-mm pipes was considered as a beyond-design-basis accident (BDBA). The scenario of the hard accident is based on the assumption of the extreme situation where a loss of primary coolant accident (LOCA) is combined with the total loss of containment cooling system. During this situation, the hermetic zone cannot be available due to the external spraying of borated water. In practice, the contribution to heating from hydrogen explosion shall also be accounted.

In addition, the extremely climatic temperature (negative or positive) impacts the external slabs and walls of the hermetic zone. The temperature boundary conditions are defined to comply with the new revision of reference temperature provided by the Slovak hydrometeorological institute (SHMU) [3] and relevant Eurocodes [22] for a return period of 104 year. Three types of the scenarios were considered (Table 2).

### 3.2. Steel and concrete material properties under high-temperature impact

The recapitulation of the research works and the standard recommendations for the steel and concrete under high-temperature effect are summarized in US NRC report [23] and Eurocodes [22].


Table 2. The scenarios of the accidents in the hermetic zone.

The recommendations for the design of the structures are described in the US standards ACI [24], CEB-FIP Model Code [25] and Eurocodes [22]. The bonding phases of concrete are from the instable substance, which can be destructed at high temperature and their microstructures are changed.

The thermal conductivity λ<sup>c</sup> of normal weight concrete may be determined between the lower and upper limits given hereafter. The upper limit has been derived from tests of steel-concrete composite structural elements. The CEB-FIP Model Code [25] and Eurocodes [22] define the stress-strain relationship for concrete and steel materials dependent on temperature θ for heating rates between 2 and 50 K/min. In the case of the concrete, the stress–strain diagram is divided into two regions. The concrete strength σc,<sup>θ</sup> increases in the first region and decreases in the second region (Figure 3).

The stress-strain relations σc,<sup>θ</sup> ≈ εc,<sup>θ</sup> in region I are defined in the following form:

In the past, the calculation of the structural reliability for containment of the type of VVER 440/213 was carried out determining the probability density function for the ultimate pressure. A basis for the calculations consisted in results of linear and non-linear analysis,

The present work analyses the impact of combination of pressure load with thermal load that can arise in extreme situations related to severe accident progression. To achieve proper results, a detailed finite-element analysis of the concrete structure was carried out using ANSYS software and the programs CRACK [3, 14–18] were employed to solve this task. The basis for the probabilistic evaluation is reviewed, considering the uncertainties connected with loads and material properties. In the "Result" section, a comparison with linear evaluation is

Most of the countries rate the safety of the NPP hermetic structures through the double-ended break guillotine test of the largest pipe in the reactor coolant system. During the design process of NPP structures with reactors VVER-440/230 type, the guillotine break of 2 500-mm pipes was considered as a beyond-design-basis accident (BDBA). The scenario of the hard accident is based on the assumption of the extreme situation where a loss of primary coolant accident (LOCA) is combined with the total loss of containment cooling system. During this situation, the hermetic zone cannot be available due to the external spraying of borated water. In

practice, the contribution to heating from hydrogen explosion shall also be accounted.

3.2. Steel and concrete material properties under high-temperature impact

Type Duration Over-pressure [kPa] Extreme temperatures [C]

I. 1 day 150 127 42 28 II. 7 days 250 150 42 28 III. 1 year – 80–120 42 28

Three types of the scenarios were considered (Table 2).

Table 2. The scenarios of the accidents in the hermetic zone.

In addition, the extremely climatic temperature (negative or positive) impacts the external slabs and walls of the hermetic zone. The temperature boundary conditions are defined to comply with the new revision of reference temperature provided by the Slovak hydrometeorological institute (SHMU) [3] and relevant Eurocodes [22] for a return period of 104 year.

The recapitulation of the research works and the standard recommendations for the steel and concrete under high-temperature effect are summarized in US NRC report [23] and Eurocodes [22].

Interior Exterior

Max. Min.

depending on whether the modelling considered the update of the properties.

also mentioned.

3.1. Scenario of the hard accident

72 Recent Improvements of Power Plants Management and Technology

$$
\sigma\_{\mathbf{c},\theta} = f\_{\mathbf{c},\theta} \left[ \mathfrak{Z} \left( \frac{\varepsilon\_{\mathbf{c},\theta}}{\varepsilon\_{\mathbf{c}\mathbf{u},\theta}} \right) \Bigg/ \left\{ 2 + \left( \frac{\varepsilon\_{\mathbf{c},\theta}}{\varepsilon\_{\mathbf{c}\mathbf{u},\theta}} \right)^{3} \right\} \right], \quad f\_{\mathbf{c},\theta} = k\_{\mathbf{c},\theta} f\_{\mathbf{c}'} \tag{1}
$$

where the strain εcu,<sup>θ</sup> corresponds to stress fc,θ, the reduction factor can be chosen according to standard [22]. The reduction factors kc,<sup>θ</sup> (kc,<sup>θ</sup> = 0.925 for θ<sup>c</sup> = 150�C) for the stress–strain relationship are considered in accordance with the standard.

The stress-strain relationships for steel (Figure 4) are considered in accordance with Eurocode [22] on dependency of temperature level for heating rates between 2 and 50 K/min. In the case of steel, the stress-strain diagram is divided into four regions.

Figure 3. Stress-strain relationship of the concrete dependent on temperature [22].

Figure 4. Stress-strain relationship of the steel dependent on temperature [22].

The stress–strain relations σa,<sup>θ</sup> ≈ εa,<sup>θ</sup> are defined in the following form in region I:

$$
\sigma\_{a,\theta} = E\_{a,\theta} \varepsilon\_{a,\theta\prime} \quad E\_{a,\theta} = k\_{\text{E},\theta} E\_a \tag{2}
$$

where the reduction factor kE,<sup>θ</sup> can be chosen according to the values of [10].

In region II:

$$\begin{aligned} \sigma\_{a,\theta} &= (f\_{ay} - c) + \frac{b}{a} \sqrt{a^2 - (\varepsilon\_{ay,\theta} - \varepsilon\_{a,\theta})^2}, \sigma^2 = (\varepsilon\_{ay,\theta} - \varepsilon\_{ap,\theta})(\varepsilon\_{ay,\theta} - \varepsilon\_{ap,\theta} + c/\mathbb{E}\_{a,\theta}), \\ b^2 &= \mathbb{E}\_{\mathfrak{a},\theta}(\varepsilon\_{ay,\theta} - \varepsilon\_{ap,\theta})c + c^2, c = \frac{(f\_{ay,\theta} - f\_{ap,\theta})^2}{\mathbb{E}\_{\mathfrak{a},\theta}(\varepsilon\_{ay,\theta} - \varepsilon\_{ap,\theta}) - 2(f\_{ay,\theta} - f\_{ap,\theta})} \end{aligned} \tag{3}$$

and in region III:

$$
\sigma\_{a,\theta} = f\_{ay,\theta} \tag{4}
$$

A graphical display of the stress-strain relationships for steel grade S235 is presented in Figure 4 up to the maximum strain of εay,<sup>θ</sup> = 2%.

The strength and deformation properties of reinforcing steels under elevating temperatures may be obtained by the same mathematical model as that presented for structural steel S235. The reduction factors kE,θ(kE,<sup>θ</sup> = 0.95 for θ<sup>a</sup> = 150�C) for the stress-strain relationship are considered in accordance with the standard.

The material properties of the concrete structures in the numerical model were considered using the experimental tests statistically evaluated during the performance of the nuclear power plant [3, 28]. The material properties of the steel structures were not changed during plant performance [3].

#### 3.3. Nonlinear model of steel and reinforced concrete structures

The theory of large strain and rate-independent plasticity was proposed during the highoverpressure loading using the SHELL181-layered shell element from (Figure 5) the ANSYS library [26].

The vector of the displacement of the <sup>l</sup>th-shell layer <sup>f</sup>ul g¼ful x,ul y,ul zg <sup>T</sup> is approximated by the quadratic polynomial [26] in the form

$$\{\boldsymbol{u}^{l}\} = \begin{Bmatrix} \boldsymbol{u}^{l}\_{\boldsymbol{x}} \\ \boldsymbol{u}^{l}\_{\boldsymbol{y}} \\ \boldsymbol{u}^{l}\_{\boldsymbol{z}} \end{Bmatrix} = \sum\_{i=1}^{4} \mathbf{N}\_{i} \cdot \begin{Bmatrix} \boldsymbol{u}\_{\boldsymbol{x}i} \\ \boldsymbol{u}\_{\boldsymbol{y}i} \\ \boldsymbol{u}\_{\boldsymbol{z}i} \end{Bmatrix} + \sum\_{i=1}^{4} \mathbf{N}\_{i} \cdot \frac{\boldsymbol{\zeta}, t\_{i}}{2} \cdot \begin{bmatrix} \boldsymbol{a}\_{1,i} & \boldsymbol{b}\_{1,i} \\ \boldsymbol{a}\_{2,i} & \boldsymbol{b}\_{2,i} \\ \boldsymbol{a}\_{3,i} & \boldsymbol{b}\_{3,i} \end{bmatrix} \cdot \begin{Bmatrix} \boldsymbol{\theta}\_{\boldsymbol{x}i} \\ \boldsymbol{\theta}\_{\boldsymbol{y}i} \end{Bmatrix} \tag{5}$$

where Ni is the shape function for i-node of the 4-node shell element, uxi, uyi, and uzi are the motions of i-node, ζ is the thickness coordinate, ti is the thickness at i-node, {a} is the unit vector in x-direction, {b} is the unit vector in the plane of element and normal to {a}, θxi or θyi are the rotations of i-node about vector {a} or {b}.

Risk Assessment of NPP Safety in Case of Emergency Situations on Technology http://dx.doi.org/10.5772/intechopen.68772 75

Figure 5. SHELL181-layered element with smeared reinforcements [26].

In the case of the elastic state, the stress–strain relations for the lth-layer are defined in the form

$$\{\sigma^l\} = [D^l\_e] \{\varepsilon^l\} \tag{6}$$

where strain and stress vectors are as follows: <sup>f</sup>ε<sup>l</sup> g <sup>T</sup> ¼ fεx, <sup>ε</sup>y, <sup>γ</sup>xy, <sup>γ</sup>yz, <sup>γ</sup>zxg, <sup>f</sup>σ<sup>l</sup> g <sup>T</sup> ¼ fσx, <sup>σ</sup>y, <sup>τ</sup>xy, τyz, τzxg and the matrix of the material stiffness.

#### 3.3.1. Geometric nonlinearity

The stress–strain relations σa,<sup>θ</sup> ≈ εa,<sup>θ</sup> are defined in the following form in region I:

where the reduction factor kE,<sup>θ</sup> can be chosen according to the values of [10].

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi a<sup>2</sup> � ðεay,<sup>θ</sup> � εa,θÞ

2

3.3. Nonlinear model of steel and reinforced concrete structures

The vector of the displacement of the <sup>l</sup>th-shell layer <sup>f</sup>ul

9 = ; <sup>¼</sup> <sup>X</sup> 4

i¼1 Ni:

ul x ul y ul z

8 < :

rotations of i-node about vector {a} or {b}.

quadratic polynomial [26] in the form

ful g ¼ 2

, c <sup>¼</sup> <sup>ð</sup><sup>f</sup> ay,<sup>θ</sup> � <sup>f</sup> ap,θ<sup>Þ</sup>

A graphical display of the stress-strain relationships for steel grade S235 is presented in

The strength and deformation properties of reinforcing steels under elevating temperatures may be obtained by the same mathematical model as that presented for structural steel S235. The reduction factors kE,θ(kE,<sup>θ</sup> = 0.95 for θ<sup>a</sup> = 150�C) for the stress-strain relationship are

The material properties of the concrete structures in the numerical model were considered using the experimental tests statistically evaluated during the performance of the nuclear power plant [3, 28]. The material properties of the steel structures were not changed during

The theory of large strain and rate-independent plasticity was proposed during the highoverpressure loading using the SHELL181-layered shell element from (Figure 5) the ANSYS

> uxi uyi uzi

9 = ;

þ<sup>X</sup> 4

where Ni is the shape function for i-node of the 4-node shell element, uxi, uyi, and uzi are the motions of i-node, ζ is the thickness coordinate, ti is the thickness at i-node, {a} is the unit vector in x-direction, {b} is the unit vector in the plane of element and normal to {a}, θxi or θyi are the

i¼1 Ni: ζ:ti 2 :

8 < : g¼ful x,ul y,ul zg

> 2 4

a1,i b1,i a2,i b2,i a3,i b3,i 3 <sup>5</sup>: <sup>θ</sup>xi θyi � �

<sup>T</sup> is approximated by the

ð5Þ

In region II:

and in region III:

plant performance [3].

library [26].

<sup>σ</sup>a,<sup>θ</sup> ¼ ð<sup>f</sup> ay � <sup>c</sup>Þ þ <sup>b</sup>

<sup>b</sup><sup>2</sup> <sup>¼</sup> Ea,θðεay,<sup>θ</sup> � <sup>ε</sup>ap,θÞ<sup>c</sup> <sup>þ</sup> <sup>c</sup>

a

Figure 4 up to the maximum strain of εay,<sup>θ</sup> = 2%.

considered in accordance with the standard.

q

74 Recent Improvements of Power Plants Management and Technology

σa,<sup>θ</sup> ¼ Ea,θεa,θ, Ea,<sup>θ</sup> ¼ kE,θEa ð2Þ

2

Ea,θðεay,<sup>θ</sup> � εap,θÞ � 2ðf ay,<sup>θ</sup> � f ap,θÞ

, a<sup>2</sup> ¼ ðεay,<sup>θ</sup> � <sup>ε</sup>ap,θÞðεay,<sup>θ</sup> � <sup>ε</sup>ap,<sup>θ</sup> <sup>þ</sup> <sup>c</sup>=Ea,θÞ,

σa,<sup>θ</sup> ¼ f ay,<sup>θ</sup> ð4Þ

ð3Þ

If the rotations are large while the mechanical strains (those that cause stresses) are small, then it is possible to use a large rotation procedure. A large rotation analysis is performed in a static analysis in the ANSYS program [26].

From the following relations, the strain in the n-step of the solution can be computed:

$$\{\varepsilon\_n\} = [B\_o][T\_n]\{\mu\_n\} \tag{7}$$

Where {un} is the displacement vector, [Bo] is the original strain–displacement relationship and [Tn] is the orthogonal transformation relating the original element coordinates to the convected (or rotated) element coordinates.

#### 3.3.2. Material nonlinearity

The technology segments on board of the hermetic zone are made from the steel. The finite element model (FEM) of these segments is based on the HMH-yield criterion for the isotropic and homogenous material properties.

Consequently the stress–strain relations are obtained from the following relations

$$\{d\sigma\} = [D\_{el}](\{d\varepsilon\} - \{d\varepsilon^{pl}\}) = [D\_{el}]\left(\{d\varepsilon\} - d\lambda \left\{\frac{\partial Q}{\partial \sigma}\right\}\right) = [D\_{ep}]\{d\varepsilon\} \tag{8}$$

where [Dep] is elastic-plastic matrix in the form

$$\begin{aligned} \left[D\_{\rm cp}\right] = \left[D\_{\rm e}\right] - \frac{\left[D\_{\rm e}\right] \left\{\frac{\partial Q}{\partial \nu}\right\} \left\{\frac{\partial F}{\partial \nu}\right\}^T \left[D\_{\rm e}\right]}{A + \left\{\frac{\partial F}{\partial \nu}\right\}^T \left[D\_{\rm e}\right] \left\{\frac{\partial Q}{\partial \nu}\right\}} \tag{9} \end{aligned} \tag{9}$$

The hardening parameter A depends on the yield function and model of hardening (isotropic or kinematic). Huber-Mises-Hencky (HMH) yield is defined in the form

$$
\sigma\_{\text{eq}} = \sigma\_T(\kappa),
\tag{10}
$$

Where σeq is the equivalent stress in the point and σo(κ) is the yield stress that depends on the hardening.

In the case of kinematic hardening by Prager (vs Ziegler) and the ideal Bauschinger's effect, it is given as

$$A = \frac{2}{9E} \sigma\_r^2 H'\tag{11}$$

The hardening modulus H' for this material is defined in the form

$$H' = \frac{d\sigma\_{eq}}{d\varepsilon\_{eq}^p} = \frac{d\sigma\_T}{d\varepsilon\_{eq}^p} \tag{12}$$

When this criterion is used with the isotropic hardening option, the yield function is given by

$$F(\sigma) = \sqrt{\{\sigma\}^T [M] \{\sigma\}} - \sigma\_o(\varepsilon\_{\text{cp}}) = 0 \tag{13}$$

where σo(εep) is the reference yield stress, εep is the equivalent plastic strain and the matrix [M] is as follows:

$$[M] = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 2 & 0 & 0 \\ 0 & 0 & 0 & 0 & 2 & 0 \\ 0 & 0 & 0 & 0 & 0 & 2 \end{bmatrix} \tag{14}$$

#### 3.3.3. Nonlinear material model of the concrete structures

The presented constitutive model is a further extension of the smeared and oriented crack model, which was developed in [3]. A new concrete cracking-layered finite shell element was developed and incorporated into the ANSYS system [3] considering the experimental tests of real reinforced concrete plate and wall structures. The layered shell elements are proposed considering the nonlinear properties of the concrete and steel depending on temperature.

The concrete compressive stress fc, the concrete tensile stress ft and the shear modulus G are reduced during the crushing or cracking of the concrete. These effects are updated on the numerical model.

In this model, the stress-strain relation is defined (Figure 6) following CEB-FIP Model Code [25]:

• Loading-compression region εcu < εeq < 0

Dep � � <sup>¼</sup> De ½ ��

76 Recent Improvements of Power Plants Management and Technology

The hardening modulus H' for this material is defined in the form

FðσÞ ¼

½M� ¼

3.3.3. Nonlinear material model of the concrete structures

q

hardening.

is given as

is as follows:

numerical model.

or kinematic). Huber-Mises-Hencky (HMH) yield is defined in the form

De ½ � <sup>∂</sup><sup>Q</sup> ∂σ n o <sup>∂</sup><sup>F</sup> ∂σ � �<sup>T</sup> De ½ �

<sup>A</sup> <sup>þ</sup> <sup>∂</sup><sup>F</sup> ∂σ � �<sup>T</sup> De ½ � <sup>∂</sup><sup>Q</sup>

The hardening parameter A depends on the yield function and model of hardening (isotropic

Where σeq is the equivalent stress in the point and σo(κ) is the yield stress that depends on the

In the case of kinematic hardening by Prager (vs Ziegler) and the ideal Bauschinger's effect, it

<sup>A</sup> <sup>¼</sup> <sup>2</sup> <sup>9</sup><sup>E</sup> <sup>σ</sup><sup>2</sup>

<sup>H</sup><sup>0</sup> <sup>¼</sup> <sup>d</sup>σeq dε p eq

When this criterion is used with the isotropic hardening option, the yield function is given by

where σo(εep) is the reference yield stress, εep is the equivalent plastic strain and the matrix [M]

The presented constitutive model is a further extension of the smeared and oriented crack model, which was developed in [3]. A new concrete cracking-layered finite shell element was developed and incorporated into the ANSYS system [3] considering the experimental tests of real reinforced concrete plate and wall structures. The layered shell elements are proposed considering the nonlinear properties of the concrete and steel depending on temperature.

The concrete compressive stress fc, the concrete tensile stress ft and the shear modulus G are reduced during the crushing or cracking of the concrete. These effects are updated on the

100000 010000 001000 000200 000020 000002

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi <sup>f</sup>σg<sup>T</sup>½M�fσ<sup>g</sup>

<sup>¼</sup> <sup>d</sup>σ<sup>T</sup> dε p eq ∂σ

σeq ¼ σTðκÞ, ð10Þ

TH<sup>0</sup> ð11Þ

� σoðεepÞ ¼ 0 ð13Þ

ð12Þ

ð14Þ

n o <sup>ð</sup>9<sup>Þ</sup>

$$
\sigma\_{\varepsilon}^{\varepsilon^{\circ}} = f\_{\varepsilon}^{\varepsilon^{\circ}} \cdot \frac{k.\eta - \eta^{2}}{1 + (k - 2).\eta}, \eta = \frac{\varepsilon^{\epsilon\eta}}{\varepsilon\_{\varepsilon}} \,, \left(\varepsilon\_{\varepsilon} \doteq -0.0022, \varepsilon\_{\varepsilon u} \doteq -0.0035\right) \tag{15}
$$

• Softening-compression region εcm < εeq < εcu

$$
\sigma\_{\varepsilon}^{\text{cf}} = f\_{\varepsilon}^{\text{cf}} \left( 1 - \frac{\varepsilon^{\text{cq}} - \varepsilon\_{\text{c}}}{\varepsilon\_{\text{cm}} - \varepsilon\_{\text{cu}}} \right) \tag{16}
$$

• The tension region ε<sup>t</sup> < εeq < ε<sup>m</sup>

$$\sigma\_c^{\circ f} = f\_{t} \cdot \exp\left(-2.(\varepsilon^{a\eta} - \varepsilon\_t)/\varepsilon\_{tm}\right), (\varepsilon\_t \dot{=} 0.0001, \varepsilon\_{tm} \dot{=} 0.002) \tag{17}$$

The equivalent values of f eq <sup>t</sup> and f eq <sup>c</sup> were considered for the plane stress state. The relation between the one and bidimensional stresses state was considered in the plane of principal stresses (σc1, σc2) of each shell layer by Kupfer (see Figure 7) [3].

The shear concrete modulus G was defined for cracking concrete by Kolmar [23] in the form

$$\mathbf{G} = r\_{\mathcal{S}}.\mathbf{G}\_{\boldsymbol{\theta}\prime}\boldsymbol{r}\_{\mathcal{S}} = \frac{1}{c\_2} \ln\left(\frac{\varepsilon\_u}{c\_1}\right), c\_1 = 7 + 333(p - 0.005), c\_2 = 10 - 167(p - 0.005) \tag{18}$$

where Go is the initial shear modulus of concrete, ε<sup>u</sup> is the strain in the normal direction to crack, c<sup>1</sup> and c<sup>2</sup> are the constants dependent on the ratio of reinforcing and p is the ratio of reinforcing transformed to the plane of the crack (0 < p < 0.02).

Figure 6. The concrete stress-strain diagram.

Figure 7. Kupfer's plasticity function.

The strain-stress relationship in the Cartesian coordinates can be defined in dependency on the direction of the crack (in the direction of principal stress, vs strain)

$$\begin{bmatrix} \sigma\_{cr} \end{bmatrix} = \begin{bmatrix} D\_{cr} \end{bmatrix} \begin{Bmatrix} \varepsilon\_{cr} \end{Bmatrix} \text{ and } \begin{bmatrix} \sigma \end{bmatrix} = \begin{bmatrix} T\_{\sigma} \end{bmatrix}^{\mathrm{T}} \begin{bmatrix} D\_{cr} \end{bmatrix} \begin{bmatrix} T\_{\varepsilon} \end{bmatrix} \{\varepsilon\} \tag{19}$$

$$\mathbb{E}\left[D\_{cr}^{l}\right] = \left[T\_{c.o}^{l}\right]^{\mathrm{T}}[D\_{cr}^{l}][T\_{c.s}^{l}] + \sum\_{s=1}^{N\_{min}} [T\_{s}^{l}]^{\mathrm{T}}[D\_{s}^{l}][T\_{s}^{l}] \tag{20}$$

where [Tc.σ], [Tc.ε] and [Ts] are the transformation matrices for the concrete and the reinforcement separately, Nrein is the number of the reinforcements in the lth–layer (Figure 8).

The stress-strain relationship for the concrete lth-layer cracked in one direction is

$$\begin{Bmatrix} \sigma\_1\\ \sigma\_2\\ \tau\_{12}\\ \tau\_{13}\\ \tau\_{23} \end{Bmatrix}\_l = \begin{bmatrix} 0 & 0 & 0 & 0 & 0\\ 0 & E & 0 & 0 & 0\\ 0 & 0 & G\_{12}^{cr} & 0 & 0\\ 0 & 0 & 0 & G\_{13}^{cr} & 0\\ 0 & 0 & 0 & 0 & G\_{23}^{cr} \end{bmatrix} \begin{Bmatrix} \varepsilon\_1\\ \varepsilon\_2\\ \mathcal{V}\_{12}\\ \mathcal{V}\_{13}\\ \mathcal{V}\_{23} \end{Bmatrix}\_l \tag{21}$$

Figure 8. SHELL181-layered element with smeared reinforcements.

The strain-stress relationship in the Cartesian coordinates can be defined in dependency on the

where [Tc.σ], [Tc.ε] and [Ts] are the transformation matrices for the concrete and the reinforce-

00 0 0 0 0 E 000

00 0 0 Gcr

<sup>12</sup> 0 0

<sup>13</sup> 0

23

8 >>>><

>>>>:

ε1 ε2 γ<sup>12</sup> γ<sup>13</sup> γ<sup>23</sup> 9 >>>>=

>>>>; l

T

c:ε� þ<sup>X</sup> Nrein

s¼1 ½Tl s� <sup>T</sup>½D<sup>l</sup> s�½T<sup>l</sup>

½Dcr�½Tε�fεg ð19Þ

<sup>s</sup>� ð20Þ

ð21Þ

½σcr�¼½Dcr�fεcrg and ½σ�¼½Tσ�

ment separately, Nrein is the number of the reinforcements in the lth–layer (Figure 8).

0 0 Gcr

00 0 Gcr

The stress-strain relationship for the concrete lth-layer cracked in one direction is

direction of the crack (in the direction of principal stress, vs strain)

½Dl cr�¼½T<sup>l</sup> <sup>c</sup>:<sup>σ</sup>� T ½Dl cr�½T<sup>l</sup>

78 Recent Improvements of Power Plants Management and Technology

Figure 7. Kupfer's plasticity function.

σ1 σ2 τ<sup>12</sup> τ<sup>13</sup> τ<sup>23</sup> 9 >>>>=

>>>>; l ¼

8 >>>><

>>>>:

When the tensile stress in the 2-directions reaches the value f 0 t , the latter cracked plane perpendicular to the first one is assumed to be formed, and the stress–strain relationship becomes

$$\begin{Bmatrix} \sigma\_1\\ \sigma\_2\\ \tau\_{12}\\ \tau\_{13}\\ \tau\_{23} \end{Bmatrix}\_l = \begin{bmatrix} 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0\\ 0 & 0 & G\_{12}^{cr}/2 & 0 & 0\\ 0 & 0 & 0 & G\_{13}^{cr} & 0\\ 0 & 0 & 0 & 0 & G\_{23}^{cr} \end{bmatrix} \begin{Bmatrix} \varepsilon\_1\\ \varepsilon\_2\\ \gamma\_{12}\\ \gamma\_{13}\\ \gamma\_{23} \end{Bmatrix}\_l \tag{22}$$

where the shear moduli are reduced by parameters rg<sup>1</sup> and rg<sup>2</sup> by Kolmar [27] as follows:

$$G\_{12}^{cr} = G\_o.r\_{\mathfrak{g}1\prime}G\_{13}^{cr} = G\_o.r\_{\mathfrak{g}1\prime}G\_{23}^{cr} = G\_o r\_{\mathfrak{g}2\prime}.$$

The stress-strain relationship defined in the direction of the principal stresses must be transformed to the reference axes XY. The simplified smeared model of the concrete cracked is more convenient for finite element formulation than the singular discrete model.

The smeared calculation model is determined by the size of the finite element, hence its characteristic dimension Lc <sup>¼</sup> ffiffiffiffi A <sup>p</sup> , where <sup>A</sup> is the element area (vs integrated point area of the element). The assumption of constant failure energies Gf = const is proposed in the form

$$\mathcal{G}\_{\!\!\!f} = \bigcap\_{\text{O}}^{\text{w}} \sigma\_n(w) dw = A\_{\text{G}}.L\_{\text{c}} \quad \text{w}\_{\text{c}} = \mathcal{e}\_{\text{w}}.L\_{\text{c}} \tag{23}$$

where wc is the width of the failure, σ<sup>n</sup> is the stress in the concrete in the normal direction and AG is the area under the stress-strain diagram of concrete in tension. The descend line of concrete stress-strain diagram can be defined on dependency on the failure energies [25] by the modulus in the form

$$E\_{c,s} = E\_c/(1 - \lambda\_c), \qquad \lambda\_c = 2\mathbf{G}\_f E\_c/(\mathbf{L}\_c \sigma\_{\text{max}}^2) \tag{24}$$

where Ec is the initial concrete modulus elasticity, σmax is the maximal stress in the concrete tension. From the condition of the real solution of relation (18), it follows that the characteristic dimension of element must satisfy the following condition:

$$\mathcal{L}\_c \le \mathcal{Z} \mathcal{G}\_f \mathcal{E}\_c / \sigma\_{\text{max}}^2 \tag{25}$$

The characteristic dimension of the element is determined by the size of the failure energy of the element. The theory of a concrete failure was implied and applied to the two-dimensional (2D)-layered shell elements SHELL181 in the ANSYS element library [26]. The CEB-FIP Model Code [25] defines the failure energies G<sup>f</sup> [N/mm] depending on the concrete grades and the aggregate size d<sup>a</sup> as follows:

$$G\_{\circ} = (0.0469d\_a^2 - 0.5d\_a + 26)(f\_c/10)^{0.7} \tag{26}$$

The limit of damage at a point is controlled by the values of the so-called crushing or total failure function Fu. The modified Kupfer's condition for the lth-layer of section is as follows:

$$F\_u^l = F\_u^l(I\_{\epsilon1}; I\_{\epsilon2}; \varepsilon\_u) = 0, \quad F\_u^l = \sqrt{\beta(\Im\_{\ell 2} + aI\_{\epsilon1})} - \varepsilon\_u = 0 \,\,\,\,\tag{27}$$

where Iε<sup>1</sup>, Iε<sup>2</sup> are the strain invariants, and ε<sup>u</sup> is the ultimate total strain extrapolated from uniaxial test results (ε<sup>u</sup> = 0.002 in the tension domain, or ε<sup>u</sup> = 0.0035 in the compression domain), and α, β are the material parameters determined from Kupfer's experiment results (β = 1.355, α = 1.355εu).

The failure function of the whole section will be obtained by the integration of the failure function through the whole section in the form

$$F\_{\mathfrak{u}} = \frac{1}{t} \cdot \int\_{0}^{t} \left( I\_{\mathfrak{u}} (I\_{\varepsilon 1}; I\_{\varepsilon 2}; \varepsilon\_{\mathfrak{u}}) \right) d\mathfrak{z} = \frac{1}{t} \sum\_{l=1}^{N\_{\mathfrak{l} \mathfrak{g}}} F\_{\mathfrak{u}}^{l} (I\_{\varepsilon 1}; I\_{\varepsilon 2}; \varepsilon\_{\mathfrak{u}}) t\_{l} \tag{28}$$

where tl is the thickness of the lth-shell layer, t is the total shell thickness and Nlay is the number of layers.

The maximum strain ε<sup>s</sup> of the reinforcement steel in the tension area (maxðεsÞ ≤ εsm ¼ 0:01) and by the maximum concrete crack width wc (maxðwcÞ ≤ wcm ¼ 0:3 mm) determine the local collapse of reinforced concrete structure.

The program CRACK based on the presented nonlinear theory of the layered reinforced concrete shell was adopted in the software ANSYS [3]. These procedures were tested in comparison with the experimental results [3, 15, 17, 28].

Risk Assessment of NPP Safety in Case of Emergency Situations on Technology http://dx.doi.org/10.5772/intechopen.68772 81

Figure 9. Comparison of experimental and nonlinear numerical analysis of plates.

The reinforced concrete plates D4 (Figure 9) with the dimensions 3590/1190/120 mm were simply supported and loaded by pressure p on the area of plate. The plate D4 was reinforced by steel grid KARI Ø8 mm, a<sup>0</sup> 150 � 150 mm at the bottom. Material characteristics of plate D4 are the following: Concrete, E<sup>c</sup> = 30.9 Gpa, μ = 0.2, f<sup>c</sup> = �34.48 Mpa, f<sup>t</sup> = 4.5 Mpa and the Reinforcement, E<sup>s</sup> = 210.7 Gpa, μ = 0.3, f<sup>s</sup> = 550.3 Mpa.

## 4. Probabilistic assessment

concrete stress-strain diagram can be defined on dependency on the failure energies [25] by the

Ec,s <sup>¼</sup> Ec=ð<sup>1</sup> � <sup>λ</sup>cÞ, <sup>λ</sup><sup>c</sup> <sup>¼</sup> <sup>2</sup>Gf Ec=ðLc:σ<sup>2</sup>

dimension of element must satisfy the following condition:

80 Recent Improvements of Power Plants Management and Technology

Gf ¼ ð0:0469d<sup>2</sup>

<sup>u</sup>ðIε<sup>1</sup>; Iε<sup>2</sup>; <sup>ε</sup>uÞ ¼ <sup>0</sup>, Fl

where Ec is the initial concrete modulus elasticity, σmax is the maximal stress in the concrete tension. From the condition of the real solution of relation (18), it follows that the characteristic

Lc ≤ 2Gf Ec=σ<sup>2</sup>

The characteristic dimension of the element is determined by the size of the failure energy of the element. The theory of a concrete failure was implied and applied to the two-dimensional (2D)-layered shell elements SHELL181 in the ANSYS element library [26]. The CEB-FIP Model Code [25] defines the failure energies G<sup>f</sup> [N/mm] depending on the concrete grades and the

The limit of damage at a point is controlled by the values of the so-called crushing or total failure function Fu. The modified Kupfer's condition for the lth-layer of section is as follows:

<sup>u</sup> ¼

where Iε<sup>1</sup>, Iε<sup>2</sup> are the strain invariants, and ε<sup>u</sup> is the ultimate total strain extrapolated from uniaxial test results (ε<sup>u</sup> = 0.002 in the tension domain, or ε<sup>u</sup> = 0.0035 in the compression domain), and α, β are the material parameters determined from Kupfer's experiment results

The failure function of the whole section will be obtained by the integration of the failure

where tl is the thickness of the lth-shell layer, t is the total shell thickness and Nlay is the number

The maximum strain ε<sup>s</sup> of the reinforcement steel in the tension area (maxðεsÞ ≤ εsm ¼ 0:01) and by the maximum concrete crack width wc (maxðwcÞ ≤ wcm ¼ 0:3 mm) determine the local col-

The program CRACK based on the presented nonlinear theory of the layered reinforced concrete shell was adopted in the software ANSYS [3]. These procedures were tested in

t X Nlay

l¼1 Fl

<sup>u</sup>ðIε<sup>1</sup>; Iε<sup>2</sup>; <sup>ε</sup>u<sup>Þ</sup> dz <sup>¼</sup> <sup>1</sup>

q

<sup>a</sup> � 0:5da þ 26Þðf <sup>c</sup>=10Þ

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi βð3Jε<sup>2</sup> þ αIε<sup>1</sup>Þ

maxÞ ð24Þ

<sup>0</sup>, <sup>7</sup> <sup>ð</sup>26<sup>Þ</sup>

� ε<sup>u</sup> ¼ 0 , ð27Þ

<sup>u</sup>ðIε<sup>1</sup>; Iε<sup>2</sup>; εuÞtl ð28Þ

max ð25Þ

modulus in the form

aggregate size d<sup>a</sup> as follows:

(β = 1.355, α = 1.355εu).

of layers.

Fl <sup>u</sup> <sup>¼</sup> Fl

function through the whole section in the form

lapse of reinforced concrete structure.

Fu <sup>¼</sup> <sup>1</sup> t : ðt

comparison with the experimental results [3, 15, 17, 28].

0 Fl Recent advances and the general accessibility of information technologies and computing techniques give rise to assumptions concerning the wider use of the probabilistic assessment of the reliability of structures through the use of simulation methods in the world [3, 14–22, 29–40]. The probabilistic definition of the reliability condition is of the form

$$\log(R, E) = R - E \ge 0 \tag{29}$$

where g(R,E) is the reliability function.

In the case of simulation methods, the failure probability is defined as the best estimation on the base of numerical simulations in the form

$$p\_f = \frac{1}{N} \sum\_{i=1}^{N} I[\mathbf{g}(X\_i) \le 0] \tag{30}$$

where N in the number of simulations, g(.) is the failure function and I[.] is the function with value 1, if the condition in the square bracket is fulfilled, otherwise is equal to 0. The semi- or full-probabilistic methods can be used for the estimation of the structure failure in the critical structural areas. In the case of the semi-probabilistic method, the probabilistic simulation in the critical areas is based on the results of the nonlinear analysis of the full FEM model for the median values of the input data. The full probabilistic method result from the nonlinear analysis of the series simulated cases considered the uncertainties of the input data.

#### 4.1. Uncertainties of the input data

The action effect E and the resistance R are calculated considering the uncertainties of the input data as follows [3]:

$$E = G\_k \mathcal{g}\_{\rm var} + Q\_k \mathcal{q}\_{\rm var} + P\_k p\_{\rm var} + T\_k t\_{\rm var} \quad \text{and} \quad R = R\_k r\_{\rm var} \tag{31}$$

The uncertainties of the input data were taken in accordance with the standard requirements [29, 30–32] (Table 3).

#### 4.2. Probabilistic simulation methods

Various simulation methods (direct, modified or approximation methods) can be used for the consideration of the influences of the uncertainty of the input data.

In the case of the nonlinear analysis of the full FEM model, the approximation method RSM (response surface method) is the most effective method [3].

The RSM method is based on the assumption that it is possible to define the dependency between the variable input and the output data through the approximation functions in the following form:

$$Y = c\_o + \sum\_{i=1}^{N} c\_i X\_i + \sum\_{i=1}^{N} c\_{i\bar{}} X\_i^2 + \sum\_{i=1}^{N-1} \sum\_{j>i}^{N} c\_{\bar{j}} X\_i X\_j \tag{32}$$

where co is the constant member; ci are the constants of the linear member and cij for the quadratic member, which are given for predetermined schemes for the optimal distribution of the variables or for using the regression analysis after calculating the response. The 'Central


Table 3. The variability of input parameters.

Composite Design Sampling' (CCD) method or the 'Box-Behnken Matrix Sampling' (BBM) method [3] can be used to determine the polynomial coefficients.

The philosophy of the RSM method is presented in Figure 10. The original system of the global surface is discretized using approximation function. The design of the experiment determines the polynomic coefficients. The efficiency of computation of the experimental design depends on the number of design points. With the increase of the number of random variables, this design approach becomes inefficient. The central composite design, developed by Box and Wilson, is more efficient.

The central CCD method is defined as follows (Figure 10):


The total number of design points is equal to N = 2<sup>k</sup> + 2k + no. The sensitivity of the variables is determined by the correlation matrices. The RSM method generates the explicit performance function for the implicit or complicated limit-state function. This method is very effective to solve robust and complicated tasks.

## 4.3. Evaluation of the fragility curve

median values of the input data. The full probabilistic method result from the nonlinear

The action effect E and the resistance R are calculated considering the uncertainties of the input

The uncertainties of the input data were taken in accordance with the standard require-

Various simulation methods (direct, modified or approximation methods) can be used for the

In the case of the nonlinear analysis of the full FEM model, the approximation method RSM

The RSM method is based on the assumption that it is possible to define the dependency between the variable input and the output data through the approximation functions in the

> ciXi <sup>þ</sup><sup>X</sup> N

Quantities Histograms

Dead load G<sup>k</sup> gvar N 1 10 Live load Qk qvar Beta 0.643 22.6 Pressure P<sup>k</sup> pvar N 1 8 Temperature T<sup>k</sup> tvar Beta 0.933 14.1 Model E<sup>k</sup> evar N 1 5 Resistance R<sup>k</sup> rvar N 1 5

i¼1

where co is the constant member; ci are the constants of the linear member and cij for the quadratic member, which are given for predetermined schemes for the optimal distribution of the variables or for using the regression analysis after calculating the response. The 'Central

Input data Charact. value Variable value Type Mean Deviation

ciiX<sup>2</sup> <sup>i</sup> þ N X�1 i¼1

X N

cijXiXj ð32Þ

μ σ [%]

j>i

consideration of the influences of the uncertainty of the input data.

(response surface method) is the most effective method [3].

<sup>Y</sup> <sup>¼</sup> co <sup>þ</sup><sup>X</sup>

N

i¼1

E ¼ Gkgvar þ Qkqvar þ Pkpvar þ Tktvar and R ¼ Rkrvar ð31Þ

analysis of the series simulated cases considered the uncertainties of the input data.

4.1. Uncertainties of the input data

82 Recent Improvements of Power Plants Management and Technology

data as follows [3]:

following form:

ments [29, 30–32] (Table 3).

4.2. Probabilistic simulation methods

Table 3. The variability of input parameters.

The PSA approach to the evaluation of probabilistic pressure capacity involves limit-state analyses. The limit states should represent possible failure modes of the confinement functions. The identification of potential failure modes is the first step of a probabilistic containment overpressure evaluation. For each failure mode, the median values were established based on the used failure criteria dependent on the applied loading consisting of temperature, pressure and dead load. Along with the pressure capacities for the leak-type failure modes, leak areas are to be estimated in a probabilistic manner. The expected leak areas are failure mode dependent. After calculation of the fragility or conditional probability of failure of containment at different locations, we must consider a combination of pressure-induced failure probabilities of different break or leak locations within containment.

Figure 10. Scheme of the RSM approximation method with the CCD design experiment.

Figure 11. Family of fragility curves showing modelling uncertainty.

Containment may fail at different locations under different failure modes (see Figure 11). Consider two failure modes A and B, each with n fragility curves and respective probabilities pi (i = 1, …, n) and qj (j = 1, …, n). Then, the union C = A∪B, the fragility FCij(x) is given by

$$F\_{\mathbb{C}\vec{\eta}}(\mathbf{x}) = F\_{Ai}(\mathbf{x}) + F\_{B\vec{\eta}}(\mathbf{x}) - F\_{Ai}(\mathbf{x}) \cap F\_{B\vec{\eta}}(\mathbf{x}) \tag{33}$$

where the subscripts i and j indicate one of the n fragility curves for the failure modes and x denotes a specific value of the pressure within the containment. The probability pij associated with fragility curve FCij(x) is given by pi. qj if the median capacities of the failure modes are independent. The result of the intersection term in Eq. (32) is FAj(x).FBj(x) when the randomness in the failure mode capacities is independent and min [FAi(x), FBj(x)] when the failure modes are perfectly dependent.

The flow is the consequence of an accident that depends on the total leak area. Multiple leaks at different locations of the containment (e.g. bellows, hatch and airlock) may contribute to the total leak area. Using the methodology described earlier, we can obtain the fragility curves for leak at each location.

For a given accident sequence, the induced accident pressure probability distribution, h(x), is known. This is convolved with the fragility curve for each leak location to obtain the probability of leak from that location (PLi). It is understood that there is no break or containment rupture at this pressure.

$$p\_{Li} = \int\_0^\infty h(\mathbf{x}) [1 - F\_b(\mathbf{x})] F\_l(\mathbf{x}) d\mathbf{x} \tag{34}$$

Here, the Fb(x) is the fragility of break at the location and Fl(x) is the fragility of the leak. The leak is for each location specified as a random variable with a probability distribution.
