Towards Data-Driven Real-Time Hybrid Simulation: Adaptive Modeling of Control Plants

We present a method for control in real-time hybrid simulation (RTHS) that relies exclusively on data processing. Our approach bypasses conventional control techniques, which presume availability of a mathematical model for the description of the control plant (e.g., the transfer system and the experimental substructure) and applies a simple plug 'n play framework for tuning of an adaptive inverse controller for use in a feedforward manner, avoiding thus any feedback loops. Our methodology involves (i) a forward adaptation part, in which a noise-free estimate of the control plant's dynamics is derived; (ii) an inverse adaptation part that performs estimation of the inverse controller; and (iii) the integration of a standard polynomial extrapolation algorithm for the compensation of the delay. One particular advantage of the method is that it requires tuning of a limited set of hyper-parameters (essentially three) for proper adaptation. The efficacy of our framework is assessed via implementation on a virtual RTHS (vRTHS) benchmark problem that was recently made available to the community. The attained results indicate that data-driven RTHS may form a competitive alternative to conventional control.


INTRODUCTION
Whilst numerical simulation methods play an ever increasing role in the analysis and design of structures, these remain insufficient in the case of complex structural systems under extreme loading conditions. As such, physical tests cannot be fully removed from the analysis and assessment process. Hybrid testing (Takanashi and Nakashima, 1987;Mahin et al., 1989) can work complimentary with numerical modeling by allowing physical testing of those regions or components of greatest interest or complexity without incurring the high costs associated with physical testing of the whole system. The presence of significant non-linearities are often a motivation for the use of hybrid simulation due to the difficulties associated in properly modeling this behavior. Non-linear components can often exhibit considerable rate dependent behavior, with this rate-dependency often not adequately compensated by scaling. This motivates the use of RTHS as the "gold standard" of hybrid simulation. Real-time in this case implies that the time scales of the numerical and physical system are the same; this allows for tests which incorporate rate dependent non-linearities in the physical component and are hence more representative of the true system (Nakashima et al., 1992;Nakashima, 2001;Benson Shing, 2008).
The real-time aspect of RTHS brings significant difficulties in comparison to an increased time scale test or pseudo-dynamic testing (Bayer et al., 2005;Pegon, 2008). Notable challenges involve the integration of the numerical system, which now must be executed in real-time, as well as enforcing robust and accurate actuator control. The ability of a hybrid simulation to recreate realistic testing conditions is reliant upon the accuracy of the applied boundary conditions. As such, the accurate recreation of the signals from the numerical substructure by the actuator are of paramount importance. With real-time testing, the accurate reproduction of these signals becomes more challenging. The dynamics of the actuation system and interaction with the testpiece can have significant effects on the reproduction of the reference signal as well as the effect of signal processing artifact. This results in two key control issues when implementing a RTHS scheme, namely, (i) the accurate reproduction (through the transfer system) of the reference signal that corresponds to the common boundaries between the experimental and the numerical substructure; and (ii) the suppression of the time delay, which is inevitably introduced by the transfer system (actuators, analog-to-digital and digital-to-analog converters, ADC and DAC, respectively, etc.) and possibly by the control scheme itself. The first issue is a typical control problem, while the second is a problem of prediction.
Initial control approaches implemented for RTHS focused on time domain approaches. The problem was first considered by Horiuchi et al. (1999), wherein a polynomial extrapolation method was used to reduce actuator delay. This method was further developed allowing for adaptation of the polynomial coefficients in Wallace et al. (2005). Further work exploited model-based control, wherein a model of the control plant is identified and used to make corrections to the reference signal for improving tracking (Carrion and B. F. Spencer, 2006;Chen and Ricles, 2009). Current state of the art methods focus on adaptive control schemes due to their superior robustness. Adaptive model based control schemes as in Najafi and Spencer (2019) demonstrate good performance, albeit present the drawback of requiring a mathematical model of the plant to be formulated. Ning et al. (2019) also demonstrate an adaptive method making use of an H ∞ filter for tracking error and a polynomial extrapolation for delay compensation. This method yields good performance but requires an offline identification stage in which a second order model form is assumed and identified to replicate the plant behavior.
Within the framework of methods for adaptive tracking control, a class of algorithms is exclusively based on data processing, rather than the integration of conventional control. As such, they require no prior knowledge on the dynamics of the transfer system, the experimental substructure and their interaction. Among other works, the adaptive time series compensator developed by Chae et al. (2013) has shown very good performance, but requires the careful prescription of userdefined limits on its hyper-parameters for ensuring the stability of test. Dertimanis et al. (2015a,b) apply adaptive signal processing concepts and succeed in reproducing reference signals with a high degree of efficacy. These algorithms can be applied online and in hard real-time, thus proving particularly beneficial in terms of robustness as the filters are adaptively optimized to the test conditions, even if these change between tests. However, although they have been experimentally validated, via a 10 Mgr linear specimen on a shaking table, they have not yet been integrated into RTHS.
The present study extends the research conducted in Dertimanis et al. (2015a,b) and proposes a data-driven, adaptive inverse control (AIC) framework for RTHS. It proceeds by first formulating a set of specifications that data-driven control schemes should meet and then splits the adaptive process into two phases, which can be executed either simultaneously, or successively, both in online and offline mode. In the first phase, the decorrelated least mean square (LMS) algorithm is applied to the identification of the control plant, while in the second phase the discrete-cosine transform LMS (DCT-LMS) takes on the tuning of the inverse controller. It is shown that the cascade of the latter and the control plant closely approximates a perfect delay, allowing thus the signal of the numerical substructure to be driven to the experimental substructure with unaltered dynamics. Accordingly, the method adopts a polynomial extrapolation method (Horiuchi et al., 1999;Wallace et al., 2005) for the compensation of the time delay, ensuring thus the stability of the RTHS loop. A particular advantage of the developed scheme is its dependence to a very small number of hyper-parameters that must be provided by the uses, an important feature that favors robustness, safety, transparency, and ease of use.
We assess our AIC method via the recently established linear vRTHS benchmark problem of Silva et al. (2020). It must be emphasized that, whilst the benchmark problem dealt with herein considers the control of a linear plant, the AIC framework has previously been successfully extended to nonlinear systems on multiple occasions. In Widrow and Plett (1996), it is discussed how the AIC framework for linear filters can be extended to take into account non-linear filtering. The demonstration of AIC to a theoretical non-linear plant then followed in Widrow et al. (1998), wherein the nonlinear filters took the form of neural networks whereby the weights are updated with a gradient descent method. A full consideration of the non-linear AIC framework is given in Plett (2003). Alternative methods have made use of other regression techniques such as support vector machines for constructing non-linear filters (Yuan et al., 2008). Physical implementations of non-linear AIC have been demonstrated on piezo-electric actuators featuring hysteresis effects (Li and Chen, 2013), on a magnetic bearing system (Jeng, 2000) and in the control of an electronic throttle system (Xiaofang et al., 2010). Furthermore, the AIC method using linear filters has been applied to physical systems that are similar to those considered in this benchmark problem, i.e., electrohydraulic shaking table, and demonstrated success both in the case of conventional real-time dynamic tests (Shen et al., 2011;Dertimanis et al., 2015a).
The contributions of this study are (i) the treatment of RTHS in a purely data-driven fashion; (ii) the establishment of a set of specifications that such an approach should meet; (iii) the derivation of a simple, plug-n-play, adaptive modeling framework for the proper estimation of a feed-forward adaptive inverse controller for RTHS; and (iv) the assessment of our methodology via the vRTHS benchmark problem of Silva et al. (2020). The context is structured as follows: in section 2 the problem is formulated, along with the list of specifications and a brief introduction to adaptive filtering is offered. Section 3 outlines the method, while section 4 contains the application study. Finally, in section 5 the outcomes of the study are summarized directions for future research are given. Figure 1A illustrates a typical RTHS loop, where a structure is split into a numerical and an experimental part and their interaction is assured via the use of a transfer system. A minimal configuration of the latter contains a set of ADC and DAC devices, a controller and an actuator (oftentimes integrating an inner control loop) that is firmly attached to the experimental substructure. Inevitably, this firm attachment causes an evolution of dynamics along two paths: this is termed controller-structure interaction (CSI) in the literature and is represented by an additional feedback loop (not shown in Figure 1, refer to Figure 1). The control plant thus consists of the transfer system and the experimental substructure.

