Optimal Sensors Placement in Dynamic Damage Detection of Beams Using a Statistical Approach

Structural monitoring plays a central role in civil engineering; in particular, optimal sensor positioning is essential for correct monitoring both in terms of usable data and for optimizing the cost of the setup sensors. In this context, we focus our attention on the identification of the dynamic response of beam-like structures with uncertain damages. In particular, the non-localized damage is described using a Gaussian distributed random damage parameter. Furthermore, a procedure for selecting an optimal number of sensor placements has been presented based on the comparison among the probability of damage occurrence and the probability to detect the damage, where the former can be evaluated from the known distribution of the random parameter, whereas the latter is evaluated exploiting the closed-form asymptotic solution provided by a perturbation approach. The presented case study shows the capability and reliability of the proposed procedure for detecting the minimum number of sensors such that the monitoring accuracy (estimated by an error function measuring the differences among the two probabilities) is not greater than a control small value.


Introduction
Strategic constructions, such as bridges and tower buildings, during their lifetime are inexorably subjected to degradation and damage of various kinds. For this reason, careful monitoring is required for profitable maintenance and for activation of safety protocols in case of danger ("early warning" systems). The main problem, that discourages the use of large-scale structural health monitoring (SHM) system, is its high cost, directly related to the cost of the instrumentation and in particular the number of sensors used.
In this framework, dynamic approaches-mainly concerned with output-only techniques-have been successfully adopted to locate and quantify structural damages for historical [1,2] and recent construction [3,4]. Other methods have been proposed for the damage detection such as data-driven [5,6] and statistical method [7].
As reported by many authors, when structural responses of dynamic systems are exploited to trace the signature of the damage, one deals with an inverse problem, where damage influence is recognized starting from its effect. On the one hand, techniques not model based, see for instance [8,9], are able to directly detect the alterations due to any variation, but their mechanical interpretation is often cumbersome. On the other hand, model-based approaches, such as the techniques based on one-dimensional elements [10,11], usually need to rely on accurate structural modeling and to select proper response signals. Damage effects have been also investigated, using non-classical continuum approaches, developed for materials with microstructure, by investigating displacement fields and wave propagation [12][13][14].
Evidently, a well position of the sensor network defines the suitability of dynamic behavior of the analyzed structure.
Therefore, a correct treatment that aims at reducing the number of sensors by limiting the cost of SHM without compromising the quality of monitoring is desirable. For this reason, optimal sensor placement (OSP) is an opening challenge in engineering [15,16] and many different strategies have been proposed by several authors. The review paper [17] reports a comprehensive overview of a number of criteria for evaluating different sensor setups such as modal assurance criterion (MAC), the information entropy (IE) [18][19][20][21], singular value decomposition ratio (SVDR) [22], Fisher information matrix (FIM) [23,24], system-realization method [25] and multi-setup modal analysis [26]. A Bayesian framework for model-based optimal sensor placement for response predictions has been presented in [27,28]. A comparison of several methods for optimal sensors placement has been proposed in [29] applied to a real case study of a suspended bridge. In a recent work [30] a study of optimal sensor placement for timber frame structure has been developed. In [31] the conditions for the invertibility of linear system models have been discussed. In these methods, an optimization process has been performed to minimize an objective function. In order to solve the optimization process, different methodologies have been adopted such as improved genetic algorithms [32,33], and discrete optimization [34].
In the above-mentioned works, attention is paid to experimental modal analysis; however, there are few works on the identifiability of damage processes suffered by a structure. The recent work [35] focused on searching the weakest part of the structure in a random context. In this work, non-localized damages are adopted, which are in general more difficult to identify than narrow damages [36,37]. Damages are modeled as a perturbation of the healthy system and therefore a perturbation approach can properly describe the effect of the damage, both for standard [38] and state-space dynamic analysis [39,40]. In [41], a perturbation method has been proposed by some of the authors to derive the asymptotic eigensolution of a vibrating beam with uncertain damage. Then, the statistics of the random damage parameter have been obtained by means of an objective function minimizing the difference among analytical and experimental fractiles of the eigenvalues. The results obtained by the authors show the importance of considering the second order terms of the perturbation approach in order to achieve a significant increase in the identifiability of stochastic damages.
Starting from this model, we here propose a technique for obtaining an optimized solution of sensors placement. The proposed approach aims at exploiting the closedform asymptotic solution of the inverse problem to compare more combinations of number and placement of the sensors. These comparisons rely on an error function measuring the deviation from the analytical probability of damage occurrence. Due to the asymptotic nature of the research pattern, the performance of the proposed placement optimization technique has a low computational effort. Some preliminary results have been shown in [42]; the peculiarities and the capabilities of the technique are here deepened.
The main novelty of the technique, with respect to the cited bibliography, concerns the closed form evaluation of the identifiability of uncertain (local or widespread) damages. Even if the technique has been so far developed only for vibrating beams, the methodology is generally suitable for all cases when the dynamics of a structural system can be described through analytical solutions.
The paper is organized as follow: In Sect. 2, the details of the proposed optimal sensor placements have been presented. In Sect. 3, a numerical example, exploiting the capabilities of the proposed method, has been shown. Finally, in Sect. 4, some final remarks and future developments are reported.

