Online Stability Analysis for Real-Time Hybrid Simulation Testing

Real-time hybrid simulation (RTHS) is an experimental technique where a critical element of a structural system is tested in the laboratory while the rest is represented through numerical simulations. A challenging aspect of this technique is the correct application of boundary conditions on the experimental substructure using actuators and sensors. The inherent dynamics of an actuator and its interaction with the physical specimen causes a time delay between commanded and measured displacements. It has been shown that delay in RTHS affects the accuracy of an experiment and even can cause instability. Therefore, to avoid stability problems, a proper partitioning choice and an appropriate compensation method for actuator dynamics should be considered. However, there will always be uncertainty in the experimental structure's behavior, so it is essential to check the system's stability during the test execution. In this paper, a stability analysis using energy methods is performed to develop an online stability indicator for the RTHS test. This indicator's goal is to detect stability problems before it can cause excessive displacements in the system, thus avoiding damage in the physical specimen or the laboratory equipment. The effectiveness of the proposed online stability indicator is demonstrated through numerical simulations taking into account the virtual RTHS benchmark problem with different compensation strategies. The proposed indicator is an excellent tool to monitor the RTHS test, improving the reliability of the experimental test while maintaining the safety of the laboratory resources.


INTRODUCTION
Laboratory tests are essential to study structural systems' behavior, calibrate mathematical models, and develop robust design methods to achieve economic and safe structures. In the case of seismic tests, the most realistic experimental technique is the shake table test. However, these experiments are very challenging for full-scale systems, not only for the laboratory requirements in terms of equipment and capacity but also for the manufacturing costs and time required to build each specimen to be tested. Alternatively, another technique called real-time hybrid simulation (RTHS) has proven to be a cost-effective and reliable approach to conducting seismic performance assessment.
Real-time hybrid simulation (RTHS) is an experimental technique for structural testing. A critical component is studied in the laboratory, while a numerical model represents the rest of the structural system (Nakashima et al., 1992). Representing the numerical part of the structure reduces the costs of each experiment considerably. The technique involves solving the equation of motion (EOM) of the numerical substructure with an integration algorithm, then imposing the calculated displacements over the experimental substructure using a transfer system (i.e., actuators). The experimental restoring forces are measured and incorporated into the equation of motion to calculate the displacement at the following time step (McCrum and Williams, 2016).
The synchronization in the boundary between numerical and experimental substructures is one of the biggest challenges of RTHS. The dynamic properties of the transfer system produce amplitude and delay errors between commanded and measured displacements. This effect not only depends on actuator properties but also on the interaction with the physical specimen (Dyke et al., 1995). The delay errors are the most harmful for RTHS because it has been proved that delay is equivalent to introducing negative damping to the hybrid system (Horiuchi et al., 1999), which affects the accuracy of the experiment and even can cause instability in the experimental setup. Several authors have studied the effect of delay in RTHS. Wallace et al. (2005a) performed a stability analysis with analytical delay differential equations. Mercan and Ricles (2008) considered pseudo-delay techniques for determining the size of delay to initiate instability. Maghareh et al. (2017) established a predictive stability indicator, demonstrating that some partitioning choices are more sensitive to delay than others. Gao and You (2019) proposed generalized EOM for the hybrid system to predict the RTHS stability limit. All of these methodologies are useful to evaluate different partitioning choices, avoiding cases susceptible to instability. However, three significant limitations are identified: (i) methodologies above were developed for linear systems; (ii) they require knowledge of the dynamic properties of experimental substructure to determine critical delay values; and (iii) these tools are used to assess the stability of the test offline before a simulation, so they cannot be implemented to check the stability online during the RTHS testing.
On the other hand, to minimize the synchronization errors and improve the accuracy of the test, several compensation methods have been proposed in the literature; for example, polynomial extrapolation (Horiuchi et al., 1999) or modelbased compensation (Carrion and Spencer, 2007;Gao et al., 2013). More sophisticated techniques are known as adaptive compensation (Wallace et al., 2005b;Bonnet et al., 2007), where the control parameters are adjusted during the RTHS to improve synchronization. Meanwhile, methods such as Tao and Mercan (2019) or Xu W. et al. (2019) are based on a frequency domain analysis to adjust parameters of a first-order transfer function. Other methods like Chae et al. (2013) or Palacio-Betancur and Gutierrez Soto (2019) estimate the plant through Taylor series expansion and adjust the parameters in the time domain. Different approaches consider polynomial extrapolation, such as Wang et al. (2019) or Xu D. et al. (2019). Adaptive modelbased compensation (Chen et al., 2015) found in this paper consists of an estimate of the plant in frequency-domain; then, compensation is implemented in time-domain using numeric derivatives of the commanded signal and adaptation based on gradient.
Even with an appropriate compensation method, it is possible to observe synchronization errors due to uncertainty or nonlinear behavior in the experimental substructure; therefore, it is necessary to monitor the performance of an RTHS test for safety purposes in a quantitative manner. Additionally, it is worth mentioning that a response of the reference structure is not always available in RTHS, so it is not trivial to evaluate the experimental response's reliability during the test.
Some studies are available in the literature where online indicators are proposed for hybrid simulation tests. Guo et al. (2014) proposed the Frequency Evaluation Index (FEI) to evaluate the tracking performance of an actuator in the frequency domain. Besides, FEI can be implemented online using a moving window Fast Fourier Transform (FFT) (Xu W. et al., 2019). Although FEI could be used to measure synchronization errors quantitatively, it does not measure the effect of this delay on the hybrid system's stability. Meanwhile, Ahmadizadeh and Mosqueda (2009) proposed the Energy Error Indicator (EEI), which can be implemented online to measure the accumulated errors in RTHS. Although experimental errors can trigger instability in a hybrid system, there is no way to determine which level of EEI is associated with instability. Thus, the EEI is not efficient if the goal is to detect instability rather than errors.
This study presents a hybrid system stability analysis using energy methods to develop a stability indicator that can be evaluated online during the RTHS test. The proposed indicator aims to detect instability before it can cause excessive displacements in the system that can damage the physical specimen or the laboratory equipment. The structure of this study is the following: Section 2 presents the methodology of the proposed indicator. Section 3 describes the virtual RTHS benchmark problem (Silva et al., 2020), which is a well-recognized and representative testbed to simulate realistic RTHS experiments by the hybrid simulation community . This benchmark problem was considered as an example of the implementation and validation of the proposed stability indicator. Adaptive model-based compensation (Chen et al., 2015), with some modifications (Fermandois et al., 2020), is employed with different control parameters to show the effectiveness of the proposed indicator in unfavorable scenarios. Section 4 provides the results of different simulations, demonstrating the capacity of the proposed indicator to detect instability early. Finally, section 5 discusses the principal findings and final remarks of this study.

