Cascading Failures and Resilience in Interdependent Critical Infrastructures: A Dynamic Approach

1 Cascading Failures and Resilience in Interdependent Critical Infrastructures: A Dynamic Approach Balaji V. Venkatasubramanian, Senior Member, IEEE, Christos Laoudias, Member, IEEE, Mathaios Panteli, Senior Member, IEEE, Abstract—Natural hazards, such as storms, increasingly threaten interconnected critical infrastructure systems (CIS), necessitating integrated resilience strategies that account for interdependencies, cascading failures, and dynamic recovery. To contextualize the research within existing developments, a bibliometric analysis is conducted, revealing key trends and gaps in hazard-resilient CIS modeling. This paper proposes a comprehensive framework that models interdependent behavior across power, telecommunication, and transportation networks under evolving storm conditions. A dedicated hazard engine captures the spatiotemporal progression of storms and simulates their direct and cascading impacts across infrastructure layers. The framework incorporates a dynamic, safety-aware restoration strategy that deploys repair crews and Mobile Energy Resources (MERs) in alignment with real-time storm conditions. The proposed framework is demonstrated through a case study inspired by the realistic infrastructure of Cyprus. Results show the model’s ability to capture cascading disruptions, illustrate the operational effectiveness of dynamic restoration, and provide insights from sensitivity and trade-off analyses. Index Terms—Cascading Failures, Dynamic Restoration, Interdependent Critical Infrastructures, Natural Hazards, Resilience. I. INTRODUCTION CRITICAL infrastructure systems (CIS)–including power, communication, transportation, and water–support essential societal functions and increasingly operate as interconnected System of Systems (SoS) [1]. Disruptions in one sector can cascade into others, amplifying impacts. Recent events such as Hurricane Beryl (2024) [2] and the 2025 blackout in Spain and Portugal [3] underscore these vulnerabilities. Addressing them requires integrated, SoS-based resilience strategies rather than isolated risk assessments. A. Research Evolution: Insights from Bibliometric Analysis This study analyzes over 800 publications retrieved from Scopus and Web of Science (2001–20251) to examine the evolution of SoS research on CIS, including power, communication, transportation, and water networks. Following the methodology in [4], titles, abstracts, and keywords were This work was supported by the European Union’s Horizon Europe research and innovation program through the project Solid Preparedness And Resilience for Robust Operations during disaster Wilderness (SPARROW) (Grant agreement ID: 101168499). Balaji V. Venkatasubramanian (email:[email protected]) and Christos Laoudias (email:[email protected]) are with KIOS Research and Innovation Center of Excellence, University of Cyprus, Cyprus. Mathaios Panteli (email:[email protected]) is with KIOS Research and Innovation Center of Excellence, Department of Electrical and Computer Engineering, University of Cyprus, Cyprus. 1For 2025, only research up to April is included. Fig. 1. Annual Research Publications related to Interdependent Critical Infrastructures within System of Systems Framework from 2001 to 2025 systematically analyzed to identify thematic evolution. Publications were collected using a comprehensive search query2 designed to capture interdependencies, cascading effects, and cross-sector interactions. CIS-related publications have grown steadily over the past two decades, particularly after 2015 (Fig.1), reflecting both stronger interconnectivity among infrastructures and increased exposure to hazards. To examine how research themes evolved, a keyword-based thematic analysis was conducted using Bibliometrix [5]. The resulting Sankey diagram (Fig. 2) highlights the dominant themes across four periods (2001–2025), illustrating a shift from conceptual modeling to hazard-driven resilience. 1) Establishing Foundations for CIS Interdependencies: In Period A (2001–2010), research emphasized CIS interdependencies through network modeling and optimization, often framed around topology and decision-making under uncertainty. Early studies linked power, telecommunication, and water systems to assess failure propagation [6], examined the computational complexity of interdependency modeling [7], and introduced multi-agent approaches [8]. While foundational, these contributions lacked hazard-oriented modeling and recovery considerations. 2(“ interdependent AND networks* OR interconnected* OR cascading OR “system of systems” OR “critical infrastructure”) AND (( power OR electricity OR “energy network” OR “electric grid” OR “power grid” AND ( communication OR “telecom network” OR ict OR “information and communication” )) OR ( power OR electricity OR “energy network” OR “electric grid” OR “power grid” AND ( water OR “water supply” OR “water network” OR “water infrastructure” )) OR ( power OR electricity OR “energy network” OR “electric grid” OR “power grid” AND ( transportation OR traffic OR “transport network” OR “transport infrastructure” )) OR ( communication OR “telecom network” OR ict OR “information and communication” AND ( water OR “water supply” OR “water network” OR “water infrastructure” )) OR ( communication OR “telecom network” OR ict OR “information and communication” AND ( transportation OR traffic OR “transport network” OR “transport infrastructure” )) OR ( water OR “water supply” OR “water network” OR “water infrastructure” AND ( transportation OR traffic OR “transport network” OR “transport infrastructure” ))). 2 Fig. 2. Thematic Evolution of Authors’ Keywords from 2001 to 2025. Flows represent transitions of key themes between periods, with thickness proportional to the number of linked studies. 2) Power-Focused Research in CIS Interdependencies: In Period B (2011–2020), research became increasingly powercentric, with keywords such as ”critical infrastructures” and ”power transmission networks” dominating. During this period, [9] examined cascading failures between power and water systems, while [10] showed how cyberattacks on communication systems affecting smart grid, and [11] introduced a multiagent framework linking energy, transport, and communication, though the analysis remained energy-centric. Although interdependencies were more explicitly recognized, most studies remained limited to power resilience, with minimal attention to cross-sector interactions or hazard-driven restoration. 3) Resilience and Initial Application of AI/ML: Recent research during Period C (2021–2024) broadened to include resilience and cross-sector dependencies. Studies addressed hurricane impacts on power–transportation [12], coordination between power and water under energy constraints [13], and governance barriers in integrating power and communication [14]. In parallel, AI/ML approaches began to emerge, mainly for vulnerability assessment and prediction; for example, [15] applied machine learning (ML) to estimate restoration times in interdependent networks, representing one of the first systematic applications of ML to CIS recovery. 4) Recent Progress in CIS Interdependencies: The most recent phase (Period D) shows a stronger hazard-driven orientation, with recurring themes of cascading failures, disasters, and resilience. Examples include earthquake-induced water–power cascades [16], power–communication failures emphasizing coupling strategies and post-disaster repair [17], and recovery models integrating community demographics [18]. ML applications also expanded, such as heuristic and graph neural network methods for resilient communication networks [19], while agent-based studies examined coordination schemes for resilience, though limited to road–water networks [20]. B. Methodological Coverage: Insights from Keyword Analysis While the thematic evolution highlights broad conceptual and sectoral trends, a complementary keyword-based analysis of over 800 abstracts examines methodological coverage. Keywords represent sector-specific tools (e.g., MATPOWER for power, SUMO for transport, EPANET for water, NS3 for telecom), restoration (post-disaster, recovery), hazard modeling (fragility curves, flood modeling), advanced approaches such as AI/ML (deep learning, reinforcement learning, clustering), and multi-agent systems (ABM, MAS). The analysis classifies papers into restoration (RS), system-specific modeling (SSM), hazard modeling (HM), or their joint occurrence (RS+SSM+HM), with RS mainly heuristic/optimizationbased, HM relying on fragility or random-failure functions, and SSM using sectoral simulators. Results are summarized in Table I, which reports coverage across infrastructure combinations. A tick mark indicates that at least one study in that category employed AI/ML or multi-agent methods. Table I shows research focuses on pairwise systems, with little on trior four-sector cases. Multi-agent methods appear sparsely since 2009, while AI/ML emerges mainly in Period C (2021–2024), both limited to pairwise combinations [15], [19]. No study jointly addresses RS, SSM, and HM across multiple infrastructures, underscoring a critical gap. C. Research Gaps and Contributions The combined thematic evolution, keyword analysis, and literature review highlight three major gaps in current CIS research. First, higher-order interdependencies involving three or more infrastructures remain largely unexplored. Second, methodological integration is limited: while RS and SSM are frequently studied in pairwise systems, HM is rarely included beyond two-sector settings, and no study simultaneously addresses all three across multiple infrastructures (Table I). Third, advanced methods show only partial uptake: multiagent approaches have appeared intermittently since 2009, and AI/ML mainly in Period C, yet both are largely confined to pairwise systems, with hazard-informed integration observed primarily in power–communication studies. Complementing these efforts, this paper proposes an integrated framework to simulate cascading failures and dynamic restoration across interdependent power, communication, and transportation networks. The framework is extensible to additional infrastructures (e.g., water), with AI/ML or multi-agent integration identified as a potential direction for future enhancement. The major contributions of this paper are: •Interdependency Modeling Framework: A detailed modeling framework that captures the spatial and functional TABLE I RESEARCH COVERAGE ACROSS INTERDEPENDENT CIS. P = POWER,C= COMMUNICATION,T=TRANSPORT,W=WATER. EACH ENTRY REPORTS THE STUDY COUNT WITH [AI/ML, MULTI-AGENT]EVIDENCE. CIS Combination RS SSM HM RS+SSM+HM P+C+T+W 0 [×,×] 1 [×,×] 0 [×,×] 0 [×,×] P+C+T 5 [×,×] 11 [×,×] 1 [×,×] 0 [×,×] P+C+W 0 [×,×] 0 [×,×] 0 [×,×] 0 [×,×] P+T+W 3 [×,×] 0 [×,×] 1 [×,×] 0 [×,×] C+T+W 1 [×,×] 0 [×,×] 0 [×,×] 0 [×,×] P+C 17 [✓,✓] 114 [✓,✓] 14 [×,✓] 1 [×,×] P+T 25 [✓,✓] 10 [×,✓] 2 [×,×] 0 [×,×] P+W 5 [×,×] 1 [×,×] 1 [×,×] 0 [×,×] C+T 5 [×,✓] 5 [×,×] 1 [×,×] 0 [×,×] T+W 2 [×,×] 0 [×,×] 1 [×,×] 0 [×,×] 3 interdependencies among power, telecommunication, and transportation systems, enabling analysis of cascading failures across sectors. •Spatiotemporal Hazard Engine: A novel engine that simulates hazard progression and its cascading impacts on interconnected infrastructure systems over space and time, allowing for more realistic, scenario-driven resilience analysis. •Dynamic and Safety-Aware Restoration Strategy: An adaptive restoration approach that dynamically deploys repair crews (for both power lines and roads) and Mobile Energy Resources (MERs) in response to evolving hazard conditions, ensuring that recovery begins only in safe areas while prioritizing responder safety. The remainder of the paper is organized as follows. Section II introduces the integrated framework and its core components: hazard modeling, interdependency analysis, and dynamic restoration. Section III presents the case study, simulation setup, and results. Section IV concludes with key insights and future directions. II. SPATIOTEMPORAL MODELING AND DYNAMIC RESILIENCE ASSESSMENT OF INTERDEPENDENT CIS This section presents the integrated framework for modeling cascading failures and assessing resilience across interdependent CIS, structured as shown in Fig. 3. The process begins with data acquisition and preparation, including infrastructure networks, spatial and terrain data, population statistics, and historical hazard records. These datasets are cleaned and spatially aligned to ensure consistency across layers. Interdependencies among power, telecommunication, and transportation systems are then modeled through structural coupling and dependency mapping. The spatiotemporal hazard engine simulates storm progression and triggers failure updates across the coupled network. Restoration and recovery modeling follows, where power line crews, road repair crews, and MERs are deployed through a dynamic restoration strategy (DRS) that adapts to evolving hazard conditions, ensuring safety-aware task execution and timely system state updates. In this context, dynamic refers to a continuous feedback process linking hazard evolution, infrastructure degradation, interdependencies, and restoration actions, such that damage propagation and recovery decisions are state-aware and updated over time. Finally, key resilience metrics are computed to evaluate performance during disruption and recovery, capturing service continuity, infrastructure failures, and societal impacts, with detailed formulations provided in the metrics subsection. A. Data Acquisition and Preparation The integrated framework begins with the acquisition and processing of spatial datasets required for hazard simulation and resilience assessment of interdependent CIS. These datasets include infrastructure topologies, geospatial boundaries, demographic information (e.g., population), terrain characteristics, vegetation distributions (e.g., tree data), and historical hazards. Each dataset type and its processing steps are described in the following subsections. Data Acquisition & Preparation Historical Hazard Events Data Cleaning Spatial Alignment Modeling Interdependencies Across Infrastructure Layers Structural Coupling, Dependency Mapping, and Layer Interactions Infrastructure Network Data Spatial Boundaries, Terrian Data, Population Statistics Spatiotemporal Hazard Engine Restoration and Recovery Modeling Key Resilience Metrics Crew and Mobile Energy Resource Modeling, Initialization, Repair Scheduling, Dynamic Routing, System State Update. Dynamic windstorm propagation across acoupled network and corresponding failure updates. DNS Assets Failed Towers Failed Population Affected Road Blockage Count System-wide Recovery Veg et at ion D at a Fig. 3. Integrated Framework for Modeling Cascading Failures and Resilience Assessment across Interdependent CIS. 1) Infrastructure Network Data: The infrastructure network data includes power network, telecommunication towers, and road network. The power network captures the topology of the grid, encompassing substations and structural components such as towers and distribution poles, depending on the voltage level and asset type. Telecommunication tower data specify geographic location, structural attributes (e.g., height), and services offered (e.g., GSM, UMTS, LTE, 5G). The road network is represented by geographic line strings describing individual road segments. All infrastructure assets are modeled as a spatial graph: key facilities (e.g., substations, towers) are represented as Point geometries, while connections (e.g., power lines, road segments) are represented as LineString geometries. These follow standard GIS conventions, commonly used in formats such as shapefiles and GeoJSON. 2) Spatial Boundaries and Terrain Data: Spatial boundary data include administrative units (e.g., districts or municipalities) within the study area. In parallel, a simulation extent is established to constrain hazard modeling. This extent is derived by computing the bounding box of the combined CIS, with an added buffer (e.g., 0.05◦≈5km) to ensure complete storm coverage. The simulation is constrained to this extent and terminates once the storm exits, enabling focused and computationally efficient hazard analysis. Terrain data is obtained from a digital elevation model (DEM) raster, providing elevation values critical for modeling storm propagation dynamics. The DEM is processed to compute slope and aspect using Sobel filters [21], influencing storm wind speed, radius, and heading deflections. These datasets are standardized to a common coordinate reference system, typically EPSG:4326 (WGS 84), to ensure spatial consistency across analyses. 3) Vegetation and Tree-specific Data: Vegetation data play a critical role in modeling road blockages caused by falling trees. To capture trees adjacent to roads, Sentinel-2 satellite imagery is used to derive the Normalized Difference Vegetation Index (NDVI) [22], which is combined with canopy height measurements from the GEDI L2A lidar dataset [23]. Treespecific data are extracted by cross-referencing GEDI canopy height points with NDVI values, retaining only points in vegetated areas with NDVI greater than 0.3. This filtering isolates canopy heights likely to represent trees. To spatially associate 4 these tree points with the road network, a spatial join operation is performed. For accuracy, both road geometries and GEDI point locations are projected into a local UTM coordinate system (e.g., EPSG:32636 for Cyprus). Each road segment is then buffered by a user-defined distance (meters), producing a polygon that captures nearby GEDI points. Within this buffer, the number of tree points is counted to represent tree density for that segment, while RH98 values (98th percentile relative height) are averaged to estimate typical tree height. Segments without overlapping GEDI points are assigned zero counts and heights, ensuring only segments adjacent to vegetated areas are considered at risk of blockage. The resulting dataset combines tree count and average height for each road segment, providing critical inputs for estimating fragility-based probabilities of road blockage under windstorm conditions. 4) Hazard Engine Parameters: The proposed framework generates synthetic storm scenarios calibrated with parameters—such as wind gust speed and storm radius—using historical data as guidance. In other words, local or regional wind climatologies [24] are used to establish baseline envelopes for calibration, while extended ranges are included as stresstest scenarios to assess resilience under more extreme but plausible conditions. Additionally, the hazard engine requires fragility curves for electrical lines and trees to represent their vulnerability to windstorms, which are adopted from published models [25], [26] and applied uniformly in the absence of region-specific calibration data. 5) Data Cleaning and Spatial Alignment: During the data cleaning stage, all spatial datasets—including telecommunication tower locations, power infrastructure elements, population distribution points, and administrative boundaries— are reviewed, corrected, and standardized for consistency in format, attributes, and spatial reference. To ensure spatial compatibility and support accurate integration, all datasets are re-projected to a common coordinate reference system, enabling cross-layer analyses without geometric distortion. B. Modeling Interdependencies Across Infrastructure Layers The resilience of interconnected CIS hinges on accurately representing the relationships across systems such as power, transportation, and telecommunications. To model these interdependencies, the framework adopts a multi-layer network structure that integrates both physical and logical connections. The following subsections detail the methods used to represent structural coupling, define cross-layer dependencies, and characterize their interactions. 1) Structural Coupling Between Infrastructure Layers: The framework uses spatially explicit methods to create topological links between infrastructure layers based on proximity and connectivity. It models two key couplings: power–telecommunication and power–transportation, establishing physical interdependencies among all CIS. A hybrid spatial-graph method links substations to telecommunication towers. Each substation associates with its nearest distribution pole, and paths from towers to substations are computed. If a path falls within a predefined threshold (e.g., 5 units), a connection is established; otherwise, the Euclidean distance Fig. 4. Conceptual Illustration of Hybrid Spatial Coupling between an Electrical Substation and Telecommunication Towers and Nearest Road. is compared to the same threshold. Nearest-neighbor queries are performed using SciPy’s cKDTree [27]. As depicted in Fig. 4, Tower A establishes a graph-based path within the 5 units threshold, Tower B meets the fallback criterion with a 3.54 units Euclidean distance, and Tower C, though having a graph-path, is excluded due to exceeding the threshold. These examples illustrate the significance of distance thresholds in determining connection establishment. Note that the 5-unit threshold included here is purely for conceptual illustration; the actual analysis applies a 2000 m threshold (see Section III). For power-transportation coupling, each substation establishes connections with its nearest road segment, facilitating physical access for restoration, crew routing, and mobility analysis. As shown in Fig. 4, the substation connected to the road network via a road access link, enabling physical access for restoration. 2) Dependency Mapping and Layer Interactions: Beyond structural coupling, dependency rules define how failures and constraints propagate across layers. Telecommunication towers are linked to designated substations. When a substation is deenergized, the associated tower switches to backup power and becomes non-operational if the outage exceeds this duration. For non-operational towers, the communication coverage radius is computed using frequency-specific signal-propagation models, and overlaid with population data to estimate affected users. This provides a spatially explicit representation of service availability during disasters, moving beyond a purely binary outage. Graded service degradation is beyond the present scope but noted for future work. The road network imposes further physical access constraints, as blocked segments delay power-line restoration and MER deployment until cleared. These dependencies capture the most immediate and widely observed storm-driven pathways that restrict restoration, while the framework remains sufficiently flexible to accommodate additional interdependencies among infrastructures in future extensions. It is noted that telecommunication towers are not modeled as directly failing from windstorm hazards or requiring dedicated repair crews; their functionality is restored once power supply to the associated substation is re-established. 5 C. Spatiotemporal Hazard Engine The spatiotemporal hazard engine simulates windstorminduced disruptions to CIS, focusing on direct impacts to power and road networks, as well as cascading effects from power outages to communication networks. By integrating geospatial datasets, terrain characteristics, tree-specific data, hazard characteristics, and fragility curves for power lines and trees, the engine generates realistic hazard scenarios capturing the dynamic interplay between storm evolution and infrastructure vulnerabilities. It is capable of generating multiple storm scenarios, each characterized by a unique storm path, wind speed profile, and impact on CIS components. The output from this engine includes time-varying matrices of CIS status, wind speeds, and failure probabilities, enabling detailed analysis of cascading impacts over a specified time horizon. 1) Windstorm Propagation Modeling: The windstorm propagation model simulates the dynamic movement of a storm across a geographic region, explicitly accounting for terrain interactions to evaluate potential impacts on CIS, such as power and transportation networks. The storm generation procedure is summarized in Algorithm 1, which initializes a storm outside the study region’s bounding box with a random direction and entry/exit points to ensure a natural approach. The path distance is divided into discrete time steps over a horizon T, and storm positions are updated using geodetic calculations that account for Earth’s curvature [28]. At each step, the algorithm retrieves elevation (et) and slope (st) from a DEM, with slope computed using Sobel filters as st= ∆e/δ. These terrain features drive the dynamic adjustment of storm parameters: wind speed (wt), radius (rt), and heading (ht). Storm dynamics are adjusted sequentially. First, wind speed is updated according to (1), where uphill slopes (st>0) reduce intensity due to drag, downhill slopes (st<0) enhance it through gravitational acceleration, and higher elevations weaken the storm due to increased surface roughness and atmospheric effects [29]. Next, the storm radius is updated by (2), expanding under stronger winds and contracting at higher elevations, reflecting the empirical coupling between storm size and intensity reported in tropical cyclone studies [28]. Finally, heading deflection is determined by (3), with variance that increases proportionally to slope magnitude, thereby capturing terrain-induced steering effects [29]. wt+1 =wt·1−0.02 max(0, st)+0.015 max(0,−st) −0.01et+1 1000 ·exp(−0.001t). (1) rt+1 =rt·1+0.01 wt−wmin wmax−wmin −0.015et+1 1200 ·exp(−0.0005t). (2) ht+1 =ht+δt, δt∼ N(0,(2 + 10|st|)◦).(3) These heuristic formulations, informed by topographic wind–terrain interaction studies, are implemented within Algorithm 1. The simulation terminates once the storm exits the region, yielding a spatiotemporal dataset of storm position, elevation, wind speed, radius, and heading for impact analysis. Algorithm 1 Windstorm Generation Require: Boundary box, time-horizon T, wind speed range, terrain data Ensure: Storm path data 1: Initialize geodetic and terrain tools 2: Select random direction; set start/end points outside boundaries 3: Compute path distance D, step size δ←D/T 4: Set initial heading h0, radius r0, wind speed w0, position, path S 5: for t= 0 to T−1do 6: Compute next position using htand δ 7: Get elevation et+1 and slope stfrom terrain analyzer 8: Update wt+1,rt+1 based on et+1 and st(see (1), (2)) 9: Deflect ht+1 randomly based on |st|(see (3)) 10: Constrain wt+1, rt+1 within bounds; apply decay 11: if outside boundaries then 12: break 13: end if 14: Append (position, et+1,wt+1,rt+1,ht+1)toS 15: Update position 16: end for 17: return Sdf (structured S) 2) Fragility-Based Failure Modeling: The fragility-based failure model quantifies the vulnerability of CIS to windstorm impacts by integrating storm trajectories with network attributes. This subsection delves into the modeling approach for these failures using fragility curves to quantify the likelihood of failure based on storm intensity. a) Power Line Failures: Power line failures are modeled by assessing the impact of storm wind speeds on line segments using fragility curves. Each power line is evaluated as a single segment for simplicity, though the model supports finer granularity. At each time step, Algorithm 2 checks whether a line segment lies within the storm’s radius, represented as a circular buffer around the storm center. For affected segments, the local wind speed at the midpoint, denoted as wseg t, is obtained by adjusting the storm gust speed with a distancebased decay factor. This value is input to a lognormal fragility curve (Fig. 5a) [25], given in (4), parameterized by mean (µ) and standard deviation (σ). A random draw against this probability determines line failure (status = 0), with failed lines remaining offline until restoration. Phw Lwseg t=Zwseg t 0 1 x σ√2πexp"−ln x−µ2 2σ2#dx, (4) b) Road Blockage: Road blockage is modeled as treefall events obstructing road segments, with fragility curves relating wind speed to tree failure probability (Fig. 5b) [26]. Leveraging tree data from Section II-A3, Algorithm 3 checks if a road segment lies within the storm radius at each time step. For exposed segments, the lognormal form in (4), parameterized with tree-specific values (µT, σT)[26], gives the tree-fall probability pt. The blockage probability is then calculated as 1−(1 −pt·b)n, where nis the segment’s tree count and 6 Fig. 5. Wind Fragility Curve(s) for (a) Power Lines [25] (b) Trees [26]. Algorithm 2 Power Line Failure Assessment Require: Storm data Sdf , power lines Glines, fragility parameters (µ, σ), time horizon T Ensure: Line status Lstatus, local wind speeds Lwind 1: Initialize Lstatus ←ones(T, |Glines|) 2: Initialize Lwind ←zeros(T, |Glines|,segments) 3: for each time step tin Sdf do 4: Compute storm center pstorm and radius rt 5: for each line iin Glines do 6: Get segment midpoint m 7: if pstorm.distance(m)< rtthen 8: Compute decay factor dand local wind wseg tat m 9: Compute failure probability Phw Lwseg tusing (4) 10: if random draw < Phw Lwseg tthen 11: Set Lstatus[t:, i]←0 12: end if 13: end if 14: end for 15: end for 16: return Lstatus,Lwind b∈[0.5,1.0] is a blockage factor accounting for uncertainties such as lean, branching, and road proximity. A road segment is marked blocked (status = 0) if a random draw falls below this probability and remains impassable until restored. Roads with zero trees are immune. In the current implementation, only tree-fall events are considered, though other obstruction mechanisms (e.g., landslides, debris, structural collapse) can be added in future extensions. 3) Telecommunication Failure Propagation: The communication failure propagation model captures how windstorminduced power disruptions cascade to affect communication services. It uses power line status outputs from Algorithm 2 to assess downstream impacts on towers. Electrically isolated substations from all sources (Pgen) are identified via graphbased connectivity analysis and marked non-operational. Dependent towers are mapped through a substation–tower dependency matrix: when a supporting substation is disconnected, towers switch to backup power and, if the outage exceeds the threshold h, are marked as non-operational (status = 0). D. Dynamic Restoration Modeling The dynamic restoration model coordinates the recovery of power lines and roads after windstorm damage, adaptively deploying repair crews and MERs based on evolving failures. Road crews clear tree-induced blockages, while line crews Algorithm 3 Road Blockage Assessment Require: Storm data Sdf , roads Groads with tree counts, fragility parameters (µT, σT), time horizon T Ensure: Road status Rstatus, blockage probabilities Rprob 1: Initialize Rstatus ←ones(T, |Groads|) 2: Initialize Rprob ←zeros(T, |Groads|) 3: for each time step tin Sdf do 4: Compute storm center pstorm and radius rt 5: for each road iin Groads do 6: if pstorm.distance(Groads[i]) < rtthen 7: ntrees ←Groads[i].tree count 8: if ntrees >0then 9: Compute ptusing (4) with (µT, σT) 10: Draw blockage factor b∈[0.5,1.0] 11: Rprob[t, i]←1−(1 −pt·b)ntrees 12: if random draw < Rprob[t, i]then 13: Rstatus[t:, i]←0 14: end if 15: end if 16: end if 17: end for 18: end for 19: return Rstatus, Rprob Algorithm 4 Failure Propagation to Telecommunication Network Require: Power line status matrix Lstatus, substation-tower dependency map, backup duration h, time horizon T Ensure: Tower status matrix Tstatus 1: Initialize Tstatus ←ones(T, |Towers|){All towers initially operational} 2: Initialize outage timers O←zeros(|Towers|) 3: for each time step tfrom 1to Tdo 4: Identify disconnected substations St Nusing graph-based traversal on Lstatus[t] 5: Determine affected towers in Tt Nbased on substationtower mapping and St N 6: for each tower iin Tt Ndo 7: O[i]←O[i]+1{Increment outage timer} 8: if O[i]> h then 9: Tstatus[t:, i]←0{Mark tower as failed due to exhausted backup} 10: end if 11: end for 12: end for 13: return Tstatus restore failed lines once access routes are cleared. MERs provide temporary supply to substations isolated by line outages, thereby sustaining critical operations until structural repairs are completed. To ensure safety, road crews are dispatched only after local storm conditions have subsided, and line crews are mobilized only once unobstructed routes are identified. The process advances in discrete time steps (e.g., one hour), synchronized with hazard progression, and continues beyond the storm horizon until all damaged assets are restored. 7 Algorithm 5 Dynamic Restoration Strategy Require: Line status Lstatus, road status Rstatus, storm data Sstorm, topology Ntopo, Crew parameters (Src, Slc, vrc, vlc, Nrc, Nlc, Mcrew, tbase), MER parameters (Smr, vmr, dmr, Nmr, Mmr, Ncrit), time-horizon T Ensure: Updated Lstatus,Rstatus, Crew and MER assignments 1: Initialize road crews (Src), line crews (Slc), and MERs (Smr) in idle state 2: Build road graph Groads and power graph Glines from Ntopo 3: Set t←1 4: while there are blocked roads, failed lines, or disconnected substations do 5: Update Groads using Rstatus[t]and Sstorm[t] 6: Identify failed lines and blocked roads from Lstatus[t] and Rstatus[t] 7: Maintain task queues for damaged roads and lines 8: for each idle road crew (≤Nrc)do 9: Select reachable failed road segment (favoring line access); dispatch if storm has passed 10: Compute travel and repair time; update Rstatus[t:, i]←1after repair 11: end for 12: for each idle line crew (≤Nlc)do 13: Select reachable failed line; dispatch only if access path is unblocked 14: Compute travel and repair time; update Lstatus[t:, j]←1after repair 15: end for 16: for each idle MER (≤Nmr)do 17: Identify disconnected substations from Lstatus[t] 18: Prioritize substations in Ncrit; select one with minimal travel time 19: if reachable then 20: Dispatch to selected substation; support for dmr or until lines are restored 21: end if 22: end for 23: Update Crew and MER task states, locations, and costs 24: t←t+ 1 25: end while 26: return Lstatus, Rstatus, Crew and MER assignments Consistent with the framework definition in Section II, the term dynamic denotes a feedback-driven process in which restoration decisions are continuously updated in response to evolving hazard conditions and changing network states. As outlined in Algorithm 5, road and line crews are dispatched from predefined depots (Src, Slc) to repair blocked roads and failed power lines identified from status matrices. The numbers of road and line crews (Nrc, Nlc) are fixed, and task queues are updated at each time step. Road crews are dispatched only after local storm clearance, with reachability assessed from the road network graph Groads in which blocked edges are removed. Dispatch decisions follow explicit prioritization rules: (i) road crews first clear blockages that unlock access to failed lines, with other blockages ranked by connectivity benefit (e.g., number of substations reconnected); (ii) line crews are assigned only to safely reachable failures, with preference for those whose repair restores greater load or more substations; and (iii) MERs are deployed first to critical substations Ncrit, then to others ranked by impact and travel time, subject to Mmr (per site) and Nmr (fleet) limits. T=D v+tbase ·ϕ(asset attributes)(5) Once a task is selected based on these priorities, its restoration time Tis determined using (5), where Dis the shortest safe-path distance in Groads,vis the resource speed (vrc, vlc, vmr), and ϕ(·)scales repair duration according to asset-specific conditions (e.g., tree density, line length, or failure severity). A maximum of Mcrew crews per asset is enforced to avoid redundancy. If no safe path exists, resources remain “Idle” rather than attempting hazardous travel. MER provides temporary support for dmr or until reconnection, and the framework advances iteratively by updating failures, redeploying resources, and recording activity, thereby capturing dynamic restoration under evolving hazard conditions. E. Key Resilience Metrics This subsection elaborates the procedure applied to quantify multi-sectoral impacts of windstorm-induced disruptions on interdependent power, telecommunication, and road networks by computing resilience metrics continuously over the entire simulation horizon, encompassing both failure progression and restoration phases. These metrics capture CIS’s operational status and societal consequences, dynamically reflecting the cascading effects of failures and recovery processes over time. At each time step t, the operational performance of the power system is evaluated by computing the DNS using a DCOptimal Power Flow (DC-OPF) model, based on the power line status matrix Lstatus and network topology Ntopo. The total number of failed power lines is computed using (6a). Similarly, for the road network, the number of blocked segments at each time step is determined by (6b). In the case of the telecommunication system, towers are considered to be failed if their power outages persist beyond the predefined backup duration h, determined by (6c). Telecommunication tower status here is inherited directly from the power-driven failure propagation model (Section II-C3), ensuring consistency between substation disconnection, backup supply, and service availability. To quantify the spatial extent of telecommunication service loss, the effective per-transceiver output power is first computed using (7). Subsequently, the total transmit power per tower is computed as: TP AW=ntrx ·PAW. Nfail,lines(t) = X i∈Glines (1 −Lstatus[t, i]) (6a) Nblocked,roads(t) = X i∈Groads (1 −Rstatus[t, i]) (6b) Nfail,towers(t) = X i∈Tt N (1 −Tstatus[t, i]) (6c) PAW=Pmax ηPA ·(1 −σfeed)(7) 8 The total power is then converted to decibels relative to dBm and input into a frequency-specific propagation model (e.g., Extended Hata for sub-3 GHz or Rapport for higher bands) to determine each tower’s coverage radius, defining a circular service-loss buffer. This step models service degradation beyond binary tower outages, consistent with interdependencybased cascading analyses in [30]. Buffers from failed towers are merged into a unified coverage polygon Ct, which is overlayed with population region maps. For each population region rintersecting Ct, the overlap area Aoverlap is computed, and the coverage ratio is given by: ρr=Aoverlap/Aregion. The affected population within each region is then estimated as: paffected r=ρr×pr; where prrepresents the total population within region r. The total affected population at each time step tis calculated by summing paffected racross all affected regions. Finally, recovery progress is tracked by monitoring restored lines and cleared roads relative to the disrupted state, providing a comprehensive view of impacts across interconnected CIS. III. CASE STUDY AND SIMULATION RESULTS This section presents simulation results analyzing the impact of windstorm scenarios on interconnected CIS. A detailed illustrative scenario first demonstrates the framework’s capabilities in modeling cascading failures and safety-aware dynamic restoration across power, telecommunication, and transportation networks. The initial analysis focuses on disruptions caused by storm progression, followed by an evaluation of coordinated restoration, including crew deployment and MER support, and their influence on key resilience metrics. To test the framework’s robustness, sensitivity and trade-off analyses are performed: the former varies wind speed ranges and power line fragility parameters, while the latter compares DRS with and without the safety-aware dispatch rule. In both cases, the same set of 100 storm scenarios is used, ensuring that differences in outcomes reflect only parameter or strategy variations. Resilience metrics are averaged across scenarios to highlight consistent trends. This broader analysis underscores the framework’s applicability to both real-time planning and long-term resilience assessment. In terms of computational performance, all simulations are implemented in Python and executed on a MacBook Pro (Apple M3 chip, 18 GB unified memory). Indicative runtimes are: storm-path and failure state generation complete within seconds, while the DRS runs within minutes. Depending on scenario complexity (e.g., number of failed CIS assets), DRS execution ranges from under 5 to about 50 minutes, showing that the framework is tractable for ensemble-based analyses such as the 100-scenario sensitivity and trade-off studies. A. Test System and Simulation Data The case study utilizes the real 62-bus power network of Cyprus’s national grid, as documented in [31]. This system is spatially integrated with georeferenced data on distribution poles, telecommunication towers, and road segments from the national open-data platform [32], and complemented with spatial layers from the Cyprus Digital Twin (CyDT3) platform 3Cyprus Digital Twin (CyDT) overview: https://youtu.be/86WFkxe5fg4 [33]. Demographic data is sourced from Eurostat’s GEOSTAT 2025 dataset [34]. Electrical substations are linked to nearby telecommunication towers—via distribution poles—within a 2 km threshold, assuming 3-hour backup power per tower. The model simulates cascading failures across interdependent power, telecommunication, and transportation systems over a 24-hour horizon. A representative windstorm scenario is generated using parameters in Table II, and dynamic restoration is executed using the configuration in Table III. All parameter values in Table III are baseline assumptions used for illustration. A subset of substations near hospitals is designated as critical and prioritized during MER deployment, as described in Algorithm 5. The modeled CIS layers, storm trajectory for the illustrative scenario, and initial locations of repair crews and MERs are presented in Fig. 6. For sensitivity analysis, the wind gust speed range is tested at 20–50 m/s, 25–55 m/s, and 30–60 m/s to examine trends under different storm intensities. In addition, the mean of the line fragility curve is varied by ±20% from the value in Table II, reflecting line hardening or aging. Each configuration is evaluated over 100 scenarios, with resilience metrics averaged to assess robustness. B. Simulation Results from Illustrative Storm Scenario A detailed illustrative scenario shown in Fig. 6 is employed to demonstrate the proposed framework’s capabilities in modeling cascading failures and dynamic restoration across CIS. Results are analyzed under two cases and compared. The first case (base case, BC) quantifies system-wide disruptions without restoration, including failures in power lines, substations, road segments, and telecommunication towers, as well as the resulting communication outages affecting the population. The second case evaluates the effect of the safety-aware dynamic restoration strategy (DRS), incorporating coordinated crew dispatch and MER support. Resilience metrics are expressed as percentages and population affected (in millions), with emphasis on peak values and their duration. 1) Base Case: Hazard Impact on CIS: The impact of the windstorm shown in Fig. 6 on interconnected CIS is quantified in this case. The storm begins affecting the systems at simulation time t= 6 (as shown in Fig. 6), with its 24-hour progression analyzed through key resilience metrics presented in Fig. 7. Figure 7a illustrates the degradation of the electrical network: power line outages (black solid line) begin at t= 6 with 3% degradation and exceeding 10% by hour 11. Substation losses (black dashed line) follow a similar trajectory, surpassing 10% by hour 7. These disruptions result in power outages, measured in terms of DNS in megawatts (MW), with an initial DNS value of 136 MW observed at TABLE II PARAMETERS FOR HAZARD ENGINE Parameter Value Range of Wind Gust Speed 30–60 m/s Line Fragility (Mean, Std) 3.81, 0.18 Tree Fragility (Mean, Std) 3.91, 0.33 Road Buffer Distance 30 m Time Horizon 24 hours 9 Fig. 6. Spatial Distribution of CIS, a Representative Windstorm Path, and Initial Locations of Restoration Crews and MERs across the Study Region. TABLE III PARAMETERS CONSIDERED FOR DYNAMIC RESTORATION Parameter Value Road crew start points Src 3 Line crew start points Slc 3 MER start points Smr 1 Number of road crews Nrc 4 Number of line crews Nlc 4 Number of MERs Nmr 4 Max MERs per substation Mmr 1 MER Support Duration (dmr ) 8 hours MER Capacity 0.5 MW Max crews per asset Mcrew 1 Crew travel speed vcrew 40 km/h MER travel speed vmr 30 km/h Base repair time tbase 2 hours MER support duration dmr 8 hours Critical substations Ncrit {10, 28, 38, 40, 44, 55, 60, 61} Source (generator) buses Pgen {5, 23, 25} t= 6 (red dotted line). The resulting substation outages induce cascading failures, leading to telecommunication tower losses once their 3-hour backup power is depleted. As shown in Fig. 7b, telecom tower outages exceed 15% by hour 14. Figure 7c depicts the affected population, which peaks at approximately 0.5 million individuals by hour 14, reflecting sustained communication disruptions in densely populated areas. Figure 7d shows the percentage of affected roads due to blockages caused by fallen trees, which progressively increases to nearly 1% (i.e., 14 road segments) by hour 12. Overall, this case study demonstrates the framework’s ability to simulate hazard-driven cascading failures in interdependent CIS. 2) Effect of Dynamic Restoration Strategy: This section presents the outcome from the safety-aware DRS that coordinates the deployment of repair crews and MERs to restore CIS impacted by the storm shown in Fig. 6. The resulting resilience metrics are illustrated in Fig. 8, which spans 40hour horizon from the onset of simulation (t= 0) through storm progression, impact, and complete restoration (t= 39). Fig. 7. Hazard Impact on CIs without Dynamic Restoration: (a) Impact on Power Network, (b) % of Lost Telecommunication Towers, (c) Affected Population due to Communication Disruptions, and (d) % of Affected Roads. Fig. 8. Hazard Impact on CIs with Dynamic Restoration: (a) Impact on Power Network, (b) % of Lost Telecommunication Towers, (c) Affected Population due to Communication Disruptions, and (d) % of Affected Roads. The corresponding repair crew and MER activity timeline is shown in Fig. 9, starting from the first impact on CIS (t= 6) and extending until complete recovery (t= 39). This timeline reflects both the storm period, when resources remain inactive due to safety constraints, and the subsequent