PROBLEM FORMULATION
In more detail, it is assumed that an external excitation, u[k] is applied to the numerical substructure and the kinematics at the boundary are calculated by applying an appropriate numerical integration scheme that solves the equations of motion (Shing, 2008). The calculated reference signal, x R [k], which can be in the form of displacement, velocity, or acceleration, is then applied though the transfer system to the experimental substructure. The response of the latter is monitored and the dynamics at the boundary, usually in the form of a force, is fed back as additional input to the numerical substructure, in order to proceed to the next step of the loop.
For the effective implementation of a RTHS loop, two fundamental problems must be solved. The first corresponds to the dynamics of the control plant and is treated by introducing an additional controller (termed the RTHS controller in Figure 1A). The second problem pertains to the delay that the transfer path inherits, which is tackled via an appropriate delay compensation method. Under this enriched configuration, the reference signal is predicted steps forward, where is an estimate of the transfer path's overall delay, and then appropriately modified via the RTHS controller to form the command signal for the transfer system. A careful tuning of all the individual blocks results in an achieved signal at the experimental substructure that is equal to the reference one, e.g., ; this equality forms the hard constraint of any RTHS loop. Under this setting, the challenge that this study aims to tackle is summarized in the following question: is it possible to tune a RTHS controller and a delay compensation algorithm without any prior information on the dynamics of the control plant, given only the availability of the reference and achieved signals?
Leaving out the delay problem for the moment, the datadriven attack to the establishment of a RTHS controller is illustrated in Figure 1B and reads as follows: under the availability of x R [k] and x A [k], a controller that forms the inverse of the control plant should be estimated, causing the cascade to perform as a perfect delay, e.g., Then, the addition of the delay compensation method should result in fulfilling the hard constraint condition of the RTHS loop. In designing such a controller, Dertimanis et al. (2015a) described a set of specifications, which are herein reformulated and enriched. According to these, a controller should be: 1. Data-driven: no need for analytical models of individual components (valves, cylinders, etc.) and identification of system parameters (stiffness, damping, oil constants, flow rates, etc.) 2. Feedforward-driven: no need for additional feedback loops. 3. Discrete-time oriented: no need for discretization of continuous-time models. Everything should be digital. 4. Minimally parametrized: the number of parameters, henceforth referred to as hyper-parameters, required for the proper tuning of the controller should be kept as small as possible. The controller should not be "too sensitive" on these hyper-parameters. 5. Of guaranteed stability: exclusive use of finite impulse response (FIR) models (refer to section 2.1), instead of infinite impulse response models in order to enforce stability and safety during test. 6. Robust: the RTHS controller should compensate for all uncertainties of the control plant, as well as the CSI problem. 7. Of minimum discrepancies: the cascade frequency response should follow the zero dB line in the maximum possible frequency band. The phase delay should be constant within this band. 8. Straightforward to implement to existing facilities: the RTHS controller should be opeational in conjunction with conventional fixed-gain controllers. 9. Applicable to a wide range of transfer systems: from smallscale actuators and light specimens, to shaking tables and specimens of several megagrams. 10. Functional for all types of command signals: both acceleration and displacement reference signals should be handled. 11. Straightforward to realize-execute: immediate implementation in commercial hardware and execution in real time. No need for sophisticated software design.
Our approach for establishing such an inverse controller utilizes the theory of adaptive signal processing (Widrow and Wallach, 2007;Diniz, 2008) and it can be implemented either online, or offline. Its effectiveness in modeling facilities for structural testing has already been demonstrated in Dertimanis et al. (2015a,b).