Substructuring in RTHS
Consider a reference structure subjected to a ground acceleration as shown in Figure 1. The reference structure is separated into a numerical substructure and a experimental substructure. The boundary conditions of the common degrees of freedom at the interface between substructures are imposed by actuators over the experimental substructure. In terms of equations, the substructuring is described as follow. The equation of motion (EOM) of the reference structure (subscript r) is expressed as: where M r , K r , and C r are the mass, stiffness, and damping matrices, respectively. x(t),ẋ(t), andẍ(t) are the displacement, velocity, and accelerations vectors, respectively, all measured relative to the ground motion.ü g (t) is the ground acceleration, and Ŵ is the seismic influence vector. The reference structure is separated into a numerical (subscript n) and a experimental (subscript e) substructures, as shown in Equation (2): Then, the EOM of the numerical substructure can be rewritten as shown in Equation (3), including the earthquake equivalent forces of the experimental substructure: where x n (t),ẋ n (t), andẍ n (t) corresponds to the displacement, velocity, and acceleration of the numerical substructure, all relative to the ground motion. Meanwhile, F e corresponds to the feedback forces from the experimental substructure described in Equation (4): where x e (t),ẋ e (t), andẍ e (t) corresponds to the experimental displacement, velocity, and acceleration of the experimental substructure, and in the ideal case, the actuators imposes x e = x n . However, in an RTHS test, the experimental displacements are not necessarily equal to the numerical displacements due to synchronization errors introduced by the transfer system and interaction with the experimental substructure (Dyke et al., 1995).

