Behaviour Investigation of SMA-Equipped Bar Hysteretic Dampers Using Machine Learning Techniques

: Most isolators have numerous displacements due to their low stiffness and damping properties. Accordingly, the supplementary damping systems have vital roles in damping enhancement and lower the isolation system displacement. Nevertheless, in many cases, even by utilising additional dampers in isolation systems, the occurrence of residual displacement is inevitable. To address this issue, in this study, a new smart type of bar hysteretic dampers equipped with shape memory alloy (SMA) bars with recentring features, as the supplementary damper, is introduced and investigated. In this regard, 630 numerical models of SMA-equipped bar hysteretic dampers (SMA-BHDs) were constructed based on experimental samples with different lengths, numbers, and cross sections of SMA bars. Furthermore, by utilising hysteresis curves and the corresponding ideal bilinear curves, the role of geometrical and mechanical parameters in the cyclic behaviour of SMA-BHDs was examined. Due to the deﬁciency of existing analytical models, proposed previously for steel bar hysteretic dampers (SBHDs), to estimate the ﬁrst yield point displacement and post-yield stiffness ratio in SMA-BHDs accurately, new models were developed by the artiﬁcial neural network (ANN) and group method of data handling (GMDH) approaches. The results showed that, although the ANN models outperform GMDH ones, both ANN-and GMDH-based models can accurately estimate the linear and nonlinear behaviour of SMA-BHDs in pre-and post-yield parts with low errors and high accuracy and consistency.


