Applicability of Geographically Distributed Simulations

Geographically distributed simulations (GDS) overcome the inherent limitations in capabilities of single research infrastructure to accurately represent large-scale complex power and energy systems within representative operating environments in real-time. The feasibility of GDS has been proven, however, there is a lack of confidence in its adoption owing to limited evidence of its stability and accuracy that ascertain its practical applicability. This paper presents detailed small signal stability models for GDS setups with two interface signals transformations. The models have been validated by empirical analysis and used for determining the boundaries for stable operation of GDS setups. For the common region of stability of the two transformations considered, accuracy analysis presented offers insights for their selection. This advancement, thereby, enables realisation of experimental setups that can cater for the growing need to design and validate operational schemes that ensure robust and resilient operation of critical national infrastructure.


I. INTRODUCTION
H ARDWARE-IN-THE-LOOP has proven to be a valuable approach to accelerate the validation and deployment of novel technologies, with an increasing number of laboratories across the world adopting the approach [1]. The rapidly increasing complexity of power and energy systems (due to increased decentralization, decarbonization and digitalization) presents a challenge for single research laboratories, as existing capabilities (computational, expertise or equipment) will have to expand to encapsulate the growing complexity to ensure rigorous validation and roll-out of novel technologies at pace. This has encouraged the development and establishment of the concept of geographically distributed simulations (GDS), allowing for complementary expertise, additional equipment, domains and computational resources to be harnessed across multiple laboratories that are geographically dispersed [2]. GDS involves the split of a power network into two or more subsystems for simulation at geographically dispersed laboratories. The subsystems are coupled over the internet, facilitated by low latency communications technologies and existing power hardware-inthe-loop (PHIL) interfaces. Some perceived benefits of GDS are summarized below, (not a comprehensive list of benefits): 1) Comprehensive Characterization: Where the real-time simulation capability is limited due to computational constraints of the digital real time simulator, the simulators across multiple laboratories can be utilized to realize more detailed large scale models for high fidelity simulations. Examples of such experiments have been reported in [3], [4]. 2) Representative Systems: With the transition of our power system towards an integrated energy system with decarbonization of heat and transport, the concept of GDS allows for domain specific laboratories (such as heat and electrical power) to be interconnected to evaluate next generation concepts to facilitate the transition. Even within a domain, the equipment at each laboratory is unique and therefore the concept of GDS allows for utilization of hardware equipment of different laboratories within an experiment enabling the realization of broader range and more representative testing scenarios [5], [6]. 3) Acceleration of Validation: GDS supports the acceleration of validation in two ways: i) Enables integration of specialist skills and services such as in [7], where the communications emulation capability of one laboratory was utilized for an experiment by another laboratory, both geographically dispersed. This saves time and effort in realizing services that might be readily available for utilization. ii) Facilitates the cooperation between industrial partners as sensitive models are no longer required to be shared, can be run at independent organization facilities with only a data link to other cooperating organizations [8]. This adds an additional level of intellectual property protection and helps accelerate the validation. In terms of industrial partners developing controls, the GDS concept further enables the testing of prototype controllers without having to be delivered to a laboratory, making it logistically more convenient [9], [10]. GDS implementations have been reported across the world, spanning from a country [11] to multiple countries [12] to continents [7] -with the longest reported link being between USA and Australia (13,570 km) [13]. A comprehensive review of the technological advancements in GDS has been reported in [14]. Although a promising approach that is gathering far reaching attention, its wide-scale adoption is hampered owing to limited understanding of its applicability.
As with PHIL simulations, the applicability of GDS is driven by the stability and accuracy of the setups. Stability and accuracy of PHIL setups have been extensively discussed in literature. PHIL setups are widely modeled as s-domain transfer functions that are subsequently utilized for stability analysis using Bode criterion [15], [16], Nyquist criterion [17] or Routh-Hurwitz criterion [18], [19], [20]. In [21], it was shown that the stability and accuracy of a PHIL setup is much influenced by the choice of interface algorithm, with a comparison of stability of commonly adopted interface algorithms presented in [22]. The impact of time delay on the stability of PHIL setups was discussed in [23], [24], [25], where deterioration in stability with increase in time delay was reported. In [19], the impact of power amplifier (switched mode) on the stability of PHIL setups was analyzed and concluded that the approximation of power amplifier as time delay in s-domain was inappropriate as the power amplification itself impacted the stability of the PHIL setup. A number of techniques to improve the stability of PHIL setups have been proposed, such as the addition of hardware inductance in [26], feedback current filtering in [27], shifting impedance method in [17], multi-rate partitioning in [28], multi-time-step transmission line interface in [29], use of open-loop inverter as power amplifier in [30], with a comparison of a few of these approaches ( [17], [26], [27], [28]) presented in [31].
From a theoretical modelling perspective, there are three differences that can be drawn out between PHIL and GDS setups: i) the GDS setups inherently incorporate larger time delays owing to the geographically dispersed nature of the setup, ii) the GDS setups may or may not incorporate a power amplifier to interconnect the two subsystems, and iii) the interface signals are transformed to DC quantities before their exchange in GDS. The s-domain transfer function models for PHIL setup can be extended for stability analysis of GDS setups -the model readily allows for integration of larger time delays and more often for simplicity the power amplifier is approximated as ideal for the analysis of stability [15], [21], [32]. However, the conventional transfer function models thus far used do not allow for the representation of the transformation of the interface signals. The interface signals in GDS are transformed as the exchange of signals as instantaneous time-domain quantities has been found to be inappropriate for GDS [4]. In addition, the transformation of the signals facilitates the incorporation of time delay compensation during the transformation of the signals back to time-domain quantities, an essential requirement for GDS due to the presence of larger inherent time delays. The transformations of signals has been more recently adopted for PHIL studies as well to realize more accurate setups incorporating time delay compensation [33], [34]. In [3], it was identified that the transformation of the interface signals impacts the stability of the GDS setup under consideration. The conventional simplified transfer function model was extended where the transformations were represented by the equivalent filters adopted for their implementation. Although a step forward, offering valuable insight of the impact of the transformations, further empirical analysis have demonstrated limited suitability of that approach over representative time delays and varied system parameters.
There are currently no representative models available for stability analysis of PHIL simulations that utilize interface signal transformations or GDS. Furthermore, with more than one transformation available for realizing GDS setups (synchronous and asynchronous, explained in more detail in Section II), there is no literature available to offer practical guidance on the choice of transformation for the implementation of GDS setups. To address these gaps, in this paper, the stability and accuracy of GDS setups is investigated to establish its practical applicability.
The key contributions of the paper are summarised below: 1) The inability of existing simplified models to characterize the stability of GDS is demonstrated through a comparison with empirical simulation analysis. With no such analysis reported in literature, this paper is the first in offering the valuable insight that warrants the development of more accurate models for stability analysis of GDS. 2) To remedy the lack of representative models for stability analysis of GDS with interface signals transformation, small signal stability models for two coupling interface transformation approaches, synchronous and asynchronous, have been developed. This will allow for accurate assessment of stability of GDS setups by the wider research community, and the developed models can be directly (or with minor modifications) adopted for accurate stability analysis of PHIL setups. 3) Boundaries for stable operation of GDS for synchronous and asynchronous coupling have been established with respect to impedance ratio, time delay and hardware composition. This will enable prospective GDS user take an informed decision on the choice of coupling transformation for given application at hand. 4) While accuracy of PHIL setups has been used as measure to communicate the confidence in simulations undertaken [35], [36], [37], this paper presents accuracy analysis of synchronous and asynchronous coupling in the common region of their stability. While the choice of coupling is obvious outside the common range of stability, accuracy analysis conducted in this paper with respect to change in voltage and frequency provides practical guidance in choice of coupling. The remainder of the paper is organized as follows: Section II presents the fundamentals of the coupling interfaces. The suitability of conventional stability analysis is evaluated by means of comparison with empirical analysis in Section III. Models for accurate stability analysis for synchronously and asynchronously coupled GDS are derived in Section IV, followed by their validation and practical application discussed in Section V. The accuracy analysis of GDS setups with the two couplings is presented in Section VI followed by two GDS implementations for real-world smart-grid applications discussed in Section VII. Future research directions have been identified in Section VIII while Section IX concludes the paper.

II. COUPLING INTERFACES
The coupling of two subsystems, split for simulation across geographically distributed research infrastructures, is subject to selection of a suitable coupling interface. A number of interface algorithms have been proposed in literature [22], this paper will adopt the voltage type ideal transformer model (V-ITM), where the voltage from Subsystem 1 PCC is reproduced at Subsystem 2 PCC by means of a controlled voltage source, while the currents in response to the reproduced voltage are fed back to Subsystem 1 PCC through a controlled current source. The V-ITM interface algorithm is a popular choice for PHIL and GDS setups owing to its ease of implementation. The choice of interface algorithm dictates the choice of interface signals, as such the chosen V-ITM requires the voltage signals from Subsystem 1 to be transferred to Subsystem 2 and the current signals to be transferred from Subsystem 2 to Subsystem 1. The exchange of interface signals as instantaneous time domain quantities has been shown to be inappropriate for packet based communications over the internet [4]. Therefore, a two step process of interface signal transformation (to DC quantities) before their transfer and signal reconstruction upon arrival is most commonly adopted for GDS implementations. This paper adopts two approaches presented in literature, namely the synchronous coupling approach using direct-quadrature (dq) transformation and the asynchronous coupling approach using root mean square (RMS).

A. Synchronous Coupling
A GDS setup adopting V-ITM interface algorithm and synchronous coupling is shown in Fig. 1(a). In this approach, the three phase quantities (voltages in the feed-forward loop and currents in the feedback loop) are transformed into DC quantities (dq components). These DC quantities are transferred to the corresponding subsystem (to Subsystem 2 in feed-forward loop and to Subsystem 1 in feedback loop), where upon their receipt three phase quantities are reconstructed before further use. In addition to the dq components, the approach requires the transfer of θ = ωt among the subsystems to enable accurate reproduction of the signals with time delay compensation. Therefore, the approach is also referred as synchronous coupling due to the fact that the two subsystems are synchronized with respect to time. The following subsections present the transformations in the feed-forward and feedback loop, followed by the details of time delay compensation employed.
1) Feed-Forward Loop: Voltage signals from Subsystem 1 are transferred to Subsystem 2 in the feed-forward loop. First, the three phase voltages (v a , v b and v c ) are transformed to v α and v β using Clarke transformation, followed by their transformation to v d and v q using Park transformation as [38] v where v d and v q are the direct and quadrature components of the voltage and θ is the phase angle of the voltage measured at Subsystem 1.
The two components with corresponding θ = ωt are sent to Subsystem 2 over the internet using User Datagram Protocol (UDP). Upon receipt of the signals at Subsystem 2, the α, β, and γ components are obtained by applying the inverse Park transformation as, [38] v where θ is the compensated phase angle calculated as is the additional phase angle added during the reconstruction of the signals to compensate for the feed-forward time delay T ff d [33]. More details on time delay compensation presented in Section II-A-3. The α, β, and γ components are then transformed to abc using the inverse Clarke transformation for injection at the PCC of Subsystem 2.
2) Feedback Loop: The current signals from Subsystem 2 are transferred to Subsystem 1 in the feedback loop. The current dq0 components (i d , i q and i 0 ) can be obtained from the instantaneous three phase currents (i a , i b and i c ) by using the Clarke-Park transformation as in the feed-forward loop. Then, upon receipt of the currents at Subsystem 1, the signals are reconstructed in a similar way as in the feed-forward loop. It should be noted that the phase compensation for the current phase angle is calculated by using the time delay of the feedback loop, i.e., T fb d .

3) Time Delay Compensation:
Time delay compensation is paramount to the realization of GDS setups with synchronous coupling. The large inherent time delay due to geographical distance between the two research infrastructures can lead to significant inaccuracies. Time delay remains a concern for the fidelity of PHIL simulations but the typically smaller delay in order of microseconds (compared to milliseconds in GDS) leads to relatively minor inaccuracies.
For an AC signal, the time delay appears as a phase shift and can be compensated by addition of phase equal to the time delay calculated as Δθ = T d · 2πf , where T d is the time delay [33]. Consider Fig. 2, where the voltage and its phase from the two Subsystems in GDS are presented. The Subsystem 2 signals appear to be phase shifted in comparison to the signals from Subsystem 1, and this phase shift is associated to the time delay. This difference is phase causes erroneous calculation of active and reactive power at the PCC leading to the inaccuracies aforementioned. The addition of phase enables reproduction of the signal at Subsystem 2 in such a way that it is in phase with the signal at Subsystem 1 as shown in Fig. 2. This eliminates any discrepancies in calculation of active and reactive powers at the PCC.
Time delay can be calculated and compensated at each Subsystem every time step, i.e., in feed-forward and feedback loop or can be compensated in entirety at one of the Subsystems. The time delay in PHIL systems can be estimated and is typically a constant (minor variability as reported in [39]). However, due to the nondeterministic nature of communications over the Internet, the time delay in GDS can vary considerably. The accurate calculation of the time delay requires GPS (Global Positioning System) clocks at each of the research infrastructures. Each UDP packet is time stamped with time before it is sent to the corresponding research infrastructure. Upon receipt of the UDP packet, the time stamp is compared with current time to determine the time delay. The synchronized measurement ensures accuracy in order of nanoseconds.

B. Asynchronous Coupling
This approach presents a much simpler alternative to synchronous coupling, as shown in Fig. 1(b) and explained in the following subsections.
1) Feed-Forward Loop: In case of asynchronous coupling, the RMS value and frequency of each phase of the three phase voltages are transferred from Subsystem 1 to Subsystem 2 in the feed-forward loop. Defining phase A voltage as v a (t) = V am (t) cos (2πf · t), its RMS value can be calculated as where T w is the window lenwhere T w is the window length of the RMS value computation at time instant t.gth of the RMS value computation at time instant t. Upon receipt of the RMS values and frequencies at Subsystem 2, the AC quantities are reconstructed as where V rms a and f are the RMS voltage and frequency received at Subsystem 2.
2) Feedback Loop: As the two subsystems are not synchronised with respect to time, the currents from Subsystem 2 cannot be utilized for injection at Subsystem 1. Instead, the active and reactive powers measured at Subsystem 2 are transferred to Subsystem 1, where upon their receipt the currents for injection are calculated with reference to the PCC voltage at Subsystem 1 as The three phase currents i a , i b , and i c can be calculated using inverse Park and inverse Clarke transformation. However, it should be noted that the feedback loop of asynchronous coupling utilises the phase angle of the voltage at Subsystem 1 PCC (θ v ) instead of requiring the transfer of the phase angle of currents from Subsystem 2 and does not include the additional phase term to compensate for the time delay.

C. Security of Communications
GDS require the transmission of data between the research infrastructures involved over the Internet. Securing the transmission of data is essential as the experimental data can be sensitive in nature or due to the fact that the data can be driving hardware equipment at the research infrastructures that can be vulnerable to manipulation of data. There are two ways reported in literature in which the risk due to cyber-attacks is mitigated: i) use of VPN (Virtual Private Network) or by utilizing software libraries that offer encryption and ii) adding data check upon receipt of signals. The incorporation of features to secure the data transmission adds additional delays, but yields comparable performance given the large inherent delays involved. A summary of such approaches has been presented in [3].

III. CONVENTIONAL STABILITY ANALYSIS
In this section, the suitability of conventional stability analysis for synchronous and asynchronous GDS is evaluated. First, the model conventionally adopted for stability analysis is presented followed by theoretical and empirical evaluation of stability.

A. Modelling
The conventional approach to modelling a PHIL setup for stability analysis involves its representation as a single input single output (SISO) system. The simplified SISO representation of a PHIL setup is presented in Fig. 3. The fundamental difference between PHIL setup and GDS setup is the lack of power amplifier as interface between the two Subsystems considered and larger time delays due to the fact that the two Subsystems are geographically apart. However, as can be observed from Fig. 3, the conventional modelling approach approximates the power amplifier behaviour as ideal (transfer function representation as unity) and therefore the same SISO is also applicable for the analysis of GDS setups. The open-loop transfer function of GDS setup can be defined as where Z 1 (s) and Z 2 (s) are the equivalent impedances of Subsystem 1 and Subsystem 2 respectively, with the time delay for the feed-forward loop as T ff d and feedback delay as T fb d . Assuming T ff d ≈ T fb d , the total time delay of the setup is considered as 2T d .

B. Theoretical Analysis
The stability of the GDS setups can be evaluated theoretically by means of Nyquist stability criterion where poles of the closedloop system are determined, equivalent to the zeros of the two characteristic equations [28] 1 + G OL (s) = 0 (8) where I is the dimension identity matrix. The impact of impedance ratio, time delay and power factor (PF) is presented in the following sub-sections. 1) Impact of Impedance Ratio Z 1 /Z 2 : The stability assessment for a range of Z 1 /Z 2 (PF = 0.95, T d = 25 ms) is shown in Fig. 4(a). The Z 1 /Z 2 boundary for stability is found as 1 with the system encircling -1 for any value of Z 1 /Z 2 > 1.

2) Impact of Time Delay:
The stability of GDS for T d in range of 10-500 ms (PF = 0.95, Z 1 /Z 2 = 0.9) is shown in Fig. 4(b). From the conventional analysis it can be concluded that the T d does not impact the stability of the GDS setups.
3) Impact of PF: The impact of PF on the stability of GDS setups with Z 1 /Z 2 = 0.4 and Z 1 /Z 2 = 0.9 (T d = 25 ms) is shown in Fig. 4(c) and 4(d). The decrease in PF tends to move the system towards instability, the impact is more significant as the Z 1 /Z 2 approaches its stability boundary.

C. Empirical Analysis
This sub-section presents the empirical evaluation of stability conducted by means of simulations in MATLAB/Simulink. The results of analysis for synchronous and asynchronous coupling are presented in Fig. 5 and Fig. 6, respectively.
The first obvious distinction is the difference in results obtained for synchronous and asynchronous couplings. The simplified model adopted for stability analysis does not cater for the impact of coupling incorporated.
1) Impact of Impedance Ratio Z 1 /Z 2 : It is evident that the impedance ratio boundary for GDS setups is different for the two couplings and much smaller than the value derived from the conventional stability analysis. The Z 1 /Z 2 boundary for synchronous coupling (T d = 25 ms, PF = 0.95) is obtained as 0.2 (as in Fig. 5(a)) and 0.25 for asynchronous coupling (as in Fig. 6(a)).
2) Impact of Time Delay: The results for stability analysis with respect to time delay for synchronous coupling for Z 1 /Z 2 = 0.2 and Z 1 /Z 2 = 0.6 are presented in Fig. 5(b) and 5(c) respectively, while the corresponding results for asynchronous coupling are presented in Fig. 6(b) and 6(c) respectively. In all the cases the PF is equal to 0.95. For Z 1 /Z 2 = 0.2, the stability of synchronous coupling tends to instability with increase in time delay (unstable for T d >25 ms) while the stability of asynchronous coupling is not impacted by the time delay. For Z 1 /Z 2 = 0.6, the stability boundary (with respect to T d ) for synchronous coupling has reduced to 4 ms, while the asynchronous coupling is unstable for any value of T d considered.
3) Impact of PF: The conventional theoretical analysis of impact of PF on stability of GDS setups complies with the empirical analysis. The results for synchronous coupling for Z 1 /Z 2 = 0.05 and Z 1 /Z 2 = 0.2 are presented in Fig. 5(d) and 5(e) respectively, while the corresponding results for asynchronous coupling are presented in Fig. 6(d) and 6(e) respectively. As can be observed, with decrease in PF the system tends towards instability with much significant impact as the Z 1 /Z 2 reaches its stability boundary.

D. Discussion
The theoretical stability analysis of the simplified models identify a larger range of impedance ratio over which the GDS setups are stable, this is in direct contradiction to the empirical results where much smaller impedance ratio boundaries have been identified for both the couplings analysed. The simplified model also fails to identify the impact of time delay on the stability of GDS when not operating at the impedance ratio   boundary of stability. Although the simplified model presents the same trend when analysing the impact of PF on the stability of GDS, the results do not match the empirical results quantitatively. The discrepancy between the theoretical and empirical analysis results is evident, highlighting the inability of the existing simplified model to accurately characterize the stability of GDS employing synchronous or asynchronous coupling. The evaluation of stability of GDS setup before its implementation is imperative to ensure a secure experimental run, particularly when hardware is involved. Erroneous results can lead in expense of time in realizing setups that are infeasible and/or damage to hardware equipment involved. This is the first of its kind evaluation, offering valuable insight that warrants the development of more accurate models that can accurately characterize the stability of GDS setups.

IV. MODELS FOR ACCURATE STABILITY ANALYSIS
In this section, accurate models for a synchronously and asynchronously coupled GDS setups are derived using small-signal stability analysis principle. Any state variable u of a dynamic system can be decomposed into two components, as u =Ū +û (10) whereŪ is the steady-state component andû is the small-signal variation component. Using this assumption, any non-linear system can be linearised around a steady-state operating point for stability analysis using linear systems theory. In the following subsections, open-loop transfer functions for synchronous and asynchronous GDS setups will be developed for accurate stability analysis.

A. Synchronous Coupling Model
The s-domain block diagram of the synchronously coupled GDS setup is illustrated in Fig. 7.
As discussed in section II-A1, the three phase voltages at the PCC in Subsystem 1 are transformed to v α and v β followed by their transformation to v d and v q . By application of small-signal theory as in (10) and neglecting the steady-state components with approximations (sinθ ≈θ, cosθ ≈ 1), the relationships of the small-signal components of voltages between αβ and dq frames can be obtained as wherev α ,V α are the small signal and steady state components of v α respectively,v β andV β are the small signal and steady-state components of v β respectively, andΘ andθ are the steady-state and small-signal components of the phase-angle of the voltage at PCC.
The positive sequence (PS) and negative sequence (NS) components can be obtained by application of Euler's formulas and Laplace transformation to (11) as wherev ± d (s) andv ± q (s) are the small-signal quantities of the PS and NS components of v d and v q in s-domain, respectively.
Upon receipt of the dq components of the voltage at Subsystem 2, v α and v β can be reconstructed using inverse Park transformation (as explained in section II-A1). Employing the small-signal principle and neglecting the steady-state quantity in q frame, the relationship can be expressed as Assuming time delay compensation of the PCC phase-angle is fully achieved, i.e. θ = θ, the α and β components of the voltage at Subsystem 2 can be represented by those received from Subsystem 1 as where T ff d is the time delay for the feed forward User Datagram Protocol (UDP) communication channel. As indicated in [40], the PS and NS components of the phase angle θ can be expressed by αβ components of voltages as where with G P LL as the transfer function of the PLL block. By substituting (16) into (15) and assuming T − P LL = T + P LL = T P LL , θ = 1 2 (θ − + θ + ), the final expression for the transfer function of the feed-forward path, denoted by G dq ff (s), can be obtained as where Same approach can be applied to derive the transfer function of the feedback path. By analysis of the block diagram in Fig. 7, the expressions describing the relationship between the α and β components of the currents at Subsystem 1 and Subsystem 2 are defined by Substitutingv α andv β byî α andî β through power relation, the final equation for the transfer function of the feedback path denoted by G dq fb (s) can be represented as where the coefficients in the equation of G dq fb (s) remain the same as of G dq ff (s). The compact form of the open loop transfer function for the synchronously coupled GDS setup depicted in Fig. 7 can be expressed as

B. Asynchronous Coupling Model
The expression for calculation of the RMS value of the phase A voltage over a moving window of length T w is presented in (4). By differentiating both sides and applying small-signal principle, the equation can be rewritten aṡ v rms Once the RMS values of the voltage arrives at Subsystem 2, the instantaneous voltage signal can be reconstructed as The small-signal quantity of v a (t) is found aŝ  (25) aŝ By substituting (16) into (28) and then the resulted expression into (29), the relationship between the α component of the voltages in the feed-forward path between two subsystems in Fig. 8 can be described aŝ where G 1 (s) = e −sT dV rms T P LL (31) Using the same procedure, the β component of the voltage reproduced at Subsystem 2 can be defined bŷ By using the Clark transformation, the equations for the threephase voltages obtained at Subsystem 2 are defined bŷ A difference between synchronous and asynchronous coupling at this stage is that the latter uses the powers for closing the loop instead of currents. The measured powers can be determined through the three-phase currents and voltages at the Subsystem 2 as where i a , i b , and i c are three-phase currents at Subsystem 2.
Employing small-signal principle for (37) and (38) and combining with (34), (35), and (36) provideŝ whereĪ rms is the steady-state component of the RMS value of the current at Subsystem 2. When the values of the measured powers reach Subsystem 1, the dq components of the current are calculated bŷ Finally, similar to (21), with some approximations, the smallsignal quantities in s-domain of the αβ components of the three-phase currents fed to the current source in Subsystem 1 is computed by using the inverse Park transformation with the phase angle derived from the phase angle of the voltage at Subsystem 1 PCC, i.e. θ in Fig. 8 î α (s) Referring to (16) and (40a,b), it can be observed that the equations of the α and β components of the reconstructed currents in (42) incorporate dynamics of two different subsystems, one coming from the current reconstructed from the powers sent from Subsystem 2 and one from the phase angle of the voltage at Subsystem 1 PCC. The use of (34) to (36) to represent the dynamics of Subsystem 2, i.e.î d (s) in (40a,b), through those of The open-loop transfer function of the asynchronously coupled GDS setup as in Fig. 8 can be represented as

V. DETAILED STABILITY ANALYSIS
The objectives of this section are two-fold: (i) to validate the small-signal models developed in Section IV, and (ii) to determine the boundaries of stability for GDS setups.

A. Small Signal Model Validation and Analysis
In this section, the small-signal models developed for GDS setups with synchronous and asynchronous coupling are validated by means of an extended analysis with respect to impedance ratio, time delay and power factor. 1) Impact of Impedance Ratio: The impedance ratio is varied from 0.1 to 0.4 in steps of 0.05 while the time delay (T d = 25ms) and power factor (PF= 0.95) are kept constant. The Nyquist plots for the two coupling interfaces, synchronous (23) and asynchronous (44) are presented in Fig. 9(a). The general trend of stability is in line with conventional analysis, i.e., with increase in impedance ratio the GDS setups tend towards instability. The stability margins derived from the proposed detailed models is in alignment with the empirical analysis, identified as 0.2 for synchronously coupled GDS setups and 0.25 for asynchronously coupled GDS setups.
2) Impact of Time Delay: The Nyquist plots for varying time delay from 10 ms to 125 ms with Z 1 /Z 2 = 0.15 and PF = 0.95 are shown in Fig. 9(b). In accordance with the empirical analysis, the setups with synchronous and asynchronous couplings are immune to variation in T d for the Z 1 /Z 2 considered. Furthermore, the suitability of the proposed models is also verified for setup with Z 1 /Z 2 = 0.6, where the results obtained are in alignment with those obtained through the empirical analysis, i.e., the stability boundary (with respect to T d ) for synchronous coupling reduced to 4 ms, while the asynchronous coupling is unstable for any value of T d considered (Fig. 5(c) and Fig. 6(c) vs Fig. 9(c)).
3) Sensitivity to Power Factor: The analysis from the developed models complies with the conventional and empirical analysis presented in Section III. The results for GDS setup with synchronous coupling for Z 1 /Z 2 equal to 0.15 and 0.2 and asynchronous coupling for Z 1 /Z 2 equal to 0.2 and 0.25 are presented in Figs. 9(d) and 9(e). As can be observed, with decrease in PF the system tends towards instability with much noticeable impact as the Z 1 /Z 2 reaches its stability boundary.

B. Boundary Conditions and Practical Applications
In this subsection, the stability boundaries for synchronously and asynchronously coupled GDS setups are derived. With the general trend of impact of impedance ratio, time delay and power factor of the system under consideration on stability presented in previous subsection, here a comprehensive stability analysis within a practical range for identified parameters has been undertaken: Z 1 /Z 2 = [0.1, 1], T d = [500 μs, 0.5 s] and PF= [0.5, 1]. The derived boundary conditions are presented in Fig. 10, with Fig 10(a) presenting boundaries for synchronously coupled GDS setups while Fig. 10(b) presents boundaries for asynchronously coupled GDS setups. The y-axis of the plot represents the impedance ratio plotted against logarithmic scale of time delay on x-axis. Each curve represents a given PF, with area under the curve representing stable region of operation. From the analysis it can be concluded that 1) For synchronously coupled GDS setups, the stability is dominated by time delay, i.e., the time delay of a given setup determines the impedance ratio of systems that can be studied using GDS. A setup with time delay of 500 μs can be utilized for studies of systems with impedance ratio less than 0.95, however, this falls drastically to an impedance ratio of ∼ 0.2 for any time delay greater than ∼ 20 ms.
2) The stability of asynchronously coupled GDS is dominated by impedance ratio of the setups and is not impacted by time delay. The impedance ratio limit for asynchronously coupled GDS is identified as ∼ 0.3. 3) A decreasing power factor generally reduces the impedance ratio limit for any given impedance ratio and time delay (compared to unity power factor). From an application point of view, a synchronously coupled GDS setup is an obvious choice for systems with impedance ratio greater than 0.3, while respecting corresponding time delay and power factor boundaries depicted. For setups with impedance ratio between 0.15 and 0.3, the choice of coupling needs careful consideration with respect to PF and time delay of the setup. For setups with impedance ratio less than 0.1, both synchronous and asynchronous coupling present good stability for the practical range of PF and time delay considered, a broader consideration of the accuracy of the setup should be in scope.

VI. ACCURACY ANALYSIS
In this section, the accuracy of interface coupling for GDS setups at the PCC is evaluated. Choosing appropriate metrics for comparison is essential. Two metrics have been defined to quantify the accuracy of the interface coupling as r Maximum Instantaneous Power Tracking Error (MIPTE) is defined as the maximum excursion of the measured apparent power S (t) of the GDS setup compared to the measured apparent power S (t) of the monolithic setup.
r Cumulative Power Tracking Error (CPTE) is defined as the sum of the tracking errors at every time step T s from the application of disturbance to the time when measured apparent power is restored within the tolerance band, i.e., for period of T r .
where N = T r /T s , with T s as the time step. T r is the error restoration time defined as the time elapsed from the disturbance initiation to when the measured apparent power returns and stays within the defined tolerance band .
T r = argmin{T r ∈ R | ∀t > T r : + < S (t) < − } (47) A smaller CPTE corresponds to better set point tracking capability. In comparison to metrics in [37], the chosen metrics present a single value (analogous to cross-sectional evaluation as opposed to time-series) encapsulating the performance of the chosen interface. The MIPTE provides a measure of maximum instantaneous error of a given interface while CPTE captures its average performance over a given period of time.
The accuracy will be evaluated with respect to the same three parameters: time delay, impedance ratio and power factor. To allow for a fair comparative accuracy analysis of the two couplings, the range of parameters have been chosen such that both coupling interfaces are stable, obtained from The remainder of this section will present the results of the comparison with respect to two disturbances, one step down in frequency (50 Hz to 49 Hz) followed by a step up in voltage (400 V to 440 V).
The results for the accuracy analysis for step change in frequency are presented in Fig. 11, with MIPTE for variation in parameters in Fig. 11(a)-11(c) and CPTE for variation in parameters in Fig. 11(d)-11(f). As can be observed, the synchronous coupling exhibits an approximately linear trend where  the accuracy deteriorates with increase in T d , Z 1 /Z 2 and decrease in PF. For all the variations in parameters considered, the asynchronous coupling presents better accuracy performance with minimal variation in accuracy with respect to the variations in parameters considered.
Likewise, the calculations of MIPTE and CPTE, which are obtained from step change in voltage, are shown in Fig. 12(a)-12(c) and Fig. 12(d)-12(f), respectively. Generally it can be stated that synchronous coupling exhibits better accuracy when compared to asynchronous coupling for a step change in voltage. MIPTE is minimally impacted by T d , higher Z 1 /Z 2 exhibits lower MIPTE while lower PF exhibits higher MIPTE. Although the MIPTE is minimally impacted with variation in T d , it can be observed that CPTE increases with increase in time delay (exhibiting exponential increase for the boundary of observation T d = 100 ms and Z 1 /Z 2 = 0.15. Similar to the analysis for step change in frequency, the accuracy for both the couplings deteriorates with decrease in PF.

VII. CASE STUDIES
This section presents two case studies where GDS implementations for real-world smart grid applications have been realized. The use of the analysis conducted and presented within the previous subsections guides the choice of interface coupling, demonstrating support in applicability and implementation of GDS setups between research infrastructures.

A. Voltage Support in a Distribution Network
This section presents an example of GDS implementation for voltage control study within a distribution network. The distribution network under consideration is shown in Fig. 13, and has been adapted based on [41]. Each bus within the network comprises a mix of loads and Photovoltaic (PV) systems, represented in aggregation for simplicity. The objective of the study is to understand the capability of the battery energy storage system (BESS) to respond to voltage variations and thereby supporting voltage regulation of the distribution network.
The two laboratories in consideration for the experiment are the Dynamic Power systems Laboratory (DPSL) at the University of Strathclyde (UoS) and the Electric Energy System Laboratory (EESL) at the National Technical University of Athens (NTUA), with a geographical distance of ∼3,500 km between the two. The average time delay between the two laboratories is ∼60 ms as shown in Fig. 14 (top), with maximum delay observed as ∼75 ms. The next step in implementation requires the evaluation of the equivalent impedance ratio at the PCC. The BESS can be integrated at any of the four buses within the distribution network, given the objective of testing the capability of the BESS to respond to variations in voltage. With the inherent time delay and impedance ratio at hand, an informed decision on choice of coupling can be made. The synchronous coupling is on the verge of stability for delays greater than ∼10 ms for Z 1 /Z 2 = 0.2, and unstable for Z 1 /Z 2 = 0.3 and Z 1 /Z 2 = 0.4 (as in Fig. 10). The stability of asynchronous coupling is not dependant upon the time delay, while the coupling can remain stable for Z 1 /Z 2 = 0.2 and Z 1 /Z 2 = 0.3. While the synchronous coupling offers better accuracy for changes in voltage (as shown in Fig. 12), the asynchronous coupling offers wider range of scenarios for evaluation and therefore is a better choice for GDS implementation.
As shown in Fig. 13, the distribution network is simulated in a DRTS at DPSL (Subsystem 1) while the components at the bus are simulated in DRTS at EESL (Subsystem 2). With asynchronous coupling adopted for implementation, the rms voltage and frequency are sent from DPSL to EESL, while the active and reactive power in response to the voltage is sent from EESL to DPSL. To test the response of the BESS, a 10 % voltage step (reduction from 1 pu to 0.9 pu as in Fig. 13 (middle)) is introduced within Subsystem 1 by means of the on-load tap changer at the grid coupling point. BESS is rated at 100 kW and is designed to respond to variations in voltage by injecting active power. The combined response of the Subsystem 2 is presented in Fig. 13 (bottom). The response emulating three different impedance ratios has been presented, practically representative of different locations within the distribution network.
The transmission rate chosen for the implementation was varied with reported results from an exchange of interface signals at 100 Hz. As is evident, a successful implementation has been realized at much lower transmission rate than required for synchronously coupled GDS setups. This presents an additional advantage that can lead to increased adoption of GDS in circumstances where hardware capabilities are limited.

B. Frequency Support in a Transmission Network
This section presents an example GDS implementation for frequency control within a transmission network. A reduced six bus dynamic model of the Great Britain (GB) power network is considered, where each bus comprises an aggregated load and an aggregated generation unit. The generation units are synchronous generators modeled to reproduce the representative dynamics of GB power system. The objective of the study is the integration of inertial response (IR) through controllable loads within the network to support the frequency regulation.
In addition to the aforementioned laboratories (DPSL and EESL), the Real-Time Laboratory (RTLab) at RWTH Aachen University, Germany, is considered for the frequency control study. RTLab and DPSL are separated over a geographical distance of ∼1250 km. The average time delay between the two research infrastructures is ∼16 ms, with maximum observed delay of ∼28 ms as shown in Fig. 16 (top). An empirical evaluation has identified an impedance ratio Z 1 /Z 2 = 0.15 when the network is split at the border of Scotland and England (Subsystem 1 -two buses in the north representing Scottish load centres, Subsystem 2 -remainder of the four buses representing the English load centres as in Fig. 15). While synchronous and  asynchronous coupling offer stability for the impedance ratio identified under the two time delays, i.e., ∼60 ms average and ∼75 ms maximum for EESL and ∼16 ms average and ∼28 ms maximum for RTLab, the evaluation in Section VI reveals better accuracy for synchronous coupling with DPSL-RTLab time delays.
The synhronous GDS setup for the frequency support study of transmission network is shown in Fig. 15. As detailed in Section II-A, the time stamped dq components of the voltage, along with phase are sent from DPSL to RTLab forming the feedforward loop while the time-stamped dq components and phase of the currents are sent back from RTLab to DPSL. To evaluate the impact of provision of IR through controllable loads on the frequency of the network, a 1GW disturbance (loss of generator emulated by increase in load) is introduced in Subsystem 1. About 20 % of the load demand is assumed to be fast acting and controllable. Inertial response in fast acting controllable loads is incorporated based on swing equation referred to as swing equation based inertial response (SEBIR) [42]. The response of the network to the disturbance with and without IR is presented in Fig. 16 (bottom). Three different levels of inertial support are considered, i.e., H = 5 H = 7.5 s and H = 10 s. The frequency response of the network improves with incorporation of IR, with better response for higher value of H. This demonstrates the successful use of synchronous GDS for frequency control study within a transmission network.
An additional simulation with different rates of data exchange is undertaken to assess the impact of transmission rate. The results for four different transmission rates (1 kHz, 2 kHz, 3 kHz and 4 kHz) have been presented in Fig. 17. As can be observed, the rate of transmission does not impact the time delay. The jitter and variability in the delay remains in the same order of magnitude. With appropriate time delay compensation implemented, the errors in frequency response are negligible.

VIII. FUTURE OUTLOOK
While the detailed models developed in this paper and the consequent establishment of the boundaries of stability will further the understanding of the applicability of GDS and support safe implementations, there are still areas that present exciting directions for further exploration to support wide-scale adoption of the concept: 1) As the concept of GDS becomes established, facilitated by adoption by research laboratories across the world, it becomes essential to maintain a database of such infrastructures and the services they offer. Furthermore, the business model or framework for collaborative use of the services needs to be explored. Most work to date is through bilateral research collaboration between research infrastructures, however commercialization opportunities can be explored. 2) This work has presented the analysis of accuracy at a single PCC, i.e., in GDS implementation where the power network is split among two research infrastructures. For large scale implementations where the power network is split across more than two infrastructures, the transference of inaccuracies from one subsystem to the other requires further consideration and evaluation. 3) Most research infrastructures undertaking GDS have reported custom implementations for the establishment of communications. A VPN tunnel is established between the two research infrastructures to facilitate secure transmission of data. A more standard approach to realization of communications is required to bolster wide scale adoption.
One option for such secure and standardized approach is the development of service through cloud solutions platforms such as Amazon Web Services or Microsoft Azure.

IX. CONCLUSION
This paper presents small-signal stability models of GDS setups with two transformations, one for synchronous coupling (dq0) and one for asynchronous coupling (RMS). The derived models account for the transformation of interface signals, in addition to typical parameters of interest for stability analysis of geographically distributed simulations (GDS) setups. The proposed models have been verified through empirical simulation based analysis. Furthermore, boundaries of stable operation for synchronous and asynchronous GDS setups have been identified to establish their practical applicability. The boundaries are supplemented with accuracy analysis to enable informed decision on choice of coupling. Two case studies of GDS implementation for real-world smart grid applications have been presented, an asynchornous GDS for voltage control of dinstribution network between the Dynamic Power Systems Laboratory (DPSL) at the University of Strathclyde (Scotland) and the Electric Energy Systems Laboratory at the National Technical University of Athens (Greece) and a synchronous GDS between DPSL and Real-Time Laboratory at RWTH Aachen (Germany). This thorough analysis will aid in understanding of the applicability of GDS and support its wide scale adoption that will enable the evaluation and subsequent realisation of experimental setups capable of encompassing the growing complexity of power systems. Future research directions to cater and support the growing interest in secure adoption of GDS have been identified.