Energy Balance
Consider that the numerical substructure is an isolated system subjected to two external forces: (i) earthquake equivalent forces; and (ii) experimental feedback forces. Assuming that an appropriate numerical integration method is implemented, so that the calculated numerical substructure response is numerically stable (Chen et al., 2009;Bas and Moustafa, 2020), the equation of motion of the numerical substructure can be expressed as an energy balance, allowing to focus the stability analysis on the interaction between hybrid substructures due to experimental errors. If numerical substructure matrices M n , C n , K n are symmetric;the energy balance can be obtained by taking the inner product with an infinitesimal numerical displacement trajectory, dx n , and then integrating both sides of Equation (3) over the displacement trajectory: ẍ T n M n dx n + ẋ T n C n dx n + x T n K n dx n Equation (3) can be later expressed as a scalar equation as shown in Equation (6): where E k , E d , and E s are the kinetic, dissipated and strain energy of numerical substructure, respectively, as described in Equation (7): E k = ẍ T n M n dx n ; E d = ẋ T n C n dx n ; E s = x T n K n dx n (7) meanwhile, W I is the work done by the earthquake forces (namely as input work), and W F is the work done by the experimental forces (namely as feedback work): It is worth mentioning that the energy balance from Equation (6) must be fulfilled during a test even if the experimental force is delayed, and the hybrid system is unstable. In this sense, the feedback force is just another external excitation supplied over the numerical substructure.

Proposed Stability Indicator
When a hybrid system becomes unstable, the mechanical energy of both substructures will grow exponentially even if the earthquake signal goes zero. This issue could happen if the experimental substructure introduces additional energy to the numerical substructure. Horiuchi et al. (1999) demonstrated that delay causes negative damping in the hybrid system, which is equivalent to adding external energy to the hybrid system. It is desirable to detect this problem before the mechanical energy grows up at levels that can damage the experimental equipment.
During an RTHS test, it is expected that the mechanical energy of the numerical substructure increases, because the earthquake forces introduces energy to the system, but the mechanical energy should not be greater than the input work. Reordering Equation (6) and grouping E k and E s in mechanical energy E mec , Equation (9) is obtained: If the feedback work, W F , gets bigger than the dissipated energy E d , the right side of Equation (9) becomes negative. Due to the energy balance, the left side becomes negative and indicates that mechanical energy gets bigger than the input energy. Thus, when the hybrid system's inherent dissipated energy is insufficient to counteract the added energy from the feedback force, the system will become unstable. Therefore, the proposed analysis consists in comparing the feedback work with the dissipated energy through an indicator described in Equation (10): where SW is called the Stability Warning Indicator, and H(·) is the Heaviside function used to compute SW only when W F is positive. Notice that the dissipated energy is always strictly positive (i.e., E d > 0); hence, E d − W F could be negative only if W F is positive. Meanwhile, C SW is a constant chosen to prevent large values of SW when the test starts (i.e., when E d is close to zero). Notice that E d = 0 at t = 0 (i.e., before the earthquake starts), so C SW should be a constant selected that is bigger enough to make the denominator to be non-zero in the SW and to prevent large values of SW if W F is dominated by noise in the force measurements. Moreover, C SW should be small enough, so the SW denominator would be dominated by the E d term during and after the earthquake. As a recommendation, both conditions could be satisfied if C SW is selected as two magnitude orders below the maximum input work over the structural system. It is worth mentioning that this input work can not be calculated before the test since the complete hybrid system's response is unknown. However, an offline estimation of the input work magnitude order is enough to select a proper C SW .
The explanation of the SW values is provided herein. When SW = 0%, due to negative values of W F , there is no risk of instability because the experimental substructure is removing energy from the numerical substructure. Next, if 0% < SW < 100%, the experimental substructure is adding energy; but, the damping of the numerical substructure is enough to counter this effect. Finally, if SW ≥ 100%, then the numerical substructure's mechanical energy overcomes the input work, causing instability.
The main difference between EEI and SW indicators is that the former evaluates simulation accuracy when affected by synchronization errors, and the latter is focused on assessing the stability of the hybrid system. For example, if a hybrid system presents a significant delay in the transfer system, the experimental errors will cause the experimental substructure to add energy to the numerical substructure. Nevertheless, if the numerical substructure has enough damping to maintain the system stable, the SW will be <100%; meanwhile, EEI will show large values due to the synchronization errors produced by the delay.