Introduction
There are several techniques that can be used for improving seismic behaviour in structures [1][2][3][4]. Adding certain elements to the structure, such as shear walls or steel braces, would enhance the seismic strength of the structure. Due to this fact, the lateral stiffness of the structure is increased, and the increased lateral stiffness will develop the force applied to major structural components [5][6][7]. Increasing the seismic behaviour of structures can be enhanced by insulating the foundations from the supports. Despite strengthening strategy of structural elements affected by lateral loads, isolated structural elements without any strengthening could be an alternative. When the seismic demand in isolated systems is reduced, the structural performance can be enhanced without affecting Steel bar hysteretic dampers (SBHDs), as added dampers to isolation systems, could withstand strong earthquakes, but the residual displacements result in higher repair costs, lower safety levels, and a significant decline in the ability of SBHD-equipped structures to withstand aftershocks. To overcome this problem, externally bonded composites can strengthen structural components [59][60][61][62][63][64]. For the near-fault areas with greater seismic risks, a substation of steel bars with shape memory alloys (SMAs) with recentring properties may be a better option. Scientists found that SMA materials could withstand significant nonlinearities and return to their original shape after unloading. As a result of their excellent corrosion resistance, fatigue resistance, and damping capacity [65], these materials reduce maintenance costs. The use of SMAs for structural purposes has been assessed through some studies throughout the past several decades, including damping devices [66], bridge supports [67][68][69], vibrational damage detection systems [70], reinforcement systems [65,[71][72][73][74], and systems of seismic isolation [75]. On the opposite of SMA advantages in improving structural behaviour, earlier studies have demonstrated that SMAs usage poses certain drawbacks as a steel bar in SBHDs. In particular, large-diameter SMAs rebar is difficult to manufacture due to their great hardness [76]. Moreover, due to the fragility of the connecting point in SMA bars, they cannot be welded to steel sections [77]. Therefore, one of the most challenging aspects of using SMAs is the high preparation costs. Despite a decline in price over the past decade, the price for these alloys remains expensive when compared with other contemporary materials. As a result, SMAs are usually used in an optimum manner only in places where more deformation is experienced.
Despite several studies conducted on various hysteretic dampers with various geometric shapes, there are still some questions regarding their behaviour, and in particular their behaviour with steel bars. Moreover, the residual displacements in steel bar hysteretic dampers (SBHDs) increase the maintenance and repair costs after earthquakes occurrence. Utilising shape memory alloys (SMAs) as a smart material with recentring properties and less residual displacement with respect to steel bars could be a proper alternative. As the prices of SMAs are more than steel materials, their usage should be conducted in an optimum manner. To address these gaps, this study proposed a novel developed bar hysteretic damper as an added damper to an isolation system, in which the steel bars are substituted optimally by SMA bars with different mechanical properties and geometrical configurations. A comprehensive study was conducted to evaluate the role of different geometrical such as the length, number, and cross section of the SMA bars on the cyclic behaviour of SMA-equipped bar hysteretic dampers (SMA-BHDs). After evaluating existing analytical models for identifying the mechanical characteristics of SMA-BHDs, new alternative models with higher accuracy and lower error values were proposed by utilising different machine learning methods.

Shape Memory Alloys
Shape memory alloys (SMAs) are considered smart and innovative materials discovered as early as 1932 and, in 1960, were manufactured from nickel and titanium [78] to exhibit better behaviour than any other element used in memory alloys. These polymorphic materials are composed of crystalline or chemical phases. SMAs are capable of enduring large stresses without causing residual strains, and it is efficient in dissipating energy [79]. In shape memory alloys, the prevailing crystalline phase is dependent on the temperature and stress. Austenite and martensite are the two crystalline phases of SMAs. Austenite and martensite phases are stable at high and low temperatures and low and high stresses, respectively [80]. Therefore, these two phases could be transformed into each other through a heating process or by applying mechanical forces. SMAs exhibit different macroscopic behaviour during the transition from one phase to another. There are certain characteristics that distinguish SMAs from other alloys and metals, specifically, shape memory and superelasticity features. SMAs' mechanical behaviour under different levels of stress and strain and various temperatures is shown in Figure 1. A shape memory feature is present in SMAs at temperatures below the austenite-transformation-to-martensitic phase during cooling, M f , and the superelasticity feature is created when the temperature is above the finish temperature of the martensite-transformation-to-austenite phase during heating, A f [79]. At temperatures below M f and in the martensitic phase, SMAs deform under stress. By load removal, the alloy could not return to its original shape and would experience residual strains. In this case, named the one-way shape memory feature, the martensitic phase is still dominant in SMAs.
would experience residual strains. In this case, named the one-way shape memory feature, the martensitic phase is still dominant in SMAs.
In a one-way shape memory feature, just heating and cooling processes could convert the martensite and austenitic phases into one another without stress being applied. This case in which SMAs could remember their different shapes at high and low temperatures is named the two-way shape memory feature. SMAs with two-way shape memory characteristics are among the few that exist and are utilised in temperature-sensitive actuators and reversible fasteners in medical research fields [78,80]. By applying stresses, the austenitic phase of SMAs would be transformed into the martensitic phase at a higher temperature than Af and lower than Md (a higher temperature than Af, at which the alloy has an elastoplastic feature). As shown in Figure 1, unloading in this situation leads to an unstable martensitic phase, which undergoes an inverse transformation that leads SMAs to return to their original state with no residual strain. This feature is named superelasticity and is utilised in many engineering fields [79]. Researchers used this smart material in many different fields due to austenitic to martensitic transformation, resulting in shape memory and superelasticity properties. The superelasticity property of SMAs is used in this study to improve the seismic behaviour of bar hysteretic dampers utilising the stress-strain model introduced by Auricchio and Sacco [81] and also used in finite element SeismoStruct software [82]. As presented in Figure 2, the most important parameters in this model are σf EA (stresses relating to the beginning transformation of the austenitic into the martensitic phase), σs SA (stresses linked to the end of the austenitic phase into the martensitic phase transformation), σs AS (stresses related to the beginning of unloading step), σf AS (stresses relating to the ending of the unloading step), εl (strain equivalent at the unloading step) and ESMA (elastic modulus in the austenitic phase) [81]. In a one-way shape memory feature, just heating and cooling processes could convert the martensite and austenitic phases into one another without stress being applied. This case in which SMAs could remember their different shapes at high and low temperatures is named the two-way shape memory feature. SMAs with two-way shape memory characteristics are among the few that exist and are utilised in temperature-sensitive actuators and reversible fasteners in medical research fields [78,80]. By applying stresses, the austenitic phase of SMAs would be transformed into the martensitic phase at a higher temperature than A f and lower than M d (a higher temperature than A f , at which the alloy has an elastoplastic feature). As shown in Figure 1, unloading in this situation leads to an unstable martensitic phase, which undergoes an inverse transformation that leads SMAs to return to their original state with no residual strain. This feature is named superelasticity and is utilised in many engineering fields [79].
Researchers used this smart material in many different fields due to austenitic to martensitic transformation, resulting in shape memory and superelasticity properties. The superelasticity property of SMAs is used in this study to improve the seismic behaviour of bar hysteretic dampers utilising the stress-strain model introduced by Auricchio and Sacco [81] and also used in finite element SeismoStruct software [82]. As presented in Figure 2, the most important parameters in this model are σ f EA (stresses relating to the beginning transformation of the austenitic into the martensitic phase), σ s SA (stresses linked to the end of the austenitic phase into the martensitic phase transformation), σ s AS (stresses related to the beginning of unloading step), σ f AS (stresses relating to the ending of the unloading step), ε l (strain equivalent at the unloading step) and E SMA (elastic modulus in the austenitic phase) [81]. Appl. Sci. 2021, 11, x FOR PEER REVIEW 5 of 29 Figure 2. Stress-strain model of SMAs with superelasticity feature [81].

Methodology
In this paper, ANN and GMDH-NN were employed to estimate the cyclic behaviour of SMA bar hysteretic dampers. The overall search flowchart in this paper is depicted in Figure 3.

The Artificial Neural Networks (ANNs)
To address complex engineering issues, artificial neural networks (ANNs) are commonly used in machine learning techniques [83,84]. It is possible to think of a neural network as a set of processes that turn inputs into outputs [85]. The ANNs are constructed by duplicating biological neural networks, which can be trained using input and output data sets [86]. To achieve a certain degree of accuracy, the neural network must be trained in order to update the correlation weights and biases. The overall structure of ANN could be seen in Figure 4, which consists of input layers (X1 to Xn), as well as an output layer (Y) and several hidden layers [87]. It is then determined whether the model is acceptable for output estimation once it has been trained, and its weights have been tested on a separated unseen dataset [88].

Methodology
In this paper, ANN and GMDH-NN were employed to estimate the cyclic behaviour of SMA bar hysteretic dampers. The overall search flowchart in this paper is depicted in Figure 3.

Methodology
In this paper, ANN and GMDH-NN were employed to estimate the cyclic behaviour of SMA bar hysteretic dampers. The overall search flowchart in this paper is depicted in Figure 3.

The Artificial Neural Networks (ANNs)
To address complex engineering issues, artificial neural networks (ANNs) are commonly used in machine learning techniques [83,84]. It is possible to think of a neural network as a set of processes that turn inputs into outputs [85]. The ANNs are constructed by duplicating biological neural networks, which can be trained using input and output data sets [86]. To achieve a certain degree of accuracy, the neural network must be trained in order to update the correlation weights and biases. The overall structure of ANN could be seen in Figure 4, which consists of input layers (X1 to Xn), as well as an output layer (Y) and several hidden layers [87]. It is then determined whether the model is acceptable for output estimation once it has been trained, and its weights have been tested on a separated unseen dataset [88].

The Artificial Neural Networks (ANNs)
To address complex engineering issues, artificial neural networks (ANNs) are commonly used in machine learning techniques [83,84]. It is possible to think of a neural network as a set of processes that turn inputs into outputs [85]. The ANNs are constructed by duplicating biological neural networks, which can be trained using input and output data sets [86]. To achieve a certain degree of accuracy, the neural network must be trained in order to update the correlation weights and biases. The overall structure of ANN could be seen in Figure 4, which consists of input layers (X 1 to X n ), as well as an output layer (Y) and several hidden layers [87]. It is then determined whether the model is acceptable for output estimation once it has been trained, and its weights have been tested on a separated unseen dataset [88]. Appl. Sci. 2021, 11, x FOR PEER REVIEW 6 of 29

The Group Method of Data Handling (GMDH)
According to Volterra functional series, the relationship between input and output variables may be estimated differently from the ANN method. The discrete equivalent of the Volterra functional series is Kolmogorov-Gabor polynomials [89].
where A vector of weights is defined as a = (a0, a1, ..., an) for input variables X = (X1, X2, ..., Xn). Despite its ability to approximate any stationary random series of observations, the Kolmogorov-Gabor polynomial has two disadvantages. A lengthy and incomplete vector of independent variables is most common, whereas the collection of observations is limited. In addition, as the input vector increases in size, the computation time for calculating all the necessary normal equations rises. Inspired by Kolmogorov-Gabor polynomials [89], Ivakhnenko [90] created the group method of data handling (GMDH) as a heuristic self-organisation method. Ivakhnenko strived to enhance the Kolmogorov polynomial using lower orders for each pair, employing heuristic and Perceptron techniques. This approach, which is based on the Perceptron-type structure, is more accurate and allows data classification into useful and harmful categories. As a result, it requires fewer observations and decreases computation time [91]. An adaptive generalised polynomial-based function model is created by GMDH-NN, which advances in complexity over time until it achieves an ideal degree of complexity, at which it is neither too simple to generalise nor too complicated to overfit, which would result in a network, as seen in Figure 5. Using the GMDH v, the number of layers, the nodes that must be selected, and the properties of those nodes are automatically determined.

The Group Method of Data Handling (GMDH)
According to Volterra functional series, the relationship between input and output variables may be estimated differently from the ANN method. The discrete equivalent of the Volterra functional series is Kolmogorov-Gabor polynomials [89].
where A vector of weights is defined as a = (a 0 , a 1 , ..., a n ) for input variables X = (X 1 , X 2 , ..., X n ). Despite its ability to approximate any stationary random series of observations, the Kolmogorov-Gabor polynomial has two disadvantages. A lengthy and incomplete vector of independent variables is most common, whereas the collection of observations is limited. In addition, as the input vector increases in size, the computation time for calculating all the necessary normal equations rises. Inspired by Kolmogorov-Gabor polynomials [89], Ivakhnenko [90] created the group method of data handling (GMDH) as a heuristic selforganisation method. Ivakhnenko strived to enhance the Kolmogorov polynomial using lower orders for each pair, employing heuristic and Perceptron techniques. This approach, which is based on the Perceptron-type structure, is more accurate and allows data classification into useful and harmful categories. As a result, it requires fewer observations and decreases computation time [91]. An adaptive generalised polynomial-based function model is created by GMDH-NN, which advances in complexity over time until it achieves an ideal degree of complexity, at which it is neither too simple to generalise nor too complicated to overfit, which would result in a network, as seen in Figure 5. Using the GMDH v, the number of layers, the nodes that must be selected, and the properties of those nodes are automatically determined.

Experimental Specimens
In this study, as shown in Figure 6, the experimental study conducted by Golzan et al. [57] on steel bar hysteretic dampers (SBHDs) as an added damper to bridge isolator systems was considered as a reference to conduct the numerical analysis. A vertical load representing gravity loads allowed forces to be applied by two vertical jacks in their test setup. At the same time, the horizontal displacement was applied to the isolation system. The geometry of the test setup was designed to accommodate added damper parallel to the isolator system. The isolation system contained a rubber isolator designed for highway bridge piers in Quebec, Canada. The added damper included six horizontal steel bars of 38 mm diameter and 1.5 m long [57]. The energy absorption mechanism in this damper occurred when the cross section of bars reached the nonlinear range in bending, and plastic hinged induced in their two ends and mid-length. More comprehensive details regarding the experimental test setup can be found in Jahangir et al. [58] research.
self-organisation method. Ivakhnenko strived to enhance the Kolmogorov polynomial using lower orders for each pair, employing heuristic and Perceptron techniques. This approach, which is based on the Perceptron-type structure, is more accurate and allows data classification into useful and harmful categories. As a result, it requires fewer observations and decreases computation time [91]. An adaptive generalised polynomial-based function model is created by GMDH-NN, which advances in complexity over time until it achieves an ideal degree of complexity, at which it is neither too simple to generalise nor too complicated to overfit, which would result in a network, as seen in Figure 5. Using the GMDH v, the number of layers, the nodes that must be selected, and the properties of those nodes are automatically determined.

Experimental Specimens
In this study, as shown in Figure 6, the experimental study conducted by Golzan et al. [57] on steel bar hysteretic dampers (SBHDs) as an added damper to bridge isolator systems was considered as a reference to conduct the numerical analysis. A vertical load representing gravity loads allowed forces to be applied by two vertical jacks in their test setup. At the same time, the horizontal displacement was applied to the isolation system. The geometry of the test setup was designed to accommodate added damper parallel to the isolator system. The isolation system contained a rubber isolator designed for highway bridge piers in Quebec, Canada. The added damper included six horizontal steel bars of 38 mm diameter and 1.5 m long [57]. The energy absorption mechanism in this damper occurred when the cross section of bars reached the nonlinear range in bending, and plastic hinged induced in their two ends and mid-length. More comprehensive details regarding the experimental test setup can be found in Jahangir et al. [58] research. The seismic behaviour of isolation systems equipped using SBHDs was investigated by applying a predefined cyclic load pattern, as shown in Figure 7. In the first step, to study the isolator behaviour, the steel bar hysteretic dampers were not added to rubber isolators [57].  The seismic behaviour of isolation systems equipped using SBHDs was investigated by applying a predefined cyclic load pattern, as shown in Figure 7. In the first step, to study the isolator behaviour, the steel bar hysteretic dampers were not added to rubber isolators [57].

Experimental Specimens
In this study, as shown in Figure 6, the experimental study conducted by Golzan et al. [57] on steel bar hysteretic dampers (SBHDs) as an added damper to bridge isolator systems was considered as a reference to conduct the numerical analysis. A vertical load representing gravity loads allowed forces to be applied by two vertical jacks in their test setup. At the same time, the horizontal displacement was applied to the isolation system. The geometry of the test setup was designed to accommodate added damper parallel to the isolator system. The isolation system contained a rubber isolator designed for highway bridge piers in Quebec, Canada. The added damper included six horizontal steel bars of 38 mm diameter and 1.5 m long [57]. The energy absorption mechanism in this damper occurred when the cross section of bars reached the nonlinear range in bending, and plastic hinged induced in their two ends and mid-length. More comprehensive details regarding the experimental test setup can be found in Jahangir et al. [58] research. The seismic behaviour of isolation systems equipped using SBHDs was investigated by applying a predefined cyclic load pattern, as shown in Figure 7. In the first step, to study the isolator behaviour, the steel bar hysteretic dampers were not added to rubber isolators [57].

Numerical Models
In this paper, considering the experimental study conducted by Golzan et al. [57], the numerical models of isolation systems and SBHDs were generated in the SeismoStruct software to determine whether the software could achieve the expected results [82]. In these models, the same cyclic loading pattern presented in Figure 7 was applied, and using the software link element, the rubber isolator was simulated by a bilinear spring. On the other hand, since this study seeks the influence of geometrical parameters on the cyclic behaviour of bar hysteretic dampers, according to Figure 8, the steel bars were precisely modelled by SeismoStruct software.

Numerical Models
In this paper, considering the experimental study conducted by Golzan et al. [57], the numerical models of isolation systems and SBHDs were generated in the SeismoStruct software to determine whether the software could achieve the expected results [82]. In these models, the same cyclic loading pattern presented in Figure 7 was applied, and using the software link element, the rubber isolator was simulated by a bilinear spring. On the other hand, since this study seeks the influence of geometrical parameters on the cyclic behaviour of bar hysteretic dampers, according to Figure 8, the steel bars were precisely modelled by SeismoStruct software.  Figure 8 shows that steel bars are integrated into the isolation system by passing through orifices of a loading plate (Loading Plate) subjected to cyclic loads. According to the experimental test setup, each free ends of the bars are restrained by two plates (Bearing Plates) from lateral movement and rotation (from the experiment supports). Based on the mechanical characteristics of steel materials reported in Golzan et al. [57] experiments, with an elastic modulus of 200 GPa, yield stress of 250 MPa, and yield strain of 0.002, the same mechanical properties of steel materials were introduced into the software for use in the numerical model.
The responses of the numerical models under cyclic loading were obtained by evaluating support reactions and mid-length displacements of bars, and the integrated results were presented as hysteresis curves. The experimental and numerical hysteretic curves shown in Figure 9 compare the bare rubber isolator and those isolators equipped by SBHDs. Figure 9 illustrates the hysteresis curves of the bare rubber isolators, and those equipped with SBHDs obtained from numerical models are fitted with the corresponding experimental results. Consequently, using SeismoStruct software and the mechanism utilised to generate this numerical model is reliable to conduct geometrical investigations.
In this paper, to compare the cyclic behaviour of utilising SMA bars as an alternative of steel bars in hysteretic dampers, the proposed SBHD sample by Golzan et al. [57] includes 6 steel bars with 38 mm diameter and 1.5 m long, which was selected as the reference. Then, all SMA-BHD numerical models were constructed by considering equal elastic stiffness (Ke) as presented in Equation (2).  Figure 8 shows that steel bars are integrated into the isolation system by passing through orifices of a loading plate (Loading Plate) subjected to cyclic loads. According to the experimental test setup, each free ends of the bars are restrained by two plates (Bearing Plates) from lateral movement and rotation (from the experiment supports). Based on the mechanical characteristics of steel materials reported in Golzan et al. [57] experiments, with an elastic modulus of 200 GPa, yield stress of 250 MPa, and yield strain of 0.002, the same mechanical properties of steel materials were introduced into the software for use in the numerical model.
The responses of the numerical models under cyclic loading were obtained by evaluating support reactions and mid-length displacements of bars, and the integrated results were presented as hysteresis curves. The experimental and numerical hysteretic curves shown in Figure 9 compare the bare rubber isolator and those isolators equipped by SBHDs. Figure 9 illustrates the hysteresis curves of the bare rubber isolators, and those equipped with SBHDs obtained from numerical models are fitted with the corresponding experimental results. Consequently, using SeismoStruct software and the mechanism utilised to generate this numerical model is reliable to conduct geometrical investigations.
In Equation (3), N, D, and L, respectively, represent the number, diameter, and length of bars, and E indicates the elasticity modulus of the utilised materials in bar hysteretic damper. In this paper, as reported in Table 1, based on an experimental study [92], three different sets of SMA materials ( Figure 2) were considered to be compared with the corresponding steel material used in bar hysteretic dampers.  For each of Set 1 to Set 3 material, one to ten (1, 2, 3, 4, 5, 6, 7, 8, 9, and 10) numbers of SMA bars (N) with cross sections ranging from 10 to 50 mm (Φ10, Φ12, Φ14, Φ16, Φ18, Φ20, Φ22, Φ24, Φ26, Φ28, Φ30, Φ32, Φ34, Φ36, Φ38, Φ40, Φ42, Φ44, Φ46, Φ48, and Φ50 rebar) in diameter (D) were selected. Then, by considering equal Ke for SMA-BHDs and corresponding reference experimental SBHD, the equivalent length (L) of SMA bars could be obtained by Equations (2) and (3). SeismoStruct software [82] was utilised to build the numerical models to investigate the influence of different geometric parameters and mechanical features on SMA-BHDs' cyclic behaviour. In total, 630 different SMA-BHD models were constructed for which the Oo_Ll_Nn was used as a notation, to label them. Here, o refers to the diameter of the SMA rebar, followed by the calculated length and numbers of SMA bars, n and l, respectively. A schematic of hysteretic damper models of SMA bars of varying numbers (1 to 10) is presented in Figure 10. In this paper, to compare the cyclic behaviour of utilising SMA bars as an alternative of steel bars in hysteretic dampers, the proposed SBHD sample by Golzan et al. [57] includes 6 steel bars with 38 mm diameter and 1.5 m long, which was selected as the reference. Then, all SMA-BHD numerical models were constructed by considering equal elastic stiffness (K e ) as presented in Equation (2).
Based on Equation (3), the elastic stiffness of bar hysteretic dampers could be calculated by their corresponding geometrical and mechanical properties.
In Equation (3), N, D, and L, respectively, represent the number, diameter, and length of bars, and E indicates the elasticity modulus of the utilised materials in bar hysteretic damper. In this paper, as reported in Table 1, based on an experimental study [92], three different sets of SMA materials ( Figure 2) were considered to be compared with the corresponding steel material used in bar hysteretic dampers. For each of Set 1 to Set 3 material, one to ten (1, 2, 3, 4, 5, 6, 7, 8, 9, and 10) numbers of SMA bars (N) with cross sections ranging from 10 to 50 mm (Φ10, Φ12, Φ14, Φ16, Φ18, Φ20, Φ22, Φ24, Φ26, Φ28, Φ30, Φ32, Φ34, Φ36, Φ38, Φ40, Φ42, Φ44, Φ46, Φ48, and Φ50 rebar) in diameter (D) were selected. Then, by considering equal K e for SMA-BHDs and corresponding reference experimental SBHD, the equivalent length (L) of SMA bars could be obtained by Equations (2) and (3). SeismoStruct software [82] was utilised to build the numerical models to investigate the influence of different geometric parameters and mechanical features on SMA-BHDs' cyclic behaviour. In total, 630 different SMA-BHD models were constructed for which the Oo_Ll_Nn was used as a notation, to label them.
Here, o refers to the diameter of the SMA rebar, followed by the calculated length and numbers of SMA bars, n and l, respectively. A schematic of hysteretic damper models of SMA bars of varying numbers (1 to 10) is presented in Figure 10. In each model, hysteresis curves were used to evaluate support reactions (Force) versus displacement (Dis.). Figures 11 and 12 show some examples of hysteresis curves for selected SMA-BHD numerical models with mechanical features of Set 1 and different numbers on bars (with 50 mm diameter) and various cross-sectional diameters (with 10 numbers of bars), respectively. It should be mentioned that the hysteretic curves of some models resulted in just elastic part, which could not obtain any information regarding the post-yield properties. Therefore, from the total number of 630 numerical models, 389 models were selected.
Hysteresis curves derived from the numerical models were evaluated to estimate the effects of length, numbers, and the cross-sectional diameter of SMA bars on the cyclic behaviour of SMA-BHDs. As shown in Figure 13, the hysteresis curves of each curve were idealised as bilinear envelope curves. The ideal bilinear curves were derived by means of displacements and forces associated with the first yield point (Dy_Hys., Fy_Hys.), and the ultimate loading point (Du_Hys. and Fu_Hys.) of hysteresis curves. This paper identified the initial yield point in the hysteresis curve based on the first slope change, and the ultimate loading point was determined based on the final loading displacement stage Du_Hys. (0.06 m). A line with a slope equal to the hysteresis curve slope before the first yield point determined the first line of the bilinear curve. The second part of the bilinear curve was depicted based on the equal slope of the hysteresis curve slope before the ultimate point was reached. The bilinear curve had a yield point (Dy_Bi and Fy_Bi) where the lines between the first and second parts intersected. Figure 13 illustrates that because of the nonlinear behaviour of hysteresis curves, the resulting displacement and the corresponding force at the yield point in bilinear curves, Dy_Bi and Fy_Bi, respectively, are greater than the displacement and force indicative of the first yield point in the hysteresis curves (Dy_Hys. and Fy_Hys., respectively). In contrast, considering the same ultimate point, the hysteresis curves and bilinear curves exhibit equal displacement (Du_Hys. = Du_Bi. = 0.06 m) and equal force (Fu_Hys. = Fu_Bi.) at the ultimate point. In each model, hysteresis curves were used to evaluate support reactions (Force) versus displacement (Dis.). Figures 11 and 12 show some examples of hysteresis curves for selected SMA-BHD numerical models with mechanical features of Set 1 and different numbers on bars (with 50 mm diameter) and various cross-sectional diameters (with 10 numbers of bars), respectively. It should be mentioned that the hysteretic curves of some models resulted in just elastic part, which could not obtain any information regarding the postyield properties. Therefore, from the total number of 630 numerical models, 389 models were selected.
Hysteresis curves derived from the numerical models were evaluated to estimate the effects of length, numbers, and the cross-sectional diameter of SMA bars on the cyclic behaviour of SMA-BHDs. As shown in Figure 13, the hysteresis curves of each curve were idealised as bilinear envelope curves. The ideal bilinear curves were derived by means of displacements and forces associated with the first yield point (D y _Hys., F y _Hys.), and the ultimate loading point (D u _Hys. and F u _Hys.) of hysteresis curves. This paper identified the initial yield point in the hysteresis curve based on the first slope change, and the ultimate loading point was determined based on the final loading displacement stage D u _Hys. (0.06 m). A line with a slope equal to the hysteresis curve slope before the first yield point determined the first line of the bilinear curve. The second part of the bilinear curve was depicted based on the equal slope of the hysteresis curve slope before the ultimate point was reached. The bilinear curve had a yield point (D y _Bi and F y _Bi) where the lines between the first and second parts intersected. Figure 13 illustrates that because of the nonlinear behaviour of hysteresis curves, the resulting displacement and the corresponding force at the yield point in bilinear curves, D y _Bi and F y _Bi, respectively, are greater than the displacement and force indicative of the first yield point in the hysteresis curves (D y _Hys. and F y _Hys., respectively). In contrast, considering the same ultimate point, the hysteresis curves and bilinear curves exhibit equal displacement (D u _Hys. = D u _Bi. = 0.06 m) and equal force (F u _Hys. = F u _Bi.) at the ultimate point.
Apart from displacement and forces at the first yield (D y and F y ) and ultimate (D u and F u ) points, there are additional parameters in hysteresis and bilinear curves that can help in assessing the cyclic behaviour. Parameters such as elastic stiffness (K e ), post-yield stiffness (K p ), and α, which is the post-yield stiffness ratio and equals to the ratio of post-yield stiffness (K p ) to the elastic stiffness (K e ), could be calculated by Equations (4)- (6): α_Hys. = K p _Hys. K e _Hys. , α_Bi. = K p _Bi.  The amount of dissipated energy in dampers is another important parameter that should be taken into consideration. Based on the adopted bilinear curve presented in Figure 14, the process of estimating dissipated energy (E) in SMA-BHD dampers are reported in Equations (7)- (9): The amount of dissipated energy in dampers is another important parameter that should be taken into consideration. Based on the adopted bilinear curve presented in Figure 14, the process of estimating dissipated energy (E) in SMA-BHD dampers are reported in Equations (7)-(9): F y _Bi. = K e _Bi. × D y _Bi. (7) where F y _Bi. and F u _Bi. represent the yield and ultimate force, and D y _Bi. and D u _Bi. show their corresponding displacements in the idealised bilinear curves, respectively. Moreover, K e _Bi. and K p _Bi. are, respectively, elastic and post-yield stiffness. Q_Bi. shows the yield force corresponding to the unloading step of the bilinear curve at the ultimate displacement and can be obtained is the stress relating to the beginning transformation of the austenitic into the martensitic phase, and σ f AS is the stress relating to the ending of the unloading step, as presented in Figure 2. Both σ f EA and σ f AS are the mechanical properties of SMA materials ( Figure 2) and, therefore, are known parameters. Moreover, F y _Bi. is equal to K e _Bi. × D y _Bi. as presented in Equation (7), where K e _Bi. (elastic stiffness of bilinear curve) is equal to K e _Hys. (elastic stuffiness of hysteretic curve), as stated in Equation (4), and D y _Bi. should be determined as an unknown parameter. Furthermore, as presented in Equation (6), K p _Bi. is equal to α × K e _Bi., where K e _Bi. is a known parameter (is equal to K e _Hys.), and α should be determined as an unknown parameter. As the ultimate point in bilinear curves is considered to be the same as hysteresis curves, their corresponding displacements are equal (D u _Bi. = D u _Hys.), and therefore, D u _Bi. is a known parameter. As a result, as shown in Figure 15, the only unknown parameters in the process of Q_Bi. determination are D y _Bi. and α which should be taken to consideration.
Instead of utilising hysteretic curves parameters, there are already existing analytical equations which used geometrical and mechanical properties of hysteretic bar dampers to identify their cyclic behaviour. The yield displacement (D y ) and an elastic stiffness (K e ) of SMA-BHDs could be described as follows [58]: In Equations (10) and (11), f y and E are, respectively, the yield stress and elasticity modulus of materials, and as reported before, N, D, and L are the numbers, diameter, and calculated length of SMA bars, respectively. In most analytical cases, the ideal bilinear curve, obtained from an experimental hysteretic curve, is considered as the analysis reference. As could be inferred from Figure 13, as a result of the high nonlinear behaviour of SMA-BHDs, although the slope of the elastic parts is equal (K e ), the first yield point of the bilinear curve is quite different from the corresponding experimental one. Therefore, the proposed analytical equation for estimating the yield point displacement (D y ) should be modified by a modification factor (β) as follows: _ .
where Fy_Bi. and Fu_Bi. represent the yield and ultimate force, and Dy_Bi. and Du_Bi. show their corresponding displacements in the idealised bilinear curves, respectively. Moreover, Ke_Bi. and Kp_Bi. are, respectively, elastic and post-yield stiffness. Q_Bi. shows the yield force corresponding to the unloading step of the bilinear curve at the ultimate displacement and can be obtained as (σf AS / σf EA ) × Fy_Bi. + Kp_Bi. × (Du_Bi. − Dy_Bi.), where σf EA is the stress relating to the beginning transformation of the austenitic into the martensitic phase, and σf AS is the stress relating to the ending of the unloading step, as presented in Figure 2. Both σf EA and σf AS are the mechanical properties of SMA materials ( Figure 2) and, therefore, are known parameters. Moreover, Fy_Bi. is equal to Ke_Bi. × Dy_Bi. as presented in Equation (7), where Ke_Bi. (elastic stiffness of bilinear curve) is equal to Ke_Hys.
(elastic stuffiness of hysteretic curve), as stated in Equation (4), and Dy_Bi. should be determined as an unknown parameter. Furthermore, as presented in Equation (6) _ .
where Fy_Bi. and Fu_Bi. represent the yield and ultimate force, and Dy_Bi. and Du_Bi. show their corresponding displacements in the idealised bilinear curves, respectively. Moreover, Ke_Bi. and Kp_Bi. are, respectively, elastic and post-yield stiffness. Q_Bi. shows the yield force corresponding to the unloading step of the bilinear curve at the ultimate displacement and can be obtained as (σf AS / σf EA ) × Fy_Bi. + Kp_Bi. × (Du_Bi. − Dy_Bi.), where σf EA is the stress relating to the beginning transformation of the austenitic into the martensitic phase, and σf AS is the stress relating to the ending of the unloading step, as presented in Figure 2. Both σf EA and σf AS are the mechanical properties of SMA materials ( Figure 2) and, therefore, are known parameters. Moreover, Fy_Bi. is equal to Ke_Bi. × Dy_Bi. as presented in Equation (7), where Ke_Bi. (elastic stiffness of bilinear curve) is equal to Ke_Hys.
(elastic stuffiness of hysteretic curve), as stated in Equation (4), and Dy_Bi. should be determined as an unknown parameter. Furthermore, as presented in Equation (6)  be determined as an unknown parameter. As the ultimate point in bilinear curves is considered to be the same as hysteresis curves, their corresponding displacements are equal (Du_Bi. = Du_Hys.), and therefore, Du_Bi. is a known parameter. As a result, as shown in Figure 15, the only unknown parameters in the process of Q_Bi. determination are Dy_Bi. and α which should be taken to consideration. where D y _Mod. is the modified version of the existing first yield point displacement (D y ) analytical equation. On the other hand, both D y and K e parameters are related to the elastic part of SMA-BHDs behaviour, and no analytical equation has been proposed to identify the post-yield behaviour, specifically the stiffness ratio α. Therefore, the purpose of this paper was to complete the analytical equations and assess their ability to predict the cyclic behaviour of SMA-BHDs in the nonlinear section via machine learning techniques such as artificial neural network (ANN) and group method of data handling (GMDH).

Numerical Database
In order to conduct machine learning approaches, the obtained database from 389 numerical models (selected from 630 numerical models based on their proper hysteretic curves including post-yield portion), with geometrical and mechanical properties of SMA-BHDs, was collected. The statistical properties of the obtained database in this paper are presented in Table 2, in which, D, L, N, E, and F y , as diameter, length and numbers, elasticity modulus, and yield stress of SMA bars, were considered as the inputs, and the β_Bi and α_Bi, as yield displacement modification factor and the post-yield stiffness ratio, resulted from ideal bilinear curves considered as outputs. As can be seen, the collected numerical data cover a wide range of each contributing parameter, and hence, these parameters can be proper inputs for estimating the exact values of modification factor of first yield point displacement (β) and the post-yield stiffness ratio (α) in SMA-BHDs.

Data Pre-Processing
With the following linear equation, the collected data from numerical models were normalised and scaled across the same ranges in order to ensure stability and convergence of weight and biases in the process of ANN and GMDH-NN development: Here, X represents the input or output variable, while X min and X max represent its minimum and maximum values, respectively. The ANN and GMDH-NN structure for the first yield point displacement modification factor (β) and also post-yield stiffness ratio (α) prediction is formed by considering inputs including the diameter (D), length (L), and numbers (N) of bars as long as elastic modulus (E) and yield strength (F y ) of SMA bars.

Proposed ANN Models
According to the proposed method by Shahin et al. [93], in ANN models, the obtained database (containing 389 data) is divided into two training (around 75% of all datasets equal to 292 data) and testing (around 25% of all datasets equal to 97 data) classes. Moreover, to assure the ANN and GMDH-NN models efficiency, the train set (around 75% of dataset equal to 292 data) itself is classified into three subsets, including train, test, and validation, with around 80% (232 data), 10% (30 data) and 10% (30 data) ratios, respectively. In order to determine the structure of the ANN models, the feed-forward back propagation model, which has a single input layer and one hidden layer with the tan-sigmoid activation function (f n 2 ), as shown in Equation (14), and an output layer having linear activation function (f n 1 , expressed as ax + b where a and b are constants) is used.
MATLAB was used to train the connecting weights of network neurons using feedforward backpropagation and the Levenberg-Marquardt method [94]. The literature has supported the usage of one hidden layer to tackle various nonlinear problems [85,95]. The optimum artificial neural network was created using the trial-and-error approach. In this study, the number of neurons in the hidden layer was assumed to be between 3 and 25. The optimal configuration was determined by using conventional statistical error and performance metrics, such as the correlation coefficient (R), the mean square error (MSE), and the mean absolute percent error (MAPE). The results of the ANN models in estimating the first yield point displacement modification factor (β) and also post-yield stiffness ratio (α) are shown in Tables 3 and 4, respectively. Figure 16 shows the R and MSE values for different neurons in the hidden layer of ANN models. As can be seen in Table 3 and Figure 16, the structure of the best ANN model for estimating β contains six neurons in the hidden layer with R values of 0.9972, 0.9958, and 0.9916 in training, testing, and validation data sets and considerably small MSE values of 0.00012 and 0.00047, and MAPE values of 2.60% and 3.31%, in training and testing data sets, respectively. On the other hand, as reported in Table 4 and Figure 16, to estimate the α values, selecting 10 neurons in the hidden layer is the most proper ANN model with the R values of 0.9986, 0.9966, and 0.9972 in training, testing, and validation data sets and MSE values of 0.00007 and 0.00008, and MAPE values of 2.34% and 2.43%, in training and testing data sets, respectively. It should also be noted that training data division in the network includes one step of verification, which uses 80%, 10%, and 10% for training, testing, and validation, respectively. However, in order to provide a more reliable and efficient model, testing data division, including 97 data (25% out of 389 data), are entirely unseen for the models and provide a double-check step for the models.  As can be seen in Table 3 and Figure 16, the structure of the best ANN model for estimating β contains six neurons in the hidden layer with R values of 0.9972, 0.9958, and 0.9916 in training, testing, and validation data sets and considerably small MSE values of 0.00012 and 0.00047, and MAPE values of 2.60% and 3.31%, in training and testing data sets, respectively. On the other hand, as reported in Table 4 and Figure 16, to estimate the α values, selecting 10 neurons in the hidden layer is the most proper ANN model with the R values of 0.9986, 0.9966, and 0.9972 in training, testing, and validation data sets and MSE values of 0.00007 and 0.00008, and MAPE values of 2.34% and 2.43%, in training and testing data sets, respectively. It should also be noted that training data division in the network includes one step of verification, which uses 80%, 10%, and 10% for training, testing, and validation, respectively. However, in order to provide a more reliable and efficient model, testing data division, including 97 data (25% out of 389 data), are entirely unseen Using the weights and bias acquired from the proposed ANN models, the following formula was developed to establish the mathematical relationship between input (D, L, N, E, and F y ) and outputs (β and α) variables: where Y denotes the values corresponding to output parameters (β and α), f n 1 and f n 2 are the transfer functions, h indicates the number of neurons in the hidden layer (in this paper equals to 6 and 10 for ANN(β) and ANN(α), respectively), X i is the input values of the network (D, L, N, E, and F y ), m is the number of the input variables (equals to 5), W ik indicates the link weight between the ith input layer and k th neuron in the hidden layer, w k is the link weight between kth neuron in the hidden layer and the independent output layer, b hk represents a bias in the kth neuron of the hidden layer, and b 0 is the bias value in the output layer. Therefore, to ensure other researchers may profit from these results, network weights and bias values were supplied in this study. The estimated β and α using ANN models, respectively, can be achieved by Equation (15), and the results are presented in Table 5.
The β and α estimations based on ANN models are expressed as following with respect to Equation (15), and the weights are presented in Table 5: where A 1 to A 6 and C 1 to C 10 are the response of the hidden neurons which feed the network output and are calculated as A and Carray elements (Equation (18)) using Equations (19) and (20).
where I is the input array, W β and B β, and W α and B α are, respectively, corresponding β and α weight and bias matrices which are presented in Equations (21)- (23).

Proposed GMDH Models
Around 25% of the database (97 data out of 389 data) was arbitrarily put aside in GMDH-NN techniques, analogous to ANN models, in order to suggest more trustworthy closed-form equations. After many tests, using the GMDH-NN methods, Equations (24) and (25) were proposed to estimate the β and α parameters, respectively. For the matter of brevity, further details of the GMDH-NN method are not included here; interested readers are encouraged to refer to the relevant published papers [85,[96][97][98][99].

ML Proposed Models Performances
To evaluate the performance of the ANN and GMDH-NN models, Equations (26)-(31) were used, which are typical criteria for measuring error and model performance [100], Appl. Sci. 2021, 11, 10057 20 of 28 including the correlation coefficient (R) and the coefficient of determination (R 2 ), mean square error (MSE), root-mean-square error (RMSE), mean absolute error (MAE) and mean absolute percentage error (MAPE). Detailed results are presented in Table 6 for each developed model. In Equations (26)- (31), A i indicates the analysed value, and F i represents the estimated value, n is the number of the considered data, A is the mean analysed values, and F is the mean estimated values.
As shown in Table 6, R and R 2 for the ANN models (ANN (β) and ANN (α)) in both training and testing steps are more than the GMDH-NN models (GMDH-NN (β) and GMDH-NN (α)). In the same manner, the evaluated error criteria demonstrate that the ANN models became more accurate than corresponding GMDH-NN models. In  16, respectively, in training and testing stages. The comparison between the measured and predicted values obtained from the ANN and the GMDH models in the training and testing phases for estimating α and β are presented in Figures 17 and 18, respectively. Figures 17 and 18 indicate that the ANN models produce more accurate and wellfitted ideal fit lines than the GMDH-NN models. Additionally, the error values of β and α proposed models for training and test data stages are presented in Figures 19 and 20, respectively. An acceptable very low error for both proposed ANN and GMDH-NN models could be inferred from Figures 19 and 20.

Sensitivity Analysis
Sensitivity analysis was performed to evaluate the contribution of each input parameter to the ANN and GMDH developed β and α models. Results for each input variable were obtained by assuming that each input variable varies between its minimum and maximum values and that other parameters are maintained at their mean values. Therefore, it was analysed how each parameter affects the correlation coefficient (R) and the root-meansquare error (RMSE). By substituting all or none of the variables with their mean values, the impact on R and RMSE could be considered to be zero and 100%, respectively. The R and RMSE impact may be defined as a percentage value using the following equation [101]:     Figures 17 and 18 indicate that the ANN models produce more accurate and wellfitted ideal fit lines than the GMDH-NN models. Additionally, the error values of β and α proposed models for training and test data stages are presented in Figures 19 and 20, respectively. An acceptable very low error for both proposed ANN and GMDH-NN models could be inferred from Figures 19 and 20.

Sensitivity Analysis
Sensitivity analysis was performed to evaluate the contribution of each input parameter to the ANN and GMDH developed β and α models. Results for each input variable In this equation, Z Var is the R and RMSE values for the considered variable, Z Ori is the effect of zero on the R and RMSE, and Z all is equal to the R and RMSE value of the proposed model when all the variables are substituted with their mean values. As shown in Figures 21 and 22, for both β and α, the calculated SMA bars length (L) is the most influential input variable in both ANN and GMDH-NN models. On the other hand, the number of SMA bars (N) has the least effect in both ANN and GMDH-NN models. were obtained by assuming that each input variable varies between its minimum and maximum values and that other parameters are maintained at their mean values. Therefore, it was analysed how each parameter affects the correlation coefficient (R) and the root-mean-square error (RMSE). By substituting all or none of the variables with their mean values, the impact on R and RMSE could be considered to be zero and 100%, respectively. The R and RMSE impact may be defined as a percentage value using the following equation [101]: In this equation, ZVar is the R and RMSE values for the considered variable, ZOri is the effect of zero on the R and RMSE, and Zall is equal to the R and RMSE value of the proposed model when all the variables are substituted with their mean values. As shown in Figures  21 and 22, for both β and α, the calculated SMA bars length (L) is the most influential input variable in both ANN and GMDH-NN models. On the other hand, the number of SMA bars (N) has the least effect in both ANN and GMDH-NN models.

Computational Costs
Various machine learning (ML) techniques could be utilised to estimate the first yield point displacement and post-yield stiffness ratio in SMA-BHDs. However, in this paper, new models were developed by the artificial neural network (ANN) and group method of data handling (GMDH) techniques. The previous sections showed that proposed ANN  were obtained by assuming that each input variable varies between its minimum and maximum values and that other parameters are maintained at their mean values. Therefore, it was analysed how each parameter affects the correlation coefficient (R) and the root-mean-square error (RMSE). By substituting all or none of the variables with their mean values, the impact on R and RMSE could be considered to be zero and 100%, respectively. The R and RMSE impact may be defined as a percentage value using the following equation [101]: In this equation, ZVar is the R and RMSE values for the considered variable, ZOri is the effect of zero on the R and RMSE, and Zall is equal to the R and RMSE value of the proposed model when all the variables are substituted with their mean values. As shown in Figures  21 and 22, for both β and α, the calculated SMA bars length (L) is the most influential input variable in both ANN and GMDH-NN models. On the other hand, the number of SMA bars (N) has the least effect in both ANN and GMDH-NN models.

Computational Costs
Various machine learning (ML) techniques could be utilised to estimate the first yield point displacement and post-yield stiffness ratio in SMA-BHDs. However, in this paper, new models were developed by the artificial neural network (ANN) and group method

Computational Costs
Various machine learning (ML) techniques could be utilised to estimate the first yield point displacement and post-yield stiffness ratio in SMA-BHDs. However, in this paper, new models were developed by the artificial neural network (ANN) and group method of data handling (GMDH) techniques. The previous sections showed that proposed ANN models outperform GMDH models with lower errors and higher with low errors and high accuracy. On the other hand, considering the computational costs related to ML models' complexity, it can be mentioned that the proposed GMDH models, with less complex closed-form equations, perform better in comparison with ANN models.

Conclusions
A novel bar hysteretic dampers equipped with shape memory alloy (SMA) bars, named SMA-BHDs, as an added damper to isolation systems, was studied in this paper. In order to predict the cyclic behaviour of these dampers, 630 numerical models including different geometrical and mechanical properties were constructed in SeismoStruct software. The obtained hysteretic curves from numerical models were idealised by bilinear curves to achieve pre-and post-yield parameters such as the first yield point displacement (D y ), elastic stiffness (K e ), and post-yield stiffness ration (α). The analyses show that the existing analytical equation for determining D y is not consistent with obtained results from hysteretic curves and needs modification. Moreover, there is no analytical equation for estimating the α parameter. Therefore, two machine learning approaches, named artificial neural network (ANN) and group method of data handling integrated by a neural network (GMDH-NN), were utilised to estimate the D y modification factor (β) and propose an equation for α parameter for the first time. An overview of the results of this paper is presented as follows: 1.
As expected, by substituting SMA bars instead of steel bar in bar hysteretic dampers, no residual displacement could be seen in hysteretic curves, which shows the excellent performance of SMA-BHDs as added dampers to isolation systems.

2.
Considering the ANN models with one hidden layer and varying neuron numbers between 3 to 25, the neural networks with 10 and 6 hidden neurons were selected as the optimised network structure for β and α parameters, respectively. In the GMDH-NN models, similarly to ANN models, around 25% of all databases (97 data from 389 data) were randomly set aside for the test stage and considered unseen data. The results show that the proposed ANN models with higher R 2 and lower error values in both the training and testing stages outperform the proposed GMDH-NN models. However, compared with the ANN model, GMDH-NN models present more user-friendly and easy-to-interpret closed-form equations.

4.
The sensitivity analysis of the input parameters in the developed ANN and GMDH-NN models for estimating both β and α parameters showed that the calculated SMA bars length (L) variable with higher R and RMSE values is the most influential input variable. Furthermore, the number of SMA bars (N) with lower impact values on the R and RMSE has the least effect.