Brief Review of Adaptive Filtering
For convenience, this section offers a brief outline of the most fundamental concepts of adaptive signal processing. The familiarized reader can safely skip this summary. A sound treatment to the topic is given in Diniz (2008), Hayes (1996), and Manolakis et al. (2005). Glentis et al. (1999) provide an excellent review on adaptive filters, including the presentation of a comprehensive review of associated algorithms.
Consider an unknown plant that can be effectively described by its impulse response g [k] in the discrete-time domain, or by its transfer function G(z) in the Z-domain. The plant is driven by a wide-sense stationary input signal u[k] and the response x[k] is measured in a noise-corrupted fashion, under the assumption of an additive disturbance d[k] at the plant's output, which is random and uncorrelated to x[k]. Using an FIR parametrization, we can represent the input-output dynamics as (1) with n corresponding to the model order. For g = (1) can be cast into a regression form Suppose now that measurements of u[k] and x[k] are acquired and the aim is to estimate a FIR model of the plant. Our model readsx wherex[k] is the model's output andĝ the vector of unknown filter weights. This is a typical parametric identification problem, which can be solved by collecting measurements over a time interval and then solving a linear least-squares problem for recovering g. Clearly, this non-recursive strategy cannot be applied in real-time, since it requires batch-processing of stored data. An alternative approach is to proceed in a recursive estimation of the weights, whenever new data becomes available. The key idea, illustrated in Figure 2, is to update the weights via FIGURE 2 | Sketch of the adaptive filtering concept for the system identification problem.
the stochastic approximation of a deterministic optimization algorithm, in which the direction is calculated using unbiased estimates on the basis of the current time index. For example, consider the steepest descent algorithm where µ is the step size and V(ĝ[k]) is the objective function, defined as the instantaneous mean square error (MSE) Plugging Equation (3) and differentiating with respect tô Approximating the expectation operator and substituting to Equation (4) yieldŝ This is the celebrated least mean square (LMS) adaptive filter developed by Widrow and Hoff (1960). Expectedly, the behavior of the algorithm depends on the step size. When d[k] is widesense stationary and the unknown plant is time-invariant, the LMS filter converges in the mean to the optimal Wiener solution, provided that µ is bounded as where λ max is the largest eigenvalue of the input's autocorrelation matrix The condition of Equation (8) does not, however, ensure stability. This is succeeded by an expression that is widely used in practice, since it is based on the energy of the input signal, which is easier to calculate than the eigenvalues. Finally, maximum convergence speed is achieved when with λ min denoting the smallest eigenvalue of Ŵ uu . Equation (11) indicates that the speed is controlled by the eigenvalue spread of Ŵ uu , e.g., From Equations (11,12) it is easy to conclude that convergence speed requires an eigenvalue spread close to one. When the input signal can be selected by the user, the best option is to let u[k] being a realization of a zero-mean Gaussian white noise process, since Ŵ uu = σ 2 uu I n+1 and the spread is equal to one.