Benchmark Problem
The virtual RTHS benchmark problem from Silva et al. (2020) is considered as a practical example to prove the capacity of the proposed indicator to detect instability. This problem consists of a three-story moment frame with three lateral degrees of freedom as the reference structure (see Figure 2). This structure is separated into a 3DOF numerical substructure and an SDOF experimental substructure.
The reference structure properties are: mass per floor m = 1, 000 kg; natural frequencies f = [3.61, 16.00, 38.09] Hz; and modal damping ζ = 3% for each mode. These properties corresponds to Case IV from Silva et al. (2020). The reference structure is partitioned as described in Equation (11): where M e = diag(m e , 0, 0); C e = diag(c e , 0, 0); and K e = diag(k e , 0, 0). The structural properties of the experimental substructure are m e = 29.1 kg; c e = 114.6 Nsec/m; and k e = 1.19 × 10 6 N/m. The reference structure and the numerical substructure are implemented in the state-space form using Simulink for direct integration. The solver utilized is a 4th Order Runge-Kutta with a fixed time step of t = 1/4, 096 s. The reference structure model is presented in Equations (12) and (13) with ground acceleration as input and reference displacements as output.
x r = I 0 The numerical substructure model is presented in Equations (14) and (15) with ground acceleration and experimental force as inputs and numerical displacements, velocities, and accelerations as outputs.
where γ = [1; 0; 0] T and f e is experimental force calculated as shown in Equation (16):  (17): This integration method is implemented in Simulink, as shown in the block diagram of Figure 3, where Discrete Finite Impulse Response (FIR) Filters and Product blocks are employed to obtain the area of each trapezoid, and a Discrete-time Integrator block is required to obtain the accumulation over time. In the case of W F the input force is P = f e and the input displacement is x = x (1) n (scalar signals); meanwhile, to evaluate E d , the input forces are P = C nẋn and the displacements are x = x n (vector signals). Additionally, to compare the proposed stability indicator, the EEI (Ahmadizadeh and Mosqueda, 2009) is calculated for each simulation. The EEI can be obtained as follow: where E * k , E * d are the kinetic and dissipated energy calculated with the first derivative of displacement instead of velocity of numerical substructure. E max E is a nonzero elastic energy term utilized to normalize EEI when W I is close to zero at the beginning of the test. It is important to declare that noncompensated delay implies added energy to the system, resulting in negative values of EEI. W exp is the work done by experimental force over the measured displacement as shown in Equation (19): In this study, the reference structure is subjected to the Kobe earthquake with a peak ground acceleration (PGA) scaled to 40%, and the maximum input work from the earthquake is around 100 Joules. This value can be computed from the integration of the reference structure. In a real RTHS experiment, the magnitude order of the input energy can be estimated from the ground motion and global properties of the structure, such as total mass and approximated damping and stiffness. Therefore, the denominator constants are selected considerably smaller, such as 1% of the maximum input work: This value is bigger enough to prevent large values for each indicator at the beginning of the earthquake. Also, these constants will have a negligible effect on each indicator during and after the earthquake, since C SW will be relatively small compared with E d in the denominator of Equation (10), just like E max E will be smaller than W I in the denominator of Equation (18).

Critical Delay Estimation
Before conducting the RTHS simulations with different delay values, an estimation of the critical delay is obtained with the method proposed by Gao and You (2019). This method consists in constructing the EOM of the hybrid system with some assumptions synchronization errors (delay and amplitude errors). The EOM of the hybrid system has the following equivalent mass, damping, and stiffness: C eq = C n −1 + M e ω 2 δ −1 + C e −1 − K e δ −1 (22) where is the modal matrix of the reference system, ω is a diagonal matrix with the natural frequencies of the reference system, δ is a diagonal matrix with the delay for each mode, and is a diagonal matrix with the amplitude error for each mode. With the equivalent properties of the hybrid system, the state matrix A eq is constructed as shown in Equation (24): The hybrid system is stable if the matrix A eq has eigenvalues with negative real part; thus, the critical delay can be obtained by searching the minimum constant delay that produces eigenvalues with positive real part. For the partitioning of this problem, assuming = I and δ = diag(τ cr ; τ cr ; τ cr ), the critical delay is: This critical delay serves as a reference for the simulations with constant delay. It is expected that simulations with delays smaller than this critical value would remain stable while running the simulations. In contrast, the opposite will result in an unstable response.

