Seismic Damage Estimation of an Actual Reinforced Concrete Structure Using Subset MCMC

To estimate seismic damage of structures under strong motions is very important to know true safety of structures. However, we have to deal with very small failure probability issue to investigate quantitatively. If we use standard MCS (Monte Calro Simulation) to discuss failure probability, e.g., 1 × 10 order, since we have to execute nonlinear dynamic response analyses of approximate 10 times, it is not realistic. Recently, a subset simulation, which reduces the computation time by replacing small failure probability into the product of conditional failure probabilities, was proposed. In this study, the subset simulation is applied to estimate failure probability of an actual reinforced concrete building with 11 stories, and discuss the safety of the structure by checking with design criteria.


Introduction
The final goal of structural design is quantitative estimation of structural safety.However, quantitative damage estimation, i.e., failure probability, subject to an individual strong motion is not obtained from the procedure based on a structural design code.Because a specific ground motion is ordinary derived from a acceleration spectrum, which reflects properties of many ground motions predicted at each site, given in the design code.Therefore, we know only whether a structure is safe for the above specific ground motion.
Recently, the subset simulation for reliability estimation of structures subject to strong motions was proposed [1].A basic idea of the subset simulation is to reduce the computation time by replacing small failure probability into the product of conditional failure probabilities.
In this paper, the above simulation is adopted for reliability estimation of an actual reinforced building with 11 stories under the conditions of uncertainty on material properties such as story stiffness and story yield strength.

Failure Probability Using Subset MCMC
The concept of subset simulation proposed by Siu-Kui Au are as follows [2].Given a failure event F , let By the definition of conditional probability, we have The above equation means that a failure probability is replaced by the production of sequence conditional probabilities.There, P(F 1 ) is obtained from standard MCS (Monte Calro Simulation), and P(F 2 ) ~ P(F m ) are obtained from subset MCMC (Markov Chain Monte Calro).Under the condition that the limit state function Z is negative, failure probability using subset methods is expressed by the following equation.

 
where t : number of samples generated in subset n i F s n : number of samples defined in subset , for example 0.1 . F be a decreasing sequence of fai- On the other hand, a conditional probability density function is expressed by the following equation. where

Ground Motions and Analytical Model
Table 1 shows input ground motions.Ground motions of No.1 -No.7 were observed at HKD180 point of K-net.
No.8 -No.11 are design ground motions which are often used for dynamic analysis in Japan.In this study, these ground motions are normalized so that PGV (Peak Ground Velocity) leads to 50 kine, then we call them the level 2 ground motions in Japan.Figure 1 shows an example of normalized ground motions, and Figure 2 shows the response spectra of No.8 -No.11 ground motions.
Analytical model is an actual reinforced concrete residential building with 11 stories and two spans in X and Y direction.Table 2 shows the structural heights and the story weights of this model.
In nonlinear dynamic response analysis, it is assumed that a skeleton curve is tri-linear, and a hysteresis rule is masing type.In Rayleigh damping ratio, it is assumed that first and second damping of natural vibration mode is 5%.Time-stepping solver for the equation of motion adopts Wilson's θ method (θ = 1.42).
As structural failure index, two indexes of MCDR (Maximum Column Drift Ratio) and MDD (Maximum Ductility Demand) are considered.In addition, it is assumed that strong sway observed on tall buildings after ground motion shake down does not occur about 11 stories building, and the following cutoff time is used.where peak-time is the latest time among occurrence times of maximum absolute amplitude on PGA (Peak Ground Acceleration), PGV (Peak Ground Velocity) and PGD (Peak Ground Displacement).
In the subset simulation, it is assumed that material properties are uncertain and input ground motions are deterministic.Concretely, first and second story stiffness, story yield strength, i.e., SK1, SK2, SSy, at each floor are treated as uncertain variables.Furthermore, it is assumed that only two variables in (SK1, SK2) or (SK2, SSy) vary at a simulation in order to simplify.Table 3 shows the mean values of these three variables specified by transforming the relationship between interstory drift and tory shear force, which is derived from pushover analy-s s, into a tri-linear skeleton curve.Each parameter is the si set of i.i.d.(independent identically distributed) according with the following Gaussian distribution: where  : coefficient of variation, it is assumed that  is equal to 5% for the combination (SK1, SK2) and 2% for the combination (SK2, SSy).Calculation flows are as follows.Firstly, standard MCS is executed.Response analysis of t n times per a simulation is repeated in MCS, and then a eakest floor is specified.In the weakest floor, frequencies that interstory drift or ductility demand shows maximum values is the biggest among all floors.
Secondly, the subset simulation is executed.Samples of w the weakest floor are generated using MCMC.t n : the number of samples generated about each subset level, is set 600 based on convergence diagnostics by Raftery and Lewis [3].And to the samples of other floors, the values generated in MCS are applied in order to reduce computation time for the parallel simulation mentioned later.

Basic Response Properties
Natural elastic period of analytical model is 0.764 (sec) in X direction, and 0.8204 (sec) in Y direction.Firstly, we describe response results subject to the 50 kine normalized 1995 Kobe NS shown in Figure 1.
Figures 3(a direction.Und and in X ame.The tendency was er the comparison of ductility dem and Y direction, the damage of X direction frame is obviously larger than Y direction fr same about results subject to other earthquake motions.Therefore, in estimation of failure probability mentioned later, only X direction frame will be described.Incidentally, the ductility demand  in X direction frame is over 1 at other floors except the top floor.appear not only at lower story but also at upper story or middle story.

Sensitivity of Parameter and Quantitative Estimation of Failure Probability
Th a paimit state function, in index, without a spe-

Sensitivity of Parameter
e subset simulation can automatically find out rameter to reduce the value of l other words to make larger failure cial knowledge.However, we know only the relative intensity on sensitivity of two parameters.In sensitivity analysis, parallel simulation is executed one time since the effect by initial value's discrepancy on parameter is small.Then, response analyses are executed 3000 times (=600 in MCS + 600 × 4 in four levels of subset) per a parallel simulation.
Ac n the ilure index is MDD, the occurrence frequency of the ed on calculations of 600 times in M ue μ 2 is No.2 in cording to the results for all input motions whe fa weakest floor is large at 1st floor: its ratio is 46.5%.The following ratios are 13% (9th floor), 12% (7th floor) and 11% (3rd floor) bas CS.While, when the failure index is MCDR, there was no occurrence of the weakest floor at 1st floor.The reason is that MCDR leads to relatively small since the structural height at 1st story is remarkably larger than other floor as shown in Table 2. Therefore, the results on failure index MDD will be hereinafter discussed mainly.
Table 4 shows that the combination of parameters is 1st and 2nd story stiffness for label "a", and is 2nd story stiffness and story yield strength for label "b".
Figure 6 shows sample parameters generated in each subset level.Here, the center of ellipse is the mean val  This tendency is same for other input motions as shown in Table 5.In the case of failure index MCDR, the results of SSy > SK2 are obtained for other ground motions except No.1, 4, 5, 6. While, sensitivity results in combination of parameter "a" are almost SK2 > SK1 as shown in Table 5.Thus, it is seemed that relative intensity of sensitivity is SSy > SK2 > SK1.

Quantitative Estimation of Failure Probability
Firstly, let's explain design criteria assuming in this study.It is general that design criteri MCDR are not valid because there is no occ eakest floor at 1 analytical model.So, it is assumed that design criteria presented by MCDR can be replaced by the criteria of MDD.Table 6 shows the relationship between MDD and MCDR at termination of subset-4 picked up the analyticcal cases with the common collapsed story.The mean value of MDD/MCDR (%) is about 2.2 as shown in this table.
Table 7 is an example of design criteria.After given criteria on MCDR as shown in the left column, MDD criteria result in 2.2 (the above mean value) × MCDR criteria.According to this table, the values of MDD, i.e., the maximum ductility demand, are over the definition of limit state shown in the right column.In the case of not high rise buildings like this analytical model, it is known that the response analysis gives fairly severer results than actual damages.It is seemed that its sign also appeared in the analytical results of this study.
Figures 7 (a)-(c) show relationship between ductility demand and failure probability and its ensemble mean value (black line) derived from parallel simulations.Each line shows the result of MCS (blue line), subset-1 (green line), subset-2 (red line), subset-3 (cyan line) and subset-4 (magenta line) in order from left side.Iteration of parallel simulation is 50 times, i.e., response analyses of 150,000 times (=3000 × 50), per a ground motion.7, the ensemble mean of MDD at termination of subset-4, which corresponds to failure probability = (1/601) × (0.1) 4 = 1.664 × 10 −7 , is over 4.4 of the life safe.Now, if the annual allowable risk for residential buildings is 1 × 10 −6 , and level 2 ground motion has 10% exceedence in 50 years, which means that occurrence probability is 2.1 × 10 −3 (=1/475, i.e., return period 475 years), the value of vertical axis of Figure 7 corresponding to the above annual risk is 4.76 × 10 −4 (=1 × 10 −6 /2.1 × 10 −3 ).The ductility demand corresponding to 4.76 × 10 −4 is 4.3 for No.2 ground motion, and is 3.9 for No.8, these values are slightly smaller than 4.4 of life safe.However, taking into account that the response analysis gives fairly severer results than actual damages as mentioned earlier, it seemed that there is still some margin for the life safe.1) It is seemed that the relative intensity of sensitivity is SSy (story yield strength) > SK2 (2nd story stiffness) > SK1 (1st story stiffness).
2) If an annual risk for residential building is 1 × 10 −6 , the ductility demand for allowable risk of No.2 ground number of Markov chain.

1 -
No.7 earthquakes are observed at HKD180 site.The earthquake names of No.5 -No.7 are temporary.

Figure 4
shows time variation of energy dissipated in X direction frame.According o this figure, yielding does not yet occur although cracking occurs in t the top floor, however yielding occurs in first floor.Figure 5 shows the ductility demands of No.1 -No.11 input motions.According to this figure, the maximum values of  , i.e., max  are over 1 for other ground motions except No.1, 4, 5, and max 

Figure 7 ( 2 R
d) shows the convergence index by Gelman and Lubin to check whether the iteration is enough.R    value is pointed out at least within 1.1 ~ 1.2 according to the convergence judgment [4].Iteration of 50 times satisfies this judgment as shown in Figure 7(d) or e largest

Table 8 .
MDD of No.2, i.e., Kushiro offshore 2004 EW, is th in all ground motions, and No.8, i.e., JMA Kobe 1995 NS, follows next.However, note that this is the results subject to level 2 ground motion.According to Figure

Table 8 . Convergence index and computation time of par- allel simulation.
Kushiro offshore 2004 EW, which gave the largest damage to the analytical building, does not yet reach to the limit state of life safe.tation time per a response analysis was 8.5267 (sec).