The Problem of Inverse Identification
Assume that the control plant can be described by the following digital rational transfer function Then, the problem of identifying an inverse controller is reduced in approximating the inverse transfer function If G(z) is strictly minimum phase (e.g., poles and zeros inside the unit circle), then the problem has a straightforward solution, since G −1 (z) will be also minimum phase, admitting a convergent expansion of the form Thus, by truncating the infinite sum up to an order n I , it is possible to derive a FIR representation for the inverse controller. If, however, G(z) is non-minimum phase, the poles of G −1 (z) [e.g., the zeros of G(z)] are located outside the unit circle and the power expansion of Equation (15) diverges. To cope with this issue, recall that a digital transfer function admits two inverse Z-transforms, a causal and a non-causal one, depending on the region of convergence (Oppenheim et al., 1999, Chapter 3). In the case of the non-minimum phase G(z), the causal inverse transform is unstable, but the non-causal is stable. A stable, non-causal expansion of G −1 (z) can be written as where, usually, the causal part is limited to very few terms (e.g., n 1 is small). If the weights of the non-causal part are significant, attempting to adapt an inverse controller would render quite poor results. However, if one would consider multiplying Equation (17) by a pure delay, z − , it would cause weights to jump from the non-causal to the causal part, increasing thus the accuracy of the inverse (Widrow and Wallach, 2007, Chapter 5).
To summarize, the adaptation of the inverse controller aims at satisfying the following equation In practice, this implies that the reference signal passes through a pure delay before being fed to the adaptation block. It is emphasized that, even if the continuous-time transfer function of a plant is minimum phase, its digital counterpart may oftentimes result non-minimum phase.

Description
The control strategy implemented herein essentially pertains to applying the adaptive filtering concepts previously outlined in identifying an appropriate inverse controller. The latter is then placed before the plant in the signal path, such that the cascade results ideally in a delayed unit impulse response. Accordingly, by integrating an appropriate delay compensation method, the achieved signal may follow the reference one with a reasonable accuracy, at a certain frequency band. The strategy used for the adaptation of this inverse controller is based on the simple architecture shown at Figure 3 and consists of two stages that are implemented in (hard) real-time: These two phases can either be performed successively or simultaneously, in the sense that it is not necessary to wait for the full convergence of the forward controller, in order to start the adaptation of the inverse one (Widrow and Wallach, 2007). The adaptation process initiates by supplying a reference signal, x R [k], of favorable properties (see below), to the control plant. This signal, together with the noise-corrupted achieved response of the plant, x A [k] are fed to the forward adaptation block, for estimating a FIR filter that describes the end-toend reference-achieved signal dynamics. The term end-to-end herein implies that the path from the reference to the achieved signal contains all individual software/hardware components, including ADCs, DACs, signal conditioners, etc., besides the control plant.
When the forward and inverse adaptation are carried out simultaneously, e.g., in synchronous mode, the instantaneous estimate of the achieved response,x A [k], together with a delayed version of the original reference signal, are fed to the inverse adaptation block and a FIR filter is tuned for describing the achieved-(delayed) reference signal dynamics. The output of this block is an estimate of x R [k − ] (see Figure 3). The delay is an algorithmic parameter (it is included in the set of method's hyperparameters, listed in Table 1) and its size depends heavily on the qualitative characteristics of the control plant.
A key detail of the presented algorithm is that the inverse filter weights are not adapted by directly feeding the output of the control plant to the inverse block, in order to avoid the propagation of the disturbance. When both stages have been successfully executed and convergence has been achieved, the inverse FIR filter can be copied before the control plant and it is implemented as an inverse controller. Ideally, the cascade of the inverse controller and the control plant form a unit impulse that reproduces a delayed version of any reference signal. It is noted that more sophisticated algorithms can be implemented as adaptive inverse controllers, as the modified filtered-X one presented in the work of Dertimanis et al. (2015a).

Forward Adaptation
Our choice for the algorithm that implements the forward adaptation path is based on (i) the reported convergence rates; and (ii) the reduction of the required tuning parameters. Based on previous results (Dertimanis et al., 2015a,b) we apply the decorrelated LMS algorithm Doherty and Porayath (1997), which belong to the class of instrumental variable methods (Glentis

Requires careful tuning
Boldface entries indicate the minimum required subset of hyper-parameters that must be supplied by the user.
et al. , 1999). At each step, the adaptive filter's output is calculated asx where and n F is the order of the filter. The coefficients are updated bŷ for a step size µ[k] In Equation (23), ǫ F is a small constant ("stability" parameter), is the error between the measured and the predicted signal and w[k] is the filter gradient, updated by with a[k] denoting the decorrelation coefficient The properties of the decorrelated LMS algorithm are studied in Doherty and Porayath (1997). See also (Douglas et al., 1999;Rørtveit and Husøy, 2009). As an alternative to the decorrelated LMS, the normalized LMS algorithm (Diniz, 2008, Chapter 4) can be also applied to the forward adaptation path, at the cost of an additional hyper-parameter that controls the adaptation step.

Inverse Adaptation
The signal that arrives at the input of the inverse adaptation block,x A [k] is, in general, highly correlated (e.g., colored), since it describes the dynamics of the control plant. Notice that when the reference signal is Gaussian-like, the inverse adaptation pertains essentially to a whitening process. This is because, the cascade of the inverse controller and the plant should give an achieved signal which is uncorrelated like the Gaussian input. As the convergence rate of the conventional LMS algorithm (including its normalized version) is adversely affected by correlated inputs, it cannot be considered as a competitive candidate. Instead, an effective solution to this problem is offered by transform domain methods, which utilize an orthogonal transformation to a space that is attributed with favorable properties (Beaufays, 1995;Diniz, 2008;Chergui and Bouguezel, 2017).
for notational convenience (see Figure 3), the output of the inverse filter reads For any orthogonal matrix S ∈ R n I , the regression form of Equation (26) can be written as where ǫ I is a small constant ("stability" parameter), and γ → 1 is the power update factor. The transformed weights are then updated by a typical LMS filtering operation with µ being the step size, e I There's a plenty of options for the selection of the orthogonal matrix S, including the discrete Fourier transform (DFT), the discrete Hartley transform (DHT) and the discrete cosine transform (DCT). Herein, we adopt the third option and we construct the orthogonal matrix as where Q p,q = 2 n I K p cos p(q + 1/2)π n I , p, q = 0, 1, . . . , n I − 1 (32) for Since the DCT transform is fully parametrized by the filter order, the orthogonal matrix S can be formulated and stored prior to the adaptive modeling process. Further details on the performance of transform domain adaptive algorithms can be found in Zhao et al. (2009), Kim and Wilde (2000), Lee and Un (1986), and Narayan et al. (1983).

Delay Compensation and Method's Hyper-Parameters
The successful adaptive modeling of the control plant and its inverse forces, ideally, the cascade to behave as a delayed unit impulse. This implies that a reference signal is driven to the control plant with unaltered dynamics, but at a delay . To compensate for this delay, we integrate a one-step ahead prediction method, which is based on polynomial extrapolation (Horiuchi et al., 1999;Darby et al., 2001;Wallace et al., 2005). The reference signal that is driven to the inverse controller is where n P is the polynomial order and p i the polynomial coefficients, which are calculated via the Lagrange basis function and are predefined for a given order. Table 1 gathers all the hyper-parameters of the method. The most important of them pertain to the orders of the forward and the inverse adaptive filters, as well as to the step size of the DCT-LMS: these are essentially the quantities that the user has to decide for, every time a new experimental substructure is attached to a transfer system, composing thus the control plant. The default order of the polynomial extrapolation is sufficient, as long as the cascade is close to a pure delay. Higher orders do not, in general improve accuracy; actually they might oftentimes lead to instability, especially when is large.

The vRTHS Benchmark Structure
In evaluating the performance of novel control schemes for RTHS in a safe manner, vRTHS offers a useful platform. VRTHS involves the implementation of a hybrid simulation fully in silico, whereby both the numerical substructure and the physical substructure are simulated. Crucially, vRTHS also involves recreation of the transfer functions associated with the actuation and signal artifacts associated to the actual physical implementation. The use of vRTHS can allow for rapid and safe assessment of various algorithms relating to the hybrid simulation procedure. The robustness of any algorithm can also be investigated by the introduction of uncertainty in the vRTHS.
The study of suitable control and delay compensation algorithms for RTHS has resulted in proposition of alternatives schemes, which have been investigated on disparate testing setups. However, hybrid simulation setups tend to be unique and hence it is difficult to fairly evaluate the relative performance of control schemes implemented on different setups. This provides the motivation for establishing benchmark problems, wherein the only variable is the control regime implemented.
Such a benchmark problem is recently proposed by Silva et al. (2020), who developed a vRTHS framework for the purpose of evaluating different control regimes on a unified system. This benchmark consists of a 2-bay, 3-story steel frame, in which the structural mass is concentrated on the floors, motion is allowed only horizontally on a single direction and damping is proportional. An originally 30-degrees of freedom (DOFs) linear elastic planar finite element model is reduced to a 3-DOF model, pinned at ground and excited on its base, which consists the reference structure (Figure 4). The benchmark is implemented in SIMULINK © and is created such that the controller block can easily be exchanged, whilst the rest of the system is left unchanged. This allows for the fair comparison of controllers. A number of uncertainties in the system parameters is also offered, to allow for analysis of controller robustness, along with a set of standardized performance metrics to aid comparison.
For the vRTHS tests the reference structure is partitioned into an "experimental" and a numerical substructure, with the former corresponding to a single bay of the ground floor and the latter consisting of the remainder of the structure (Figure 4). The structural properties of the experimental substructure are kept fixed, while the ones of the numerical substructure vary in accordance to the partition schemes of Table 2 The dynamics of the control plant (transfer system plus experimental substructure) are described as shown in Figure 5. The open-loop transfer function between the command and the achieved signals in continuous-time reads for   The numerical values of all associated parameters are listed in Table 3. To incorporate a degree of uncertainty, for testing the robustness of a proposed controller, some parameters are defined as normally distributed random variables. All vRTHS tests are conducted in SIMULINK © through the integration of the structural equations of motion via an explicit 4th-order Runge-Kutta numerical framework, at a fixed sampling rate F s = 4, 096 Hz. It is emphasized that the choice of the integration scheme has detrimental effects on the evolution of the dynamics and, consequently, on the behavior on the proposed adaptive modeling method. To demonstrate this, Figure 6 displays the discrete-time impulse response of the control plant under two different discretization schemes: in the upper plot, g[k] is calculated via the aforementioned integration, by applying a unit impulse excitation to the part of the SIMULINK © model that corresponds to Figure 5. In the lower part, the impulse response is obtained by the impulseinvariance transformation of Equation (35), e.g., g[k] = T s g(t = kT s ), where g(t) is the continuous-time impulse response (the inverse Laplace transform of Equation 35). A quite different "damped" behavior is observed, which renders the integrationbased impulse response attaining a much slower decay rate, necessitating the use of very high order FIR models for effective adaptive modeling of the control plant.

Adaptive Modeling
For generating the reference signal, a Markov-1 process of the form is adopted, where e[k] is a zero-mean Gaussian white noise stochastic process of variance σ 2 ee . In deciding for the values of φ 1 and σ 2 ee , recall that the maximum and minimum eigenvalues of Ŵ uu are provided by the corresponding maximum and minimum of the power spectrum of u[k], which is given by Thus, λ max = (1 − φ 1 ) 2 σ 2 ee λ min = (1 + φ 1 ) 2 σ 2 ee and the eigenvalue spread for Ŵ uu is In order to maintain this number close to one, our choice is φ 1 = −0.01 (the minus sign is adopted for attributing S uu (f ) with low pass characteristics). The variance of e[k] is then determined by first setting 2σ uu = 0.01, for constraining 95% of the input's amplitude within the ±10 mm range. From the theoretical Markov-1 process variance solving for σ 2 ee leads to σ 2 ee ≈ 2.5 × 10 −5 . The adaptation process is carried out simultaneously, that is, the forward and the inverse adaptation counterparts operate in synchronous mode. Having established the statistics of the reference signal, several trials are performed for the rest of the hyper-parameters of the method (e.g., the forward and inverse filter orders and the step size of the DCT-LMS), while the total adaptation time frame is fixed to 60 s. Some critical observations during the evolution of the whole adaptation process can be summarized by the following remarks: -In all trials the convergence in synchronous mode is quite fast and robust against the DCT-LMS step size. The forward adaptation phase requires at least around 10-15 s for proper convergence. -Within this time the DCT-LMS makes a rapid improvement, indicating that full convergence of the forward filter is not necessary for adapting the inverse. -A proper selection of the DCT-LMS step size is important.
In many cases, a quite fast convergence is observed, yet, the final result is not optimal in terms of cascade dynamics, as the corresponding impulse response returns "noisy." -It is noted that a fair number of trials is conducted by applying non-Markovian inputs, as well as other adaptive algorithms. Regarding the former, signals composed from filtering Gaussian white noise at cut-off frequencies equal to 100 Hz and F s /2 are also tested. The result is a severe distortion in the filters' weights, followed by inability to converge within the specified time frame. Very slow convergence is also observed in the case where the normalized LMS algorithm replaces the DCT-LMS one for the adaptation of the inverse filter. This is, though, an expected result, attributed to the inability of the algorithm to cope with colored inputs.
Aspects of these remarks can be visually validated from Figures 7-10, which display the results for n F = 10, 000 and n I = 100 ( = 50, µ = 5 × 10 −3 ), our final choice for the forward and inverse FIR filters. The fact thatĝ F [k] is very long (Figure 7, top) is explained by the high sampling rate of the simulation (4,096 Hz): in absolute times this impulse response evolves over 2.14 s (the settling time of the control plant's step response is around 1 s). It is noted that when lower orders are considered (for example 8,000 and 9,000 weights), an apparent distortion around the last weights of the cascade is caused: in fact, this is a practical way of detecting that the forward filter requires higher orders. The absolute error (Figure 8, upper part) confirms that convergence has been succeeded already after 15-20 s. In contrast, the inverse filter dynamics are expanded in a considerably smaller time frame (around 24 ms, Figure 7, middle). The dominant part ofĝ I [k] is located around the chosen delay (which is typical for inverse filters) and there exist many leading and trailing weights that could, potentially, be treated as zeros. However, adaptation at lower order/delay pairs (indicatively 80/40, 60/30, and 40/20), or at fixed order (n I = 100) and lower delays (e.g., = 40, 30, and 20) does not improve the inverse controller. Thus, the 100/50 pair results as the lowest possible accurate choice for this plant.
The absolute adaptation error (Figure 8, bottom) exhibits fast convergence (after approximately 25-30 s), with longer time needed however for a finer tuning of the weights; it is observed that by reducing the adaptation time to half the fixed frame (e.g., 30 s) and keeping all other hyper-parameters unaltered leads to sub-par results. It is also worth mentioning that the converged values for |e I [k]| are higher in comparison to the ones of |e F [k]|. This is also expected since, as already mentioned in section 3.3, the inverse adaptation performs as a whitening filter. Indicatively, Figure 9 shows the synchronization plot between the delayed reference and the achieved signals for the last 0.25 s of the adaptation process. A normalized root-mean-square error, defined as of approximately 7% is estimated. The cascade ofĝ F [k] andĝ I [k], the first 1,024 weights of which are plotted in Figure 7, bottom, depicts a very good approximation to a pure delay. This is verified in the frequency domain against the perfect digital impulse with = 50; as Figure 10 displays, there is a quite close resemblance between the frequency responses and the associated phase delays between the estimated and the perfect delayed impulses. The highest distortion in amplitude, located around 40 Hz, is not more than −0.8 dB, whereas the largest difference in phase delay is approximately 0.5 ms.
In order to obtain a better insight on the effects of the observed distortions in the frequency domain, a number of open loop simulations is conducted with the adaptive inverse controller  placed before the control plant. Different signals are selected as references, including sinusoids at different frequencies, a chirp, low pass filtered Gaussian white noise, as well as the structural response of DOF x 1 for the El Centro earthquake (40% scaling). The results of the simulations are again interpreted in terms of synchronization plots and are shown in Figure 11. Apart from Quantity J 1 (ms) J 2 (%) J 3 (%) J 4 (%) J 5 (%) J 6 (%) J 7 (%) J 8 (%) J 9 (%) the case of 40 Hz (Figure 11C), which returns a NRMSE value of approximately 8%, all other simulations produce achieved signals of good quality in a wide frequency range.

vRTHS
We are now almost ready to integrate the estimated adaptive inverse controller to the vRTHS closed loop. To this end, the delay compensation scheme outlined in section 3.4 is applied and the reference signal that is sent to the adaptive inverse controller is calculated from Equation (34) as is the displacement of the interface DOF calculated at time k and = 50, which corresponds to a cascade delay of 12.2 ms. The numerical substructure is excited by the El Centro earthquake at 40% intensity, while for each case of Table 2 20 individual simulations are executed. The results are presented in terms of the criteria J 1 to J 9 that are reported in Silva et al. (2020) and listed in Table 4.
The following points summarize the performance of the proposed data-driven method: -The delay between the reference and the achieved displacement (e.g., J 1 ) is exactly zero in all nominal simulations, and very low in the perturbed simulations of all partitioning case studies. The statistics of the perturbed simulations are comparable, indicating robustness of the controller against uncertainty. It is noted that the total number of non-zero delay simulations are 3, 2, 5 and 5 in partitioning case studies 1 to 4, respectively. This is illustrated in Supplementary Tables 1-4. -The NRMSE between the reference and the achieved displacement (e.g., J 2 ) is also kept quite low, ranging from 2.85 to 3.29% in the nominal simulations. These numbers are within the range of the NRMSEs shown at the synchronization plots of Figure 11. The mean values of the perturbed simulations follow closely the ones of their nominal counterparts, whereas the corresponding standard deviations are one order of magnitude lower. -Expectedly, the normalized peak tracking error (e.g., J 3 ) returns a bit increased, compared to J 2 , but again is kept in very low levels. -Criteria J 4 to J 6 calculate the NRMSEs between the relative displacements of the storys during vRTHSs and the ones during the simulation of the reference structure (e.g., without substructuring), while J 7 to J 9 the normalized peak tracking errors of the same quantities. The associated errors never exceed 6% in the nominal cases, whereas the mean values of the perturbed simulations are higher (reaching up to 8.52% in partitioning case study 4), followed by increased dispersions. It is, however, noted that in the offered benchmark problem, the simulation of the reference structure is accomplished via a different discretization scheme (either zero, or first order hold), to the one implemented for vRTHS (e.g., the explicit 4th-order Runge-Kutta integration) and, as shown in Figure 6 there might be significant differences in the calculated structural responses.
We conclude that the efficacy and robustness of the adaptive inverse controller is quite competitive in all partitioning case studies. In comparing the method presented herein, the reader may refer to Silva et al. (2020)

CONCLUSIONS
In this study, we explore the possibility of conducting RTHS from a purely data-driven perspective. By setting a number of specifications that such a framework should fulfill, we demonstrate that this is indeed feasible and it is characterized by several potential advantages, including transparency, robustness and minimum tuning. Further improvements are possible, such as a more rigorous methodology for automating/adapting the selection of the method's hyper-parameters to the examined case study, the integration of the delay compensation to the adaptation process, as well as techniques for reducing the delay of the cascade's impulse and increasing the effective frequency band of RTHS. These are left as objectives of future work together with the design and conduction of an extensive experimental campaign for the validation of the method.

DATA AVAILABILITY STATEMENT
The software and data developed for this study are available upon request to the corresponding author.

AUTHOR CONTRIBUTIONS
TS: methodology, conduction of simulations, and manuscript preparation. VD: conceptualization, methodology, conduction of simulations, and manuscript preparation. EC: methodology, supervision, manuscript preparation, and funding acquisition.

FUNDING
This work was carried out as part of the ITN project DyVirt and has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement no. 764547.