RTHS With Constant Delay
The Simulink block diagram for the simulations with constant delay is presented in Figure 4, where the output of the numerical substructure, denoted as "target vector, " is defined as x t = x (1) n ;ẋ (1) n ;ẍ (1) n T (displacement, velocity and acceleration of first degree of freedom). Then, the Delay block generates a discrete delay of the input signal x t by a specified number of samples d, so the time delay τ is a multiple of the simulation time step. Then, the delayed boundary conditions are: Finally, the outputs of experimental substructure are described in Equations (27) and (28): f e = k e ; c e ; m e x d Four simulations are performed with different delay values from 0 to 9 ms. The structure is subjected to the Kobe earthquake with a peak ground acceleration (PGA) scaled to 40%. The measured displacement of experimental substructure x m is compared with the first degree of freedom displacement x r1 from the reference structure, as shown in Figure 5. Measured displacements are very close to reference displacements for lower values of delay (τ = 0 ms and τ = 0.98 ms). Then, for higher delay τ = 7.1 ms, a notorious synchronization error in displacements appear, but at least the measured displacement decays at the end of the simulation. This does not happen with τ = 9 msec > τ cr , where the measured displacement grows exponentially after the earthquake reaches zero. In Figure 6, the error of displacement x m respect to x r1 for each delay is presented, including the normalized root mean square error (NRMSE) for each case. Clearly, for τ = 0 ms and τ = 0.98 ms, the system is stable; while, for τ = 9 ms is unstable (i.e., error grows exponentially). However, for τ = 7.1 ms, the error is considerable, with NRSME = 144%, but the system is still stable. It is easy to classify the case with τ = 9 ms as unstable after  the test. Still, it would be beneficial to detect unstable behavior before the system reaches a large response.
In Figure 7, the SW indicator is calculated for each model scenario. For τ = 0 ms, τ = 0.98 ms, and τ = 7.1 ms, the indicator stays under 100% during the entire test, showing the stable behavior of each simulation. In contrast, for the unstable case with τ = 9 ms, the indicator exceeds 100% near 10 s of simulation, much earlier than the instant when the system reaches the large displacements presented in Figure 5. Thus, the SW indicator could help stop the simulation in time, avoiding any dangerous behavior in the RTHS test.
For the sake of comparison, Figure 8 shows the absolute value of EEI for each case scenario. The EEI grows up with higher values of delay but is not evident in which case is stable and in which it is not. For τ = 0 ms, the EEI values are relatively small. Moreover, for τ = 0.98 ms, the EEI stays under 15% during the test; both cases are indeed stable. However, for τ = 7.1 ms, the behavior of EEI is similar to τ = 9 ms before 15 s of simulation, so it is not possible to detect instability until the earthquake signal dies out. If the EEI is utilized to stop unstable testing, it is quite challenging to define a threshold that relates EEI with instability. For example, for a threshold of EEI = 100%, the simulations  with τ = 7.1 ms and τ = 9 ms should be stopped before 10 sec, even when the simulation with τ = 7.1 ms is stable and does not present any risk of large displacements. Nevertheless, it must be recognized that EEI is an excellent tool to assess the accuracy of the results from an RTHS test through a post-processing analysis.

