{'fraction_non_alphanumeric': 0.029908043924649585, 'fraction_numerical': 0.009329524149629497, 'mean_word_length': 4.977321237993596, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
On the simulation of the seismic energy transmission mechanisms

In recent years, considerable attention has been paid to research and development methods able to assess the seismic energy propagation on the territory. The seismic energy propagation is strongly related to the complexity of the source and it is affected by the attenuation and the scattering effects along the path. Thus, the effect of the earthquake is the result of a complex interaction between the signal emitted by the source and the propagation effects. The purpose of this work is to develop a methodology able to reproduce the propagation law of seismic energy, hypothesizing the"transmission"mechanisms that preside over the distribution of seismic effects on the territory, by means of a structural optimization process with a predetermined energy distribution. Briefly, the approach, based on a deterministic physical model, determines an objective correction of the detected distributions of seismic intensity on the soil, forcing the compatibility of the observed data with the physical-mechanical model. It is based on two hypotheses: (1) the earthquake at the epicentre is simulated by means of a system of distortions split into three parameters; (2) the intensity is considered coincident to the density of elastic energy. The optimal distribution of the beams stiffness is achieved, by reducing the difference between the values of intensity distribution computed on the mesh and those observed during four regional events historically reported concerning the Campania region (Italy).


Introduction
After several seismic events, anomalous distributions of buildings damage, with strong differences between adjacent areas, can be observed. It can be caused by different vulnerability of buildings, as well as by changes in the characteristics of ground motion due to local geological and geomorphological conditions. It is well known that the seismic wave propagation is influenced by: the characteristics of the earthquake excitation (magnitude, type of fault breaking and source-site distance) and, the 2 local response in function of very complex interactions between seismic waves and local conditions, generally known as "site effects". More specifically, spatial variation in seismic ground motions is manifested as measurable differences in amplitude and phase of seismic motions recorded over extended areas. It has an important effect on the response of lifelines, such as bridges, because these structures extend over long distances parallel to the ground and their supports undergo differential motions during an earthquake (Fontara et al, 2015;Fontara et al, 2017). A challenge for the research networks of different areas has been placed worldwide with regard to evaluate the seismic attenuation laws, in order to determine the seismic motion in a given area, once the motion in a reference site is known. The formulation of physicalmathematical models that represent the actual complexity of the phenomenon and their resolution, has considerable difficulties. With the exception of a very limited number of simple situations, the physical-mathematical approach does not lead to solutions in a closed form, that is, integrally analytical (Douglas, 2014 and2017). Therefore, in general it is necessary to adopt numerical procedures, often with approaches based on continuous discretization processes.
Various researchers have published reviews of ground-motion simulation techniques (e.g. Aki, 1982;Shinozuka, 1988;Anderson, 1991;Erdik and Durukal, 2003;Douglas and Aochi, 2008). The models for the prediction of earthquake ground motions are based on two approaches (Ólafsson et al., 2001): (1) the 'mathematical' ones, where the model is analytically based on physical principles; and (2) the 'experimental' ones, where the mathematical model, which is not necessarily based on physical insight, is fitted to experimental data. In addition, there are hybrid approaches combining elements of both approaches. Among bidimensional (2D) numerical procedures, the Finite Element Method (FEM) is frequently used to assess the seismic wave propagation (Lysmer and Drake, 1972;Bao et al., 1998;Ma et al., 2007). The continuous domain (or reference site") is dividing into an equivalent system of smaller subdomains (of volume, surface and beam),.taking into account the soil heterogeneity, the boundary conditions and the not linear behaviour of materials. In this perspective, the present research activity is developed. The main objective is to furnish a methodology able to define the propagation modalities of the seismic energy, through a structural optimization procedure with a predetermined energy distribution (Baratta and Zuccaro, 1990). The approach, based on a physical model of a deterministic type, determines the seismic intensity distribution on the territory, based on objective correction, forcing the compatibility of the observed data with the physical-mechanical model.

Problem definition
The seismic event is identified with the elastic energy E(e,s), symmetric, accumulated in proximity of the site "s" caused by a seismic event at the epicentre "e"; the measure of the seismic event is determined by suitably correlating said energy with an intensity grading, according to any macro-seismic scale (MCS, MSK, EMS). The investigated soil, assimilated to continuous elastic, is discretized into a plane reticular structure formed by braced square meshes with side of approximately 10Km (Figure 1).
The solution to the problem is pursued by researching the function E(e,s), applying the laws of reciprocity which govern the theory of elasticity, with the object of obtaining a prefixed distribution of energy, coincident with the aforesaid, deriving from macro-seismic intensities observed during previous seismic events. In practice, we must resolve a problem of structural optimizing at known geometry. In the present study, the Campania-Lucania territory is investigated. From the Isoseistic Atlas of the P.F.G. -C.N.R. (1985) four seismic events registered during the 20th century in the Campania-Lucania Appennine region have been located ( Figure 2). For each seismic event, from the quoted plane of the intensities, the values of the nodal intensities I oi of the equivalent reticular structure have been calculated, averaging the distance of the node i-th from the closer isoseismic curves, that is, by interpolating linearly. As may be noted, the calculated values of the observed nodal intensities I oi are not expressed by integer numbers -as defined by any macro-seismic intensity scale (MCS, MSK, and EMS). However such an assumption turns out to be less arbitrary if we consider the intensity as a measure of the propagated seismic energy. Such observed intensities I oi are correlated to the elastic energy  oi stored in proximity of the generic node of the mesh, through the empiric relation: (1)

Analytical treatment
For each seismic event, we first locate within the mesh the Cell nearest to the epicentral area, and we simulate the seismic event by a self-balanced system of distorting parameters, applied to the 4 nodes defining the cell (Figure 3). Detailed description of the terms that define the problem is given in Table 1.   Vector of the free nodal displacements f are unknown. At solution, the previous variables have to fulfil equilibrium and compatibility of the deformed configuration, compatibly with a suitable simulation of the assigned energy distribution.

Equilibrium and compatibility
Equilibrium is enforced by setting: Which, as already defined in the table 1 is equals to: Eq.
(2) represents a system of "n" equations in the "N" unknown n. However, in general, the number of the unknowns is greater than the number of the equilibrium equations, (N>n) so 6 that the problem is hyperstatic; in order to determine the hyperstatic unknowns. A partition of the compatibility matrix can be done as follows: in which the first "n x n" lines the matrix   T n n, 1 B is square and non-singular, while the second The same operation can be performed on the axial strain vector n collecting the first "n" entries in a vector n 1 and the additional "N-n" entries in a vector n 2 : Thus eq. (2) is partitioned as follows: The matrix   T n n, 1 B is square and usually invertible. Hence it is possible, to express the previous system of equations in a staggered from by inverting the   T n n, 1 B . In fact, solving eq. (6) with respect to n 1 and denoting for conciseness the vector of the axial hyperstatic strains n 2 (N-n) as Setting further: is the identity matrix and 0 (N-n) is the null vector, The vector n(n) can me expressed as: In this way the axial strain vector is expressed as the sum of a known vector n 0(N), in balance with the nodal and an additional distribution of axial strains, depending the forces hyperstatic unknowns. The known axial strain n 0(N) may be conveniently regarded as the sum of the axial strains in the beans associated with for each distortive parameter: and eq. (12) eq. (9) can be written as: In conclusion the axial strain n i in the i-th bean can be expressed as: 8 In order to enforce collimation between the observed energy and the computed one we need to introduce. The elastic elongation in the i-th beam is introduced and defined as: Where bij denotes the entries of the compatibility matrix B.
(n) n) (N, (N) The "object" energy distribution Introducing the matrix of connection matrix ) , ( N n C whose generic entry is defined by: Assuming that the elastic energy of each beam equally splits between the two nodes to which it is connected, the energy at each node is given by the well known Clapeyron formula: The energetic collimation between the observed energy ) n (

 and the computed one )
(N x is expressed as: In particular, the energetic comparison at the i-th node reads: Eq. (19) of energetic collimation is a system of " n "equations in the "N" unknowns the vector (N) x , in general the number of the unknowns is greater than the number of the equations making the problem hyperstatic. In order to compute the hyperstatic unknowns it is convenient to partition the connection matrix as follows: Thus, the energetic collimation equation can be reformulated (19) by writing: The solution of the previous systems of equation is searched by pivoting on ) ( 1 n x by introducing the vector t of the free energetic collimation variables: Thus the first " n " energy unknowns can be expressed as function of the remaining "Nn " unknowns. In particular, setting: This operation has the purpose of isolating the hyperstatic variables, in the sense that the vector of the beam energies is expressed as the sum of a known vector x 0(N) , obtained from the observed energies, and from a further energetic distribution dependent on the free variable unknowns of energetic collimation t.

Objective function
The unknowns of the problem are actually evaluated by searching the minimum of the average square error E between the calculated energies and those observed in each beam; on account of the equilibrium equation (2) and compatibility one (15) the average square can be expressed: Whose minimum reached is under the condition ensuring the positivity of the rigidity of the ith beam: As a matter of fact, the solution is obtained by solving a structural optimization problem (Non-Linear bounded Programming, NLP) for a given geometry. In conclusions, the mathematical problem is to solve amounts to evaluate the minimum of non-linear function of several variables subject to nonlinear constraints. Among several procedures available in the literature, the Penalization method is selected. Since it amounts to search the constrained minimum by solving a unconstrained optimization problem based on suitably modified objective function that incorporates the constraints as a penalization term.
For the case at hand thus term reads so that the unconstrained minimum of the functional: Therefore, researching the free minimum of the function: retrieves the constrained minimum of the functional E when   r 7

Formulization of the inverse problem
Having indentified the elastic proprieties of the equivalent reticular structure (mesh), it is possible to run the inverse procedure. Once derived the distortional parameters f at the epicentre, an elastic analysis of the mesh is done, on the basis of the values of the calculated optimal stiffness: where f (n) is the vector of the active loads (known): is the vector containing the unknown nodal displacements, with which the examined area has been discretized, matrix D (N,N) contains the constitutive low of the material; K (N,n) is the stiffness matrix and B (N,n) is the compatibility matrix eq. (16). The solution to the eq. (31) is achieved by inverting the matrix K (N,n) arriving at the following relation: By which, the elastic energy stored is obtained in each node of the equivalent reticular structure through: In which: Assuming that, the vector of the nodal elastic energy is known, nodal intensities can be evaluated through the inverse of the eq. (1):

Optimization
The problem is developed assuming the following Mono Objective Optimization of the four events considered: With e=event=1,…, 4; by which the distortive parameters of the four epicentre cells are assumed known, referred to the events considered. The optimization of the eq. (37) is carried out by the Descent Methods, iterative methods, which, as known, employ information on the inclination f(x) of the objective function f(x) to determine a direction of research, d k . Such methods have a major computation weight, because the necessary operations during each iteration are more complex than those carried out by methods that do not utilize derivations. Nevertheless, the latter are unable to guarantee the same convergence characteristics, and normally a greater number of iterations are 12 necessary to obtain a solution (the number of iterations is usually very high and sometimes, in order to reach the optimum, said number tends towards the infinite). What distinguishes one method of descent from another, are the criteria used to choose the direction of descent, d k and the length of the step, k . The direction of research may be determined by considering a suitable approximation of the objective function regarded as function of the vector direction of optimizing, d k . In the numeric calculation code, developed ad-hoc, the function approximation is of quadratic type with continuous Hessian, symmetric and defined positive. By the Taylor formula, terminated at the 2nd order terms for sufficiently low values of the vector increment rule d k : The optimal solution is obtained when the partial derivatives are null: The optimal point in solution may be written as: The method consists in fixing, for the direction, precisely d k , and along this direction  k is moved by one step . The problem regarding the determination of the step,  k, in direction of the descent, d k , takes the name of line search (as it happens along a "line", or the direction of descent). The problem is to find a new iteration of the form Where x k denotes the current iteration, d k is the direction of descent, and  k * is a parameter of the step's length (scaling), which is the distance from the minimum. The minimum along the line formed by the descent direction is generally estimated by using a research procedure, or by a polynomial method which involves extrapolation or interpolation. In the implemented calculation code, the minimum has been estimated by using the polynomial cubic interpolation method. The basic idea of interpolation is that of utilizing an approximated third-degree polynomial representation of the step in the interval of the values of interest, which are those included between 0 and the current step. It is known that the Quasi-Newton method avoids the numeric calculation of the Hessian matrix (Newton Method), by utilising the observed behaviour of the objective function f(x) and of the gradient f(x) in order to develop information on the curvature to make an approximation of the Hessian matrix, using appropriate updating techniques. In literature, a great number of updating-inversion methods of the Hessian matrix have been developed. Among these, two methods have been chosen: BFGS and DFP. In the BFGS method (Broyden 1980, Fletcher 1970, Godfarb 1970and Shanno 1970, the matrix represent by the equation (42): while the DFP method (Davidon 1959, Fletcher andPowell 1963) is similar with the previous one by replacing the qk by sk. The research direction is determined in automatic by code -by a choice between the BFGS method and the DFP method. The Hessian matrix always remains positive, while the research direction, dk, is always in the direction of descent. The positivity of the Hessian is realized by assuring that it is initialized to be defined positive and that, thereafter, the following term will always be positive: In the polynomial cubic method proposed at each iteration the gradient and the function are calculated, the updating of the Hessian is accomplished when the new point found: meets the condition: That is, when the method is actually a method of descent. If such a condition doesn't occur, then the step  k * is reduced to form a new step  k+1 .

9
Numeric results: mono-objective optimization Shown below are the results obtained by the proposed modelling, focalized on the four events under exam and illustrated in Table 2. The map of the isosistes superimposed on the equivalent reticular structure (mesh) is represented in Figure 4a. The comparison between the values of the calculated intensities and 14 the observed intensities in the various knots of the mesh is shown in Figure 4b. The nodal displacements obtained from the optimization is shown in Figure 4c.
We have compared the values of the classes of calculated intensity in each cell of the mesh, obtained from the average of the intensity values of the four bordering knots, with the observed macro-seismic field, Figure 4d, and with the observed isosistes, Figure 4e. In addition, Figure 4f shows the trend of the intensity residues.
 The 23 July 1930 earthquake Figure 5a shows the Isosistes map, overlaid on the equivalent reticular structure (mesh); Figure 5b represents the comparison between the values of the calculated intensities and the observed intensities in the various knots of the mesh; in Figure 5c the nodal displacements obtained by optimization is shown. We have compared the values of the classes of intensity calculated within each cell of the mesh, obtained as an average of the values of the intensities of the four surrounding knots, with the observed macro-seismic field, in Figure 5d, and with the observed Isoseists, in Figure 5e. Figure 5f shows the intensity residue trends.
 The 21 August 1962 earthquake We have represented in Figure 6a the Isosistes map, overlaid on the equivalent reticular structure (mesh); in Figure 6b the comparison between the calculated intensities and the observed intensities in the various mesh knots; in Figure 6c the nodal shifting obtained by optimization. We compared the values of the intensity classes calculated in each meshcell, obtained as the average of the intensity values of the four adjacent knots, with the observed macro-seismic field, Figure 6d, and with the observed Isoseistes, Figure 6e. Figure 6f shows the trend of the intensity residues.
 The 23 November 1980 earthquake We represent in Figure 7a the Isosistes map overlaid on the equivalent reticular structure (mesh); in Figure 7b the comparison between the calculated intensities and the observed intensities in the various mesh knots; in Figure 7c the nodal shifting obtained by optimizing. We compared the values of the classes of intensities calculated in each mesh cell, obtained as the average of the intensity values of the four adjacent knots with the map of the observed Isosistes, Figure 7d, and with the observed macro-seismic field, Figure 7e. Figure 7f shows the trend of the intensity residues.
The research is mainly aimed at assessing seismic risk at territorial scale (Zuccaro and De Gregorio, 2013;Zuccaro et al., 2017), but it is also useful for the development of new construction technologies and materials in seismic engineering (Cavalieri et al., 2017;Pingue et al., 2015).

Conclusions
The purpose of this work is to develop a novel methodology able to reproduce the propagation law of seismic energy, hypothesizing the "transmission" mechanisms that preside over the distribution of seismic effects on the territory, by means of a structural optimization process with a predetermined energy distribution. In order to achieve this goal, the present paper proposes a physical mechanical prototype of seismic energy transmission in the Campania territory.
The results obtained have shown promising results; the nodal values of the intensities tend to show a good accordance, the intensities observed less than 5-6 are those that differ most from the intensities generated by the model, due probably to the bigger uncertainties tied to the low intensities.
Through the comparison of the classes of intensity, calculated according to the map of observed Isosistes and the field of macro-seismic observations, we note the close approximation of the intensities generated by the model and therefore, a satisfactory grade of approximation of the simulation technique proposed in each single event.
We also note, from the trend of the intensity residuals, that there is a close grade of approximation between the observed values and the calculated ones; in fact, there isn't a tendency to growth of the residuals and they are all considerably contained.
The results of the proposed model have particular interest, both as an aid in the defining of local laws of the attenuation of seismic energy let out at the epicentre and as an eventual analysis of the evaluation of local geological effect at macro-scale. The latter is of great actuality concerning the topics of seismic reclassification on a regional scale.
For the future developments of the present research, the authors propose a procedure based on multi-objective analysis, having the capacity of averaging, for the four selected earthquakes, the characteristics of stiffness of the soil according to a star of directions.