Beam and Damage Model
We consider in this work the classical beam model of the Euler-Bernoulli in which the shear deformation effect is assumed negligible. The governing equation of the free flexural undamped vibration reads: where B is the flexural stiffness and ρ = ρ V A is the mass density per unit length (ρ V and A being the volume mass density of the material and the area of the cross section, respectively). x is the abscissa along the beam, of length L, v is the transversal displacement, and ε is a random variable related to damage (Fig. 1). In addition, the following assumptions were considered. Negligible changes of the mass, i.e., ρ(x, ε) = ρ. That is a fairly common hypothesis when one deals with damage identification. 2. Gaussian distribution of the random parameter ε with mean value μ ε and standard deviation σ ε , i.e., ε ∼ N (μ ε , σ ε ). The random parameter ε defines the damage state of the beam and the related damage probability P d can be evaluated as: where f (ε) is the probability density function of ε. In detail, since the random parameter ε has been assumed Gaussian distributed, the probability of damage P d depends only on the ratio among the standard deviation σ ε and the mean value μ ε , that is, P d is a function of the so-called coefficient of variation C V ε = σ ε /μ ε . 3. Linear decomposition of the flexural stiffness as follow: where μ b characterizes the position around which the damage is centred and σ b is the relevant extension. The 'small' perturbation parameter ε linearly scales the stiffness variation B 1 (x) of an initial (uniform) bending stiffness B 0 . In our applications we pose a negative function B 1 (x), so that the random parameter ε can be regarded as a random damage, where ε = 0 implies an healthy beam, whereas a positive value of the parameter indicates a damage state (the higher the parameter, the greater the damage). We remark that B 1 (x) is a function establishing where the damage is located and how wide it is, and that the resolutive equations shown below are still true even for functions other than the bell-shaped one proposed here.
These positions reduce the entirety of the problem, they provide a useful simplification, allowing for a proper description of stochastic reductions (damage) and increases (reinforcement) of the beam stiffness [42].
With these assumptions and expanding all the terms in Eq. (1) one gets: For the sake of simplicity, we adopt the following notation for the spatial (•) = ∂ (•) /∂ x and time( •) = ∂ (•) /∂t derivatives, respectively.
Introducing of the modal projection: where q i (t, ε) = C i sin(ω i (ε)t +θ i ) and φ i (x, ε) are the mode shapes (eigenfunctions) and ω i (ε) the relevant (circular) frequencies related to the eigenvalues λ i : The governing equation in the modal space reads: Since ε is a small parameter, Eq. (6) can be solved through a perturbation approach. Moreover, the dynamic system under investigation is non-defective, and therefore, the eigensolution admits a perturbation series of integer powers of the small parameter: Defining the operator L[•] = −B 0 λ 4 0 (•) + B 0 (•) and balancing the same powers in ε, Eqs. (6) and (7) can be written as a series of equations as follows (where the dependence on i has been omitted): The solution of Eq. (8) 1 -so-called generating equation-is the well-known one: The boundary conditions posed (BC) allow us to compute the eigenvalues λ 0i and three of the four constants C j , j = 1, . . . , 4 (the fourth one is a generic constant scaling the mode). Once the pair (λ 0i ; φ 0i (x)) are known, Eq. (8) 2 gives (λ 1i ; φ 1i (x)), Eq. (8) 3 gives (λ 2i ; φ 2i (x)) and so on. After some algebra, the expressions of the eigenvalues and the eigenfunctions for the first and second order terms can be calculated: for the first and second order terms of the eigenvalues, and: for the first-and second-order terms of the eigenfunctions.

Optimal Sensor Placement (OSP)
Let us focus on a beam defined by its stiffness B(x, ε), linear mass density ρ(x, ε) = ρ and length L, such as described in Sect. 2.1. The solution of Eq. (7) gives a closed-form expression of the eigenvalues and eigenfunctions here adopted to estimate the identifiability of damages in the vibrating uncertain beam. Moreover, the damage state of the beam is governed by the random parameter ε, that has been assumed gaussian distributed (see Sect. 2.1), and the probability of damage P d , defined in Eq. (2) is a function of the coefficient of variation C V ε = σ ε /μ ε . The probability to detect the damage resorting to a statistical dynamic analysis is related to the reduction of the eigenvalues: being n λ the number of the first (measured) eigenvalues (related to sampling frequency of the sensors). Using Eq. (7) 1 , the differences λ i (ε) − λ 0i can be written as: , whereas the term (ε/σ ε ) 2 has a noncentral chi-squared distribution X 2 (k χ , λ χ ), with a unitary degree of freedom k χ and non-centrality parameter λ χ equals to (μ ε /σ ε ) 2 . Therefore, the probability density function of the term λ i (ε) − λ 0i is a mixture of two known distributions. For the sake of brevity, the relevant probability P ddi is here explicitly evaluated only in the case when the two derivatives λ 1i and λ 2i are both negative (which is the most usual case): with i = 1, . . . , n λ and where Erf and Erfc are the error function and the complementary error function, respectively. The two weights w 1i and w 2i can be expressed as: The penalty indexes affecting the first and the second order term are: respectively. These two indexes are introduced to describe in a synthetic way the accuracy in the reconstruction of the variations encountered by the mode shapes. Several definitions can be adopted for this purpose. In this work, penalty is related to the discretization error of experimental mode shapes, i.e., γ 1i and γ 2i are evaluated comparing the continuous (analytical) mode shapes and the discretized (measurable by the network) ones. A 1i and A 2i are the analytical absolute values of the areas under the graph of |φ 1i (x)| and |φ 2i (x)|, respectively, whereasÃ 1i andÃ 2i are the discretized counterparts: where x k and n s are the positions of sensors along the beam axis and the number of sensors, respectively. An example of the areas under the graph of |φ 1i (x k )| for different values of the number of sensors is reported in Fig. 2, which shows that γ 1i and γ 2i tend to 1 (which implies no penalty) when the number n s of the sensors tends to infinity.
The discrepancy between the probability of damage occurrence, P d defined in Eq. (2), and the probability to detect the damage, P ddi defined Eq. (13), measures the capabilities of the sensor network to detect the changes in eigenvalues and eigenfunctions due to structural damages (that is, the damage): Eq. (17) is a function of the number n λ of measured eigenvalues, the positions x k of sensors along the beam axis and the number n s of sensors, i.e., ΔP = ΔP(n λ , n s , x k ).
Assuming that an optimal sensors placement is the one requiring less sensors, the best sensors network for a given number n λ of measured eigenvalues can be defined through the following constrained optimization problem, where the objective function is not-strictly convex and constraint is an "hard-like" nonlinear inequality: determine min F (n s ) = n s subjected to: with the bounds on the unknowns: where T ol = ΔP is a chosen threshold value for the difference among the probability of damage and the probability of detecting the damage. The proposed procedure is resumed in the flowchart of Fig. 3, where E 0 is the initial value of the Young's modulus and I 0 the initial value of the cross-sectional moment of inertia (so that the undamaged bending stiffness B 0 is equal to E 0 I 0 ).

Case Study
The case study is a simply supported beam, which geometry and material in the With these values we intend to represent usual real-life cases: when one needs to set a sensor network to monitor the health of a structure (beam-like in our case), both the position in which the damage is likely to occur and its magnitude can be postulated or estimated quite easily (for instance, with nonlinear models in ultimate limit state conditions). From this standpoint, the central section is certainly the most prone to damage for a simply supported beam (hence the value L/2 for μ b ); at the same time, a 20 % reduction in stiffness (μ ε of 0.20) is also a typical and reasonable value.   Conversely, it is difficult to model and calculate the spatial amplitude σ b of the damage and the scattering (uncertainty) σ ε of its magnitude. For this reason, here we assume a wide range of possible values for σ b (from L/60 to L/15) and σ ε (from 0.04 to 0.36).
To give an idea of the damage that is being considered, Fig. 5 shows the trend of the bending stiffness B(x, ε) with the variation of the spatial amplitude σ b (considering for the random parameter intensity ε a value of 0.20), whereas in the case where σ b is equal to L/60, the damage is practically concentrated at the center section, for σ b equal to L/15 the damage extends for almost a third of the beam.
As regards the chosen values of the standard deviation σ ε of the random parameter ε, the relevant probability of damage P d is shown with the help of Fig. 6: when the coefficient of variation attains the values C V ε = {0.2, 1.0, 1.8}, the probability of damage P d is equal to 100, 84, and 71 %, respectively. It clearly turns out that the first case can be assimilated to a deterministic damage.
To run the flowchart of Fig. 3, two other parameters must be set: -the tolerance value T ol for the probability discrepancy ΔP. Here we use a value of 10 %, which seems to be a good compromise between the wish of prediction accuracy and the need to contain the number of sensors; -the number of measured eigenvalues n λ . Three values are here adopted, n λ = {2, 4, 8}, well representative of the number of eigenvalues typically measured in real beams.

Results
In this section, we show the results of the constrained optimization problem in Eq. (18).
Since three different values have been assumed for the coefficient of variation C V ε , for the spatial amplitude of the damage σ b and for the number of measured eigenvalues n λ , a total amount of 27 runs are here performed. As also shown in the flow-chart of Fig. 3, the results of the technique are provided in terms of the minimum number of sensors n s and their position x k , k = 1, . . . , n s along the beam axis. To pose a computationally simpler problem, we assume hereafter that networks are composed of equally spaced sensors: under this hypothesis, the solution is just given as n s , i.e., the minimum number of sensors.
With reference to the case of four measured eigensolutions (n λ = 4) and coefficient of variation of 0.20 (C V ε = 0.20), Fig. 7 shows the variation of the percentage probability difference ΔP (Eq. (17)) when n s ranges from 1 to 10. When the number of sensors is increased, there is a "physically sound" reduction of the percentage difference; moreover, as even suggested by the mechanical meaning of the spatial amplitude σ b , the reduction with n s is more pronounced when σ b is greater, that is, when the damaged area is increased. Recalling the adopted threshold T ol of 10 %, for all these three cases, we obtain as optimal (equally spaced) sensors placement the one considering eight sensors. Table 1 shows the results obtained for the 27 runs of the constrained optimization problem. In detail, the table shows the minimum number of sensors n s needed to ensure a percentage probability difference ΔP less than the adopted threshold T ol of 10 %. As a general comment, within the ranges of variation considered in this application, the minimum number of sensors n s changes from 6 to 10; the plots shown in the previous Fig. 7 are related to the values highlighted in the table through a gray background. Some significant remarks can be collected: -following the previous comments on Fig. 7, for a given number n λ of eigenvalues and a given coefficient of variation C V ε , the minimum number n s tends to decrease if the damage is more widespread (i.e., when σ b is increased); -for a given number n λ of eigenvalues and a given spatial amplitude σ b of the damage, the minimum number n s tends to decrease if the damage is more uncertain (i.e., when the coefficient of variation C V ε is increased). This evidence is to be connected to a simple circumstance: as the coefficient of variation C V ε increases, the probability of damage P d decreases, therefore, having imposed a constant tolerance T ol (between the probability P d and the probability to detect the damage P ddi ), the minimum number n s must necessarily decrease; -the last condition to be analyzed is the one revealed considering a given spatial amplitude σ b of the damage and a given coefficient of variation C V ε of the random parameter ε. The results show that when the number of eigenvalues n λ increases, also the required minimum number of sensors n s increases. At a first glance, this result may seem counterintuitive, because one would expect that as the number of information (the number of eigenvalues) of the single sensor increases, the network consequently needs fewer sensors to achieve a certain performance. In truth, the result here obtained is fully consistent with the proposed procedure. Indeed, our technique provides as criterion for selecting the minimum number of sensors n s , the evaluation of a difference in probability which is an average value considering all the contributions of the measured eigenvalues (see Eq. (13) and Eq. (17)). Taking into account that these contributions of the individual eigensolution are Table 1 Optimal (minimum) number of sensors n s for different combination of: number of eigensolutions n λ , coefficient of variation C V ε of the random parameter ε and spatial amplitude σ b of the damage weighted on the discretization error of the perturbed mode shapes (see the penalty indexes of Eq. (13)), it is evident that as the number of eigenvalues is increased, even the required minimum number of sensors n s increases (because it is required to describe higher, "more distorted", modes). In other words, to see a benefit linked to a greater number of measured eigenvalues, a different relation should be used for the difference in probability Eq. (17), for example, instead of the average value, the minimum value of the difference between the damage probability P d and the probability to detect the damage P ddi could be considered.

Conclusions
The main goal of this paper was the damage identification in transversely vibrating beams with uncertain stiffness. The beam model is non-deformable for any shear deformation and changes in mass due to structural damages are assumed negligible. Damage is then introduced as a "small" perturbation of the uniform initial bending stiffness. In detail, the damaged scenario is described with two terms: a bell-shaped function, ruling the position around which the damage is centred and its extension; and the "small" random parameter, which controls the intensity of the damage. Then, the adoption of a modal projection in conjunction with a perturbation technique allows us to achieve the eigensolutions in closed-form.
Moving towards the sensor network, the hypothesis of Gaussian distribution of the damage random parameter leads to a robust procedure for an optimal placement of sensors. It is based on the comparison among the probability of damage occurrence, evaluated from the known distribution of the random parameter, and the probability to detect the damage, evaluated exploiting the closed-form asymptotic solution described above. This comparison brings us to a constrained optimization problem, where the constraint is ruled by an error function encompassing three parameters: the number of measured eigenvalues, the positions of the sensors along the beam axis and the number of sensors.
A simple and paradigmatic case study, whose behavior is well known and clear, was considered to highlight the advantages and disadvantages of the proposed technique. The parametric investigation of a simply supported beam with rectangular cross section showed how the proposed approach is a suitable tool to collect and compare the performance of different sensors configurations. On the other hand, the main limitation of the proposed technique is the adoption of analytical solutions. This approach excludes its generalization to structures more general or complex than beam-like, for which closed form solutions are not pervasive.
Starting from the procedure here developed, there are some future developments capable of improving and generalizing the procedure itself. In particular: -here the optimization refers to a given, expected, damaged scenario, but more plausible damage models would be desirable; -a penalty index taking into account the sensitivity of the sensors should be introduced to describe in the procedure the lower frequency shift detectable by the sensors; -the definition used to measure the accuracy of the sensor networks, based on an average difference (Eq. (17)), could be "weaken" or "strengthen" (looking at the lower or at the greater difference of the summation); the relevant effects on the optimal sensors placements are being analyzed; -an extension of the proposed methodology for 2D continua is under investigation in order to consider structural elements more general than beam-like (slabs, plates, shells, and so on).