RTHS With Actuator Model and Dynamic Compensation
In this subsection, the same partitioning as in the previous subsection is utilized, but in this case a transfer system is modeled as presented in Silva et al. (2020), and adaptive compensation with different design parameters is implemented to minimize the errors produced by the transfer system. The block diagram utilized to model this problem in Simulink is presented in Figure 9, where the output of the numerical substructure is the target displacement x t = x (1) n (displacement of first degree of freedom). Then, the target displacement is sent to the adaptive controller to determine the commanded displacement for the actuator. The plant (i.e., actuator connected to the experimental substructure)  is modeled with the transfer function presented in Equation (29), which includes actuators parameters and interaction with the experimental substructure, and relates the commanded displacement x c with the measured displacement x m : x m = 2.19 × 10 13 (29.1)s 5 + (1.26 × 10 4 )s 4 + (8.42 × 10 6 )s 3 + (2.33 × 10 9 )s 2 + (5.44 × 10 11 )s + (2.17 × 10 13 ) x c where s = iω, is the Laplace variable, i is the complex number and ω is the circular frequency. This transfer function has a frequency-dependent phase. The phase is another way to state that the measured and commanded displacements are delayed. Thus, the delay is frequency-dependent for stable and underdamped linear systems. The interaction forces are obtained with Equation (30): Besides, sensor noise has been added to the measured displacement x m and the experimental force f e using bandlimited white noise (BLWN) blocks to model the physical sensors. Then, the signal x m is sent to the controller for the adaptation process. Consequently, the signal f e is fed back to the numerical substructure to be considered in the integration of the numerical EOM. Notice that force measurement noise could impact the experimental force power, but it has a small effect on the work done by this force due to the numerical integration process to calculate work. However, measurement noises can affect the hybrid loop, depending on the structural properties, transfer system, and dynamic compensation. The adaptive compensation implemented in the controller consist in adaptive model-based compensation (Chen et al., 2015)  with the modifications presented in Fermandois et al. (2020). A first-order adaptive feedforward determines the commanded displacement as shown in Equation (31): where a 0 and a 1 are the adaptive control parameters, andẋ t is approximated using the backwards difference method. The adaptation laws of a 0 and a 1 are described in Equation (32): where Ŵ = diag(Ŵ 0 , Ŵ 1 ) is the adaptive gain matrix associated with the adaptation rate of parameters a 0 and a 1 , respectively, and e is an estimation error of the adaptive parameters described in Equation (33): Both signals x c and x m are filtered with a fourth-order Butterworth low-pass filter with a cutoff frequency of 20 Hz, anḋ x m is obtained with the backward difference method. This subsection aims to demonstrate the performance of the SW stability indicator for different compensation scenarios. Notice that the goal of the SW indicator is to detect instability during a test with a predefined compensation design, and not to tune the compensation/adaptation parameters. In these simulations, the initial adaptive parameters are selected as a 0 (t = 0) = 1 and a 1 (t = 0) = 10/1, 000 s. With these parameters, FIGURE 10 | Reference and measured displacements for different adaptive gains.
Frontiers in Built Environment | www.frontiersin.org the controller cannot compensate for the control plant delay successfully; therefore, parameter adaptation is necessary.
Four simulations are performed with different adaptive gains presented in Table 1, where Ŵ a corresponds to a non-adaptive case, Ŵ d to an optimally calibrated case (Fermandois et al., 2020), and Ŵ b and Ŵ c are intermediate cases of adaptation. The structure is subjected to the Kobe earthquake scaled to 40% PGA.
The measured displacement of experimental substructure x m is compared with the first degree of freedom displacement x r1 from the reference structure in Figure 10, where the case FIGURE 11 | Displacement error for different adaptive gains.
FIGURE 12 | Adaptation of control parameters a 0 and a 1 for different adaptive gains.
Frontiers in Built Environment | www.frontiersin.org with Ŵ a is unstable, while the case with Ŵ b presents large measured displacements, but that decreases after the earthquake. Considering that in true RTHS tests, the reference displacement is not available, it would be difficult to distinguish if the measured displacements are correct or are associated with an unstable response. On the other hand, for Ŵ c and Ŵ d simulations are stable.
In Figure 11, the error of displacement x m respect to x r1 for each adaptive gains is presented, including the normalized root mean square error (NRMSE) for each case. The cases with Ŵ a present unacceptable errors, while for Ŵ b , the errors are substantial, but at least decrease after the earthquake. For the cases with Ŵ c and Ŵ d , the errors are relatively small.  The difference between each case is explained by the adaptation process presented in Figure 12, where there is no adaptation for Ŵ a , slow adaptation for Ŵ b , and the same convergence for Ŵ c and Ŵ d but with different rates. The parameter a 1 is directly related to the delay compensation, so it is expected that for Ŵ a and Ŵ b , the non-compensated delay produces large errors in the simulations. However, the adaptation for Ŵ b is apparently enough to stabilize the system during the test.
In Figure 13, the SW indicator for each case is presented. For Ŵ c and Ŵ d , the indicator is growing at the beginning of the earthquake. Then, the SW indicator decays due to the adaptation and good compensation, and stays under 100% during the entire simulation, confirming the stability of these simulations. For Ŵ a and Ŵ b , the indicator reaches 100% near 10 s of simulation, so much before each simulation presents large displacements. So, the stability warning indicator could help stop the simulation in time before any damage occurs in the test setup. In the particular case with Ŵ b , the stability warning indicator exceeds 100% between 10 and 20 s, consistent with the large response of the system in that interval. After 20 s of simulation, when the controller has reached enough adaptation to compensate for the delay, the SW indicator decays below 100%. This result is consistent with the stable behavior of this case at the end of the simulation. However, this simulation could be stopped when the indicator reaches 100% to avoid the behavior shown between 10 and 20 s.
Moreover, notice that for Ŵ b and Ŵ c cases, the SW does not decay to zero at the end of the simulation. This situation occurs because of uncompensated delay during the earthquake, which causes a surplus of energy inserted into the numerical substructure. Later, when the delay is compensated, displacements and forces are too small, so the work done by the experimental forces remains constant. This effect does not happen for Ŵ d , where the delay is compensated early when the earthquake strikes. The experimental substructure dissipates energy during the rest of the earthquake, leading W F to zero and consequently, SW to zero as well.
Furthermore, the control plant in the benchmark problem is modeled as a linear time-invariant system, but in these simulations, the adaptive controller produces time-varying compensation. Therefore, the delay between target and measured displacements is a time-varying property. This scenario is a complicated situation to predict the stability of a system before the test. However, with the online analysis of the SW it is possible to detect which case results in stable behavior and which are not.
For comparison, Figure 14 shows the absolute value of EEI for each adaptive case scenario. For Ŵ a and Ŵ b , the EEI presents large values, while for Ŵ d , the EEI keeps low during the test. So, the EEI could allow lab technicians to stop the simulations with Ŵ a and Ŵ b , and accept the simulation with Ŵ d . However, for Ŵ c case, the EEI reaches 100% before 10 s where the displacements are small, and the controller is adapting, and the test still has good chance to give acceptable results without the risk of instability. Thus, the EEI is an excellent tool to analyze the errors after the RTHS test, but it is not reliable for detecting instability online.
Additionally, in these simulations, displacement and force measurement noises do not have a notorious impact on the hybrid system response. The numerical substructure acts like a low-pass filter for the force measurement, and the adaptive compensator is not affected by the displacement noise.

