Determination of damping ratio using Monte-Carlo method

This paper highlights an innovant method for calculating the damping ratio using Monte-Carlo method. Compared with other approaches reported in the literature, the clear benefit of the method is the ability to directly determinate the damping ratio upon irregular stress-stain curves through Python. By using Monte-Carlo method, it is found that (i) the computed damping ratios are consistent with the undrained cyclic triaxial result on a dense reference sand specimen; (ii) while total random samplings number for each hysteresis loop greater than 5000 times, the method could provide a satisfactory analysis.


Introduction
In the field of materials science, damping ratio is an important indicator directly reflecting the capacity of energy dissipation (e.g., concrete [1], steel [2]). In soil mechanics, when subjected to cyclic shearing in undrained condition, sand could undergo a transient loss of shear resistance owing to the rapid buildup of excess pore water pressure. During this process, sand flows almost like a liquid [3], which could be directly identified by the apparent rise of damping ratio.
The damping ratio [4] is defined as ∆ , where W and W stand for energy dissipated and peak energy in each single cycle, respectively. In practice, the damping ratio (D) is commonly determined in stress-strain curves (see Figure. 1), and its value is proportional to area ratio between the dissipated energy (blue hysteresis loop) and the peak energy (red triangle). In this sense, a traditional method is to sequentially calculate D one after another by assuming that hysteresis loop could be possibly viewed as an ellipse. Unfortunately, this method is extremely time-consuming since the cycle number N cyc is usually huge. In order to solve this problem, Wei [5] deduced that D could be roughly calculated by a simple formula: , where  d ,  m stand for intersect stress value with y-axis and maximum stress in one cycle, respectively. Later, Chen [6] further pointed out that the polygon had a better fitting performance for hysteresis loops, as compared with ellipse. However, numerous experimental studies [7] indicate that the hysteresis loop is definitely irregular, especially in large

Monte-Carlo method
Originated in the 1940s, Monte-Carlo method (also referred to as "Monte-Carlo experiments") relies on a great amount of repeated random samplings for the purpose of obtaining a probability as the solution to a given event that might be deterministic. Nowadays, with the rapid development of calculation power, the mimicking of random samplings could be quite easily realized in many commercially available platforms. In particular, Python has been widely proved to have greater advantages in this field since it is totally Open-Source and object-orienting. Thereafter, Monte-Carlo method based on Python has gradually become one of the standard tools in physical or mathematical problems when it is extremely difficult or seemingly at least impossible to apply other approaches. In image processing, this method is commonly employed in calculating the surface of irregular graphics. Thus, the area of irregular hysteresis loops (S loop ) linking to damping ratio of sand response is considered as a reference instance in the present work, and the main processes could be divided into the following steps ( Figure. 2): 1. Split the target loop from the complete stress-stain curve with a constant external frame. The boundary values of external frame are respectively x min , x max and y min , y max for x-and y-axes, implying that the surface of external frame S tot is equal to (y max -y min )×(x max -x min ). During this process, it should be noticed that the after-splitting single loop needs to be strictly enclosed. In other words, the presence of intersection point (represented by solid green point in Figure. 2) is imperative. 2. Generate random points (x rand , y rand ) with a total samplings number T ot in the coupled range of [x min , x max ] and [y min , y max ]. 3. For each random point, interpolate (e.g., linear manner) the upper and lower boundaries of the hysteresis loop in terms of the x-coordinate of random point (x rand ) to find the corresponding y intmin and y int-max . 4. Estimate the relative position of random point to hysteresis loop. If y rand <y int-max and y rand >y int-min could be simultaneously satisfied, it could be logically deduced that this random point is an  Figure. 2. Note that, in order to render all judgements regarding the relative position achievable, the hysteresis loop needs to be geometrically enclosed, as previously mentioned in the 1 st step. 5. Repeat the 3 rd and 4 th steps until all generated random points are treated. 6. Count the number of inside random points T ins . Then, the probability of these points is equal to T ins /T ot and, at same time, equal to S loop /S tot . Accordingly, the desired hysteresis loop surface (S loop ) could then be calculated: (S tot ×T ins )/T ot .
It is worth noting that although the probability is dependent on the used external frame (e.g., size and shape), the final result (S loop ) is still irrespective of the frame. Besides, the more random points there are, the more accurate the final result becomes since one might expect that the probability tends to gradually approach the true value with the increasing number in random points.

Computation results
The stress-stain curve (q/ a : deviator stress/axial strain) of undrained cyclic triaxial tests on a clean dense HN31 sand specimen (density index of sand matrix I Dmat =0.70, cyclic stress ratio T cc =0.22 and effective confining pressure  ' c =200 kPa) is depicted in Figure. 3. Note that detailed index properties of the used HN31 sand and experimental program are reported by Zhu [8] and Jradi [9], respectively. Three distinct behaviours can be distinguished as follows: (a) At the very beginning of test (e.g. | a |≤0.5%), the hysteresis loops are almost linear and coincide with each other. This suggests that in this minor range of  a , the sand response could be viewed as "pseudo-elastic" without significant energy dissipation.
(b) With the further development of  a , the plasticity of sand becomes much more important, which can be embodied by the increasing discrepancy in the subsequent loops.
(c) In the extremely large range of  a (e.g., | a |≥5%), the hysteresis loops become quasi-"Z" shape, exhibiting a strong dissipation capacity. This is caused by the triggering of liquefaction in which the sand specimen behaves almost like a liquid due to the rapid build-up of excess pore water pressure (see the subplot of Figure.  Due to the constraints of limited space, only several treated hysteresis loops are schematized in Figure. 4 where the damping ratio D (determined by Monte-Carlo method of total random samplings number T ot =5000) and the corresponding cycle number N cyc are marked at the top-left corner, respectively. It can be seen that (i) all loops are well enclosed; (ii) for N cyc altering from 1 to 9, the damping ratios are very low since their surfaces are extremely minor (pseudo-elastic); (iii) from N cyc equalling to 13, the hysteresis loops become more obvious with the significant increase in D; (iv) in the stage of large deformation (e.g., N cyc =24, 27), the two corresponding damping ratios attain the high level of 55.79% and 59.66%, reconfirming the occurrence of liquefaction softening.  v =( a +2 r )/3=0, where  a / r stand for axial/radial principal stresses and  v represents volumetric strain.
In the case of T ot =25, the damping curve is characterized by a fake peak value of about 82%, which is contradictory to the sand behavior identified in Figure. 3. This misleading could be attributed to the fact that the random points are too few (only 25) to gain an accurate probability related to the damping ratio.
In the cases of T ot =500, 5000, 10000 and 50000, it is found that (i) the fake peak damping ratio disappears; (ii) with the increase in T ot , the curves gradually converge to a single one, especially for T ot =5000, 10000 and 50000; (iii) three distinct responses of sand specimens could be clearly identified, which is consistent with experimental results in Figure. 3. While assessing damping ratio using Monte-Carlo method, the above phenomena further suggest that (i) T ot is the most prominent factor and (ii) a value higher than 5000 times could, in general, ensure a satisfactory accuracy.

Conclusion
This paper investigates an innovant method for determining the damping ratio with Monte-Carlo method. Compared with other approaches presented in the literature, the great advantage of the method is directly treating irregular hysteresis loops without any geometrical assumption. In addition, it is found that the computed result is highly sensitive to Tot and a value of Tot higher than 5000 times could, in general, provide an accurate result.