CONCLUSIONS
Energy analysis for real-time hybrid simulation tests is conducted to establish a stability indicator for online evaluation. The features of the proposed indicator to detect instability is demonstrated through numerical simulations with constant delay, including an actuator model with adaptive compensation for different design scenarios. The stability warning indicator (SW) allows the detection of unstable behavior before the system reaches large displacements that can damage the experimental substructure or laboratory equipment. The proposed indicator only requires the measured force from the experimental substructure and information from the numerical substructure. Hence, it is independent of the compensation method and does not require a model of the experimental substructure or the transfer system (i.e., actuation system). This indicator can be implemented to stop an RTHS test that is rendered unstable automatically. The test should be stopped if the SW reaches 100% (i.e., it is unnecessary to stop a test with a SW < 100). So, this indicator is an excellent complement to any compensation method. It provides further safety guarantees to the test, especially for adaptive compensation if there is considerable uncertainty in the plant or the controller's adaptation capacity. It is worth mentioning that the proposed indicator provides safety for the test, but it does not ensure accuracy compared to the reference structure's true response. However, priority should always be given to conducting a safe experiment. Then, the reliability of the results could be studied in detail after the test has ended.
Finally, the effectiveness of the SW indicator is demonstrated for linear systems. Future work will consider non-linear systems and the corresponding experimental validation of the proposed indicator, together with the implementation of existing compensation techniques in RTHS testing and its synergism.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/supplementary material.

AUTHOR CONTRIBUTIONS
CG and GF conceived the presented idea for this study. CG developed the methodology and performed the numerical simulations. GF supervised the findings of this work. CG wrote the manuscript with support from GF. All authors discussed the results and contributed to the final manuscript.