This paper presents a method that allows urban planners and municipal engineers to identify critical components of interdependent infrastructure systems. The in- tentofthemethodistoprovideameansofmodelingtheimpactofcapacity-related changes (e.g., population growth, component degradation) on a city’s ability to deliver resources to critical locations. Infrastructure systems are modeled as ﬂow networks in which capacities, demands, and supply constraints vary over time; demand nodes also have criticality ratings that allow a user to model levels of importance. Interconnections between infrastructure systems are represented by physicalandgeospatialdependenciesatacomponentlevel. Aﬂow-basedcentral- ity measure is used to rank components according to their role in the delivery of resources to critical locations. A simple instantiation of the method is presented andevaluatedonadistrict-scalemodelofacitythatcontainsinterconnectedwater and electricity networks. Finally, two forms of reliability analysis are demon- strated: a composite measure incorporating edge reliability, and a variation on standard component failure/degradation analysis. Keywords: Component importance measures, Centrality measures, Complex systems, Network science, Infrastructure reliability 1. Introduction 1 This paper presents a method for identifying critical components in interde- 2 pendent,urbaninfrastructuresystems,whereacomponentisdeemedcriticaltothe 3 extentthatitisreliedupontodeliverresourcestocriticallocations(e.g,. hospitals). 4 The goal of the work is to provide urban planners and municipal engineers with 5 a means of visualizing the impact of capacity-related changes. This goal is met 6 by combining a ﬂow-based centrality measure with a network-based approach to 7 Preprint submitted to Int. Journal of Critical Infrastructure Protection April 3, 2021 infrastructuresystemsthatpermitsmodelersto: (1)deﬁnephysicalandgeospatial 8 interdependencies;(2)representtemporaldynamicsthroughtheuseoftimeseries 9 for key system variables, and; (3) label locations with criticality values in order 10 to represent asset importance. While useful in its own right, this basic method 11 may be combined with standard reliability methods in order to model the eﬀects 12 of component failure. 13 The perspective in this work is resource-based. Residents of cities depend on 14 infrastructuresystemstodelivernotonlyphysicalresourcessuchaswaterandgas, 15 butalsoarangeofsocialgoodsrangingfromeducationtohealthcare. Disruptions 16 in the delivery of resources and/or services can have extremely deleterious con- 17 sequences, particularly for critical locations such as hospitals and transportation 18 hubs. Methods for identifying infrastructure components that are relied upon to 19 supply critical locations with resources can be used for a variety of purposes, 20 including maintenance scheduling, disaster recovery, critical infrastructure pro- 21 tection, and urban planning. 22 In keeping with the focus on resource delivery, the method presented in this 23 work represents infrastructure systems as ﬂow networks –graphs that allow com- 24 modities to ﬂow from a set of supply nodes to a set of demand nodes. Edges of 25 thenetwork(i.e.,pipes)arelabeledwithcapacityvaluesthatlimittheirmaximum 26 ﬂow. Demand nodes represent consumers of resources such as hospitals, laun- 27 dromats, and homes. The modeler supplies the various system variables, such as 28 capacities, supply constraints, and demands. 29 Unlike typical ﬂow networks, however, the demand nodes are labeled with 30 criticality ratings that indicate the relative importance of the structure at that 31 location. Thesecriticalityratings,whicharealsodeﬁnedbythemodeler,arethen 32 used to calculate a centrality measure–the critical ﬂow centrality (“CFC”). CFC 33 valuesarepropagatedfromthedemandnodesthroughtherestoftheinfrastructure 34 system’s components. The user can then visualize the system’s components (e.g., 35 withaheatmap)accordingtotheirroleindeliveringresourcestocriticallocations. 36 Since urban infrastructure systems are interconnected in various ways, users 37 may deﬁne physical and geospatial dependencies . Physical dependencies exist 38 whencomponentsofsystem requireresourcessuppliedbysystem . Geospatial 39 dependencies exist when components of andare co-located. The inclusion 40 of dependencies allows the user to model the potential for cascading failure, 41 simultaneous disruption due to weather events, and other common scenarios. In 42 particular, criticality ratings are propagated across interdependent systems so that 43 the user may immediately determine which components of System play an 44 important role in supplying the critical components of System with resources. 45 2 Apartfromdependencymodelingandcriticalityratings,thethirddistinguish- 46 ingcharacteristicofthisworkistheuseof timeseries . Keysystemvariables(e.g., 47 capacityconstraints,supplyconstraints,demandlevels,andcriticalityratings)are 48 dynamic rather than static, allowing the user to model phenomena at a variety of 49 time scales. For instance, a geographer may be interested in studying the impact 50 ofpopulationgrowthonaparticulardistrictoverdecades,whileanurbanplanner 51 might use the same model with diﬀerent data in order to visualize the short-term 52 consequences of planning interventions. 53 Finally,theCFCmeasurebecomesevenmoreinterestingwhencombinedwith 54 standard reliability methods. The CFC was designed to identify components that 55 play a substantial role in delivering resources to critical locations. However, the 56 failure of a component with a high CFC value will not necessarily result in a loss 57 of service, as there may be alternative (fallback) routes available. Combining the 58 CFC with reliability methods allows a modeler to reason explicitly about fallback 59 routes and the consequences of component failure. 60 Thestructureofthispaperisasfollows. Section2providesusefulbackground 61 information, while Section 3 introduces the methodology used in this paper, in- 62 cludingtheCFCmeasure. Section4providesanevaluationofthemethodologyon 63 adistrict-scalemodelofacity. Section5discussestwoformsofreliabilityanalysis 64 that can be combined with the CFC measure to model potential service loss due 65 to component failure. The paper closes with suggestions for future research. 66 2. Background 67 The method in this paper can be viewed as a combination of techniques from 68 network science and critical infrastructure protection. The fundamental build- 69 ing block is a component importance measure (“CIM”) (e.g,. [1], [2], [3]) that 70 estimates the degree to which a given component participates in the delivery of 71 resources to critical locations. Before discussing the method in detail, a quick 72 discussion of relevant background material is required. 73 2.1. Network Science and Centrality Measures 74 Networks are a common choice of modeling mechanism in many ﬁelds, and 75 critical infrastructure protection is no exception (see [4]). For example, many ap- 76 proachestoinfrastructurevulnerabilityandresiliencemakeuseoftechniquesfrom 77 networkscience. Fromtheperspectiveofthecurrentpaper,themostimportantof 78 these techniques are the centrality measures , which are used to identify the most 79 central components in a network (see [5, 6, 7]). 80 3 Numerous centrality measures exist [8], the most intuitive of which are: (1) 81 nearness measures , which determine a given component’s centrality by means 82 of its proximity to other components, and; (2) betweenness measures , which 83 deem components to be central to the extent to which they stand between other 84 components as intermediaries. While nearness and betweenness measures focus 85 on network topology, dynamical measures are based on processes (e.g., ﬂows, 86 diﬀusion) taking place on the network. 87 Theprogenitorofthemethodusedinthispaperis ﬂowcentrality [9]. Consider 88 asimplenetworkwithnodes +andlinks. AnodeEisconsideredtobe between 89 othernodesDandFtotheextentthatthemaximumﬂowbetween DandFdepends 90 onE. Nodes are deemed centralto the extent that they facilitate maximum ﬂow. 91 Stated formally, for DEF2+, let<DFbe the maximum ﬂow between Dand F,andlet<DF¹Eºbethemaximumﬂowbetween DandFthatdependson E. Then theﬂow centrality (“FC”) of a node E2+is the degree to which the maximum ﬂow between all unordered pairs of nodes depends on E: ¹Eº=Õ D<F<E<DF¹Eº (1) 2.2. Interdependent Infrastructures 92 Infrastructure systems are typically coupled to the extent that the failure of 93 components in one system can cause failures in connected systems [10]. These 94 interdependent systemsaretypicallymorefragilethansolitarysystems[11],with 95 additionalfailuremodes(e.g., cascadingfailures [12])thatcanbe quitecomplex. 96 For example, water distribution systems impose much greater cascading damage 97 onothersystemsthantheyreceiveinreturn[13],andtheyseemtodisplayagreater 98 propensity to initiate cascading failure in other systems [14]. 99 Infrastructuresystemscanbedisruptedinnumerousways,includingdeliberate 100 attacks, component failures, and natural disasters. Much of the existing research 101 oncriticalinfrastructureprotection,forinstance,hasfocusedonprotectinginfras- 102 tructures against damage due to extreme weather or deliberate attacks [4, 15, 16]. 103 Component failure has been studied extensively in the ﬁeld of reliability engi- 104 neering (e.g., [17]) and in the various engineering disciplines (e.g., water [18], 105 drainage[19],electricity[20],telecommunications[21],andtransportation[22]). 106 Disruptionofnetworkshasalsobeenconsideredinoperationsresearch(e.g.,[23]), 107 computerscience(e.g.,[24,25]),networkreliability(e.g.,[26,27]),graphtheory 108 (e.g., [28, 29]), and network science (e.g., [30]). 109 However, disruptions are not the only phenomena to have an impact on urban 110 infrastructure. Infrastructure systems are inﬂuenced by a variety of factors, in- 111 4 cluding: (1) populationgrowth ,whichcanpushthecapacityconstraintsoflegacy 112 infrastructure;(2) componentdegradation ,whichcanreducecomponentcapacity; 113 (3)maintenanceactivities ,whichcaninterdictﬂowthroughselectedcomponents, 114 and; (4)planning interventions (e.g., the development of new residential subdivi- 115 sions), which can have eﬀects both on system topology and on demand patterns. 116 In some cases, these factors can alter ﬂow distributions so that resources needed 117 by critical locations are dependent on unreliable infrastructure. 118 Various research communities have advocated an integrated view of infras- 119 tructure systems, and a growing body of work is available on interdependencies 120 (e.g., [31, 32]). For instance, homeland security initiatives following the Septem- 121 ber11thterroristattacksintheUnitedStatesspurrednumerouseﬀortsaddressing 122 infrastructureinterdependencies(e.g.,[33]). Overviewsoftechniquesforthemod- 123 eling and simulating interdependent critical infrastructure systems may be found 124 in several places, including [34]. 125 2.3. Modeling Interdependent Infrastructures with Networks 126 Oneapproachtoanalyzinginterdependentinfrastructuresystemsinvolvesmod- 127 eling them as interdependent networks [32]. Interdependent (or multilayer [35]) 128 networks have received increasing amounts of attention of late, particularly from 129 thephysicsandnetworksciencecommunities. Arecentsurveypapercanbefound 130 in[36],whilebooksonthetopicarereadilyavailable(e.g.,[37,38,32,39,40,35]). 131 Networkisdependent onnetworkifthestateof caninﬂuencethestate 132 of[41] (see also [42]). Dependencies can be classiﬁed as follows [14]:1 133 1.Physical dependencies , in which the state of is aﬀected by the material 134 outputs/ﬂows of ; 135 2.Geospatial dependencies , in which certain components of andare in 136 suchclosespatialproximitysuchthatlocaleventscanaﬀectbothnetworks; 137 3.Informationaldependencies ,inwhichandareconnectedbyinforma- 138 tion and communications technology; 139 4.Social dependencies , in whichaﬀectsalong social dimensions; 140 5.Proceduraldependencies ,whereaﬀectsonthebasisoforganizational 141 or regulatory structures, and; 142 6.Financialdependencies ,wheremarketconditions,ﬁnancialcrisesandother 143 economic events allow one network to aﬀect another. 144 1Alternative classiﬁcations appear in [43, 44]. 5 There are many ways to represent these dependencies in network models, a dis- 145 cussion of which is beyond the scope of the paper. 146 2.4. Finding Critical Components in Interdependent Networks 147 Numerous researchers have proposed methods for identifying critical compo- 148 nents in interdependent networks. Typical examples are described below. 149 Apostolakis and Lemon [45] evaluate the vulnerability of geospatially in- 150 terdependent infrastructure systems (gas, water, electric) by identifying critical 151 locations—geographical points that are susceptible to attack. Each system is rep- 152 resentedasadirectednetworkinwhichverticescanrepresentnotjustjunctionsbut 153 also physical features (e.g., manhole covers). Co-location of assets (e.g., shared 154 service tunnels) is modeled by allowing vertices from one graph to appear in an- 155 other. (Physical dependencies, such as the use of electricity by the water system, 156 are not modeled). In their approach, a set of attack scenarios is identiﬁed and 157 the networks are analyzed in order to identify minimal cut sets (see [20]). The 158 resulting vulnerabilities are prioritized by: (1) the degree to which the targets are 159 accessible to the attacker (i.e., susceptibility), and; (2) the value of the targets 160 from the standpoint of the decision-maker, calculated by summing their expected 161 disutilities. The susceptibility and value are combined to yield a vulnerability 162 category — one of ﬁve colors ranging from green to red. 163 Lee et al. [43] provide a method for prioritizing service restoration activities 164 in an interdependent system-of-systems. Each independent system is represented 165 as a ﬂow network that carries commodities, composed of edges and vertices that 166 may both have capacity constraints. Dependencies are modeled as additional 167 constraints in a mixed integer network ﬂow model. In addition to geospatial and 168 physicaldependencies,theyallow shareddependencies (i.e.,formulti-commodity 169 ﬂow networks) and exclusive-or dependencies (i.e., to allow ﬂow on a multi- 170 commodity network to be restricted to one type of commodity at a time). 171 Duenas-Osorio et al. [46] study the interdependency of electricity and water 172 systemsfromatopologicalstandpoint. Bothgeospatialandphysicaldependencies 173 are modeled, with the water system requiring electricity for pumps, lift stations, 174 and control units. Conditional probability distributions are used to model po- 175 tential failures of water system components given failure of electricity system 176 components. Three types of vertex removal strategies are used to model disrup- 177 tions; for each such disruption, a set of metrics are calculated: (1) nodal degree; 178 (2)characteristicpathlength[47];(3)clusteringcoeﬃcient[48],and;redundancy 179 ratio. Flows of water or electricity are not modeled. 180 6 Buldyrev et al. [10] examine the impact of electricity system disruptions on 181 the internet. Geospatial dependencies are modeled by assigning each internet 182 server to the closest power station. Disruptions are initiated by removing power 183 stationsandtrackingresultingnodalfailures—inparticular,anode Eisconsidered 184 to be failed if: (1) all of E’s neighboring nodes are failed, or; (2) the geospatially 185 coupled node in the electricity network is failed. Nodes are ranked according to 186 theconsequencesofremoval. Theauthorsarguethatdisruptionofasmallnumber 187 of nodes in the electricity system is suﬃcient to provide cascading failures in the 188 internet network. 189 Galvan and Agarwal [49] perform vulnerability analysis on interdependent 190 infrastructures by examining the impact of disruptions. Each infrastructure is 191 representedasaﬂownetworkwithauniqueresourcetype. Ineachiterationofthe 192 analysis, a single node is selected for failure (disruption). After recomputing the 193 ﬂow solution, the algorithm identiﬁes every node that is in violation of capacity 194 constraints. These latter nodes are then disabled and the process repeats itself 195 until no more failures occur. The authors introduce a new vulnerability metric 196 -1, deﬁned as the fraction of nodes that fail after the ﬁrst step of the cascading 197 failure process. After using -1to rank nodes, they compare the results against 198 traditional centrality measures (i.e., nodal degree, the ﬂow value for the non- 199 disrupted solution, and network eﬃciency). 200 Svendsen and Wolthusen examine interdependent critical infrastructures in a 201 series of papers [50, 51, 52, 53]. Their models represent multiple concurrent 202 types of dependencies, categorized at a high level into storableandnon-storable 203 types. Each vertex Ein a network can act as a producer or consumer of up to 204 <diﬀerent resources, and for each such resource Ehas a corresponding buﬀer. 205 Theauthorsinvestigatenumerousissues,includingthebehaviourofsystemswith 206 cyclic interdependencies. 207 Kotzanikolaouetal[54]provideamethodforidentifyingthreatstoinfrastruc- 208 ture systems that may have a signiﬁcant cumulative eﬀect. From a risk analysis 209 table detailing the risks to particular infrastructures, the authors construct a risk 210 dependencygraph (“RDG”)inwhich: (1)anoderepresentsaninfrastructuresys- 211 tem, and; (2) a directed edge -!.from system -to system.represents a 212 risk dependency (i.e., -poses a risk to .). Edges are labeled with risk values 213 thatrepresentthelikelihoodofdisruption,aswellassocietalimpact. Theauthors 214 provide a method for computing higher order risks from such a graph. 215 Stergiopoulos et al [55] use an RDG in combination with centrality measures 216 (e.g.,betweenness,eigenvectorcentrality,nodedegree)toidentifycriticalsystems. 217 Thecentralitymeasuresarecomputedinordertoidentifynodesthataﬀectcritical 218 7 risk paths in the RDG. A decision-tree algorithm is then used to select a subset 219 of these nodes for risk mitigation. Testing their approach on empirical data, the 220 authors make a number of observations about the relationship between centrality 221 and risk (e.g., that the most critical paths in the RDG tend to involve nodes with 222 high centrality). 223 ShahraeiniandKotzanikolaou[56]provideamethodtoaidinthedesignof wide 224 areameasurementsystems (“WAMS”),whicharecomposedofthemeasuringand 225 communicationslayersofasmartgrid. Theauthorsprovideamodelthatcaptures 226 theinternaldependenciesbetweentheselayersinadependencygraph. Givensuch 227 a graph, the importance matrix of the bus components is determined, and then 228 the importance metrics of all of the rest of the WAMS components are computed 229 using a centrality measure. An optimization algorithm is used to distribute the 230 importance of the WAMS elements in order to avoid single points of failure. 231 Interestingly, while the authors acknowledge that the importance of power grid 232 componentsisdynamic,changingovertime,theyuseanon-temporalapproachfor 233 the initial design of a WAMS. 234 The present work is diﬀerent from these examples of prior art. As discussed 235 below,theentiresetofinterdependentsystemsisrepresentedasagraph. Similarly 236 toRDGs,eachnodeisanindividualinfrastructuresystem. However,unlikeRDGs, 237 the edges in the graph do not represent risk dependencies; rather, they represent 238 actual component-level interdependencies. Furthermore, the graph formed from 239 theinterdependentsystemsisusedonlyforcontrollingtheorderofexecution,and 240 not for analysis. 241 Of all the aforementioned works, the method in this paper is closest to Apos- 242 tolakis and Lemon [45] in overall intent. In both cases the computational model 243 allows for co-location. However, the present work also allows for physical depen- 244 dencies,aswellasdynamicbehaviorthroughtheuseoftimeseriesforkeysystem 245 variables. Thetemporalaspectsofthecurrentworkalsodistinguishitfrommany 246 of its predecessors. 247 8 3. Methodology 248 Thegoalofthisworkistoexploremeansbywhichurbanplanners,municipal 249 engineers and other decision makers can identify critical components of interde- 250 pendent infrastructure networks. When embodied in software, such methods can 251 be used to support decision makers engaged in maintenance scheduling, zoning, 252 capacity planning, or other activities related to municipal infrastructure. 253 3.1. Overview 254 Thepaperprovidesanexampleofsuchamethod,basedonacentralitymeasure 255 thatcombinesclassicalﬂowcentrality[9]withconceptsfromcriticalinfrastructure 256 systems (e.g., [45]). The perspective in the paper is resource-based , focusing on 257 theroutesbywhichresourcesaredeliveredtoconsumers. Componentsaredeemed 258 critical to the extent that they are involved in facilitating the ﬂow of resources to 259 critical locations. 260 Computation of the centrality measure, critical ﬂow centrality (“CFC”), can 261 beaccomplishedinseveralways(see[57]). Inthecurrentpaper,adiscrete-valued 262 approach is taken in which: (1) an infrastructure system is represented as a ﬂow- 263 network; (2) demands, capacities, and supply limits are given as integers, and; 264 (3) each demand node in the network is assigned a real-valued criticality rating. 265 Network ﬂows are simulated with a standard maximum ﬂow algorithm; once a 266 ﬂowhasbeendeﬁned,asearchalgorithmcomputesexpectedcontributionofeach 267 component to the critical ﬂow within the network. 268 Since infrastructure networks are not independent of each other, physical and 269 geospatial dependencies may be introduced between individual infrastructures. 270 The most important of these for the present paper are physical dependencies in 271 which resources provided by one system (e.g., electricity) are used by another 272 system(e.g.,waterpumps). Oneofthemaincontributionsofthepaperistoshow 273 how CFC values can be propagated from one infrastructure system to another. 274 The method is demonstrated by applying it to a district-level model of a city. 275 Each structure has a type, a criticality rating, and a set of demand curves (time 276 series) for resources. For reasons of brevity, only two infrastructure systems 277 (electricity, water) are shown. The simple method provided in this paper also 278 assumes that the physical dependencies between individual infrastructures are 279 acyclic. 280 The main thrust of the demonstration is to show that: (1) the computation 281 of CFC values can be performed eﬃciently, enabling their use in interactive GIS 282 applications; (2) CFC values can correctly propagate between system models, 283 9 and;(3)CFCcomputationscanbeintegratedwithstandardreliabilitymeasuresto 284 provide a composite view of a system. 285 The CFC measure itself is completely general, requiring only a ﬂow solution 286 and a network topology. The method presented in this paper uses the discrete 287 algorithms to compute values for each individual infrastructure system—namely, 288 (1) an integer-valued, maximum-ﬂow algorithm to approximate resource ﬂow 289 withininfrastructures,and;(2)amodiﬁedgraph-searchalgorithmtocomputeCFC 290 values. These design choices are for ease of explanation, and more sophisticated, 291 heterogeneous systems can be accommodated. One can model a water system 292 using hydraulic techniques [58], for example, coupling it to an electricity system 293 that is simulated using its own domain-speciﬁc methods. Given a ﬂow solution 294 andnetworktopology,CFCvaluescanbecomputedbyusingMarkov-chainMonte 295 Carlo or random walks (see [57] for details). 296 3.1.1. Integration with GIS 297 This work was motivated by the problem of providing adequate decision sup- 298 portforurbanplanning. Forinstance,densiﬁcationofurbanareasisaccompanied 299 bygreaterdemandforresources;theincreaseddemandcould: (1)violatecapacity 300 constraints, as in the case of the London sewer systems [59, 60], or; (2) threaten 301 the ability of a legacy infrastructure system to reliably deliver services to critical 302 locations such as hospitals and transportation hubs. Urban planners could beneﬁt 303 from tools that allow them to visualize the impacts of land-use decisions on the 304 provision of critical resources and/or services. 305 Eﬀective modeling of integrated infrastructure systems requires more than a 306 static,single-perspectiveapproach. Managementofdisruption(andpreventionof 307 cascading failures) requires an understanding of system dynamics [61]. Further- 308 more, any model used to study the disruption of interdependent infrastructures 309 needs to support two diﬀerent perspectives [43]: (1) a ‘system-of-systems’ view 310 thatfocusesondependencies,and;(2)atraditionalviewofeachindividualsystem 311 that is familiar to managers/specialists. 312 One means of providing infrastructure models that support multiple perspec- 313 tives is through the use of geographical information systems (“GIS”) software. 314 In fact, the critical information protection community has begun to use GIS as a 315 platformforresilienceandvulnerabilityanalysis[62]. Forthisreason,themethod 316 describedinthispaperwasexplicitlydesignedforintegrationwithinGISsoftware. 317 10 3.1.2. Data Sources 318 Two major challenges arise when data sources are considered. First, data on 319 infrastructure systems does not always exist, and particularly not in a form that 320 permits detailed analysis of interdependencies. Second, infrastructure systems 321 in many countries (e.g., the United States power grid) are not under the control 322 of a single entity [61], making the data collection process diﬃcult. The lack of 323 information on infrastructure assets has motivated some researchers to develop 324 techniques for inferring asset locations from proxy data sources (e.g., [63, 64]). 325 The model used in this paper is a mixture of synthetic and empirical com- 326 ponents. The basic topology (i.e., road and parcel structure) was taken from 327 downtown Toronto, albeit the boundaries were simpliﬁed in order to make dia- 328 gramsfeasibleandtoconveythebasicmethodclearly. Resourcedemanddatawas 329 obtained from published studies (e.g., [65]) and from municipal utilities. 330 3.1.3. Implementation 331 ThesamplemethodwasimplementeddirectlyinC++andOpenGL.Roadand 332 building information was obtained from OpenStreetMaps, imported into ESRI 333 CityEngine,andeditedmanuallytoremoveartifacts. Custompythonscriptswere 334 usedtoexporttheroadnetworktopology,block/lotgeometry,andbuildingshapes 335 from CityEngine to Extensible Markup Language (“XML”) ﬁles. Infrastructure 336 systemswerecreatedmanuallyusingtheapplication’seditingfunctionality. Lastly, 337 the diagrams shown in this paper were generated by exporting model geometry 338 directly to Scalable Vector Graphics (“SVG”) format. 339 3.2. Modeling Approach 340 This section discusses the building blocks of the simpliﬁed model, including: 341 (1)thenetworkrepresentation;(2)theuseoftimeseriesforkeysystemvariables; 342 (3) criticality ratings, and; (4) interdependencies. 343 3.2.1. Network Representation 344 An infrastructure system is modeled as a weighted, capacitated, ﬂow network 345 =h+iwhere is a set of nodes, ++is a set of edges: 346 •each nodeE2 +has Euclidean coordinate®F¹Eº=¹EGEHEIº2R3, as 347 well as an (optional) capacity constraint 2¹Eº2N. 348 •each edge4=¹E8E9º2 has acapacity2¹4º2N, aﬂow5¹4º2N, 349 and alength;¹4º2Rdeﬁned ask®F¹E8º®F¹E9ºk2. 350 11 Each network is amulti-graph in which multiple edges may connect a given 351 pair of nodes, allowing for redundant (fallback) connections. Bi-directional rela- 352 tionships, cycles, and self-loops are all permitted. 353 Asaﬂownetwork, containsbothsource(supply)andsink(demand)nodes. 354 The set of source nodes is+(=fB1B2B?g+, and the set of demand 355 nodesis+=f31323:g +with+(\+=. All other nodes are 356 calledtransmissionnodes . Multi-functionalnodesaresupportedusingastandard 357 maximum ﬂow reduction (as described in Section 3.4.1). 358 Aﬂowon isareal-valuedfunction 5:!Ron ’sedgesthatobeysthree 359 ﬂow properties: 360 1.Capacity Constraints : for all4=¹E8E9º2, we have5¹4º2¹4º. 361 2.SkewSymmetry : forall4=¹E8E9º2,wehave5¹¹E8E9ºº=5¹¹E9E8ºº. 362 3.FlowConservation : foralltransmissionnodes EC2+¹+[+(º,wehave 363Í E2+5¹¹ECEºº=0. 364 Agivennetwork supportsasingletypeofresource/commodity,unlikethemulti- 365 commodity approach in [52]. The value of a ﬂow is deﬁned as the ﬂow exiting 366 the source nodes:j5j=Í E2+Í B2(5¹BEº. 367 3.2.2. Supply Constraints and Demand Distributions 368 Supply constraints and resource demands are represented as discrete, integer- 369 valuedtime series (see [66]). (While capacities can also be represented as time 370 series,thedemonstrationassumesnodeandedgecapacitiesarestatic.) Forsimplic- 371 ity, each time series is assumed to be regularly sampled at times C82)=»01¼. 372 They can be interpreted as the output of functions: 373 •Each supply node E2+(may be assigned an optional supply constraint 374 function5B E¹Cº:)!N¸that gives the maximum amount of resource that 375 may be supplied from Eat timeC. 376 •Eachdemandnode 32+hasamandatory demandfunction X3¹Cº:)! 377 N¸that gives the amount of ﬂow required by node 3at timeC. 378 Anassignment toanetworkinvolvesspecifyingfunctions(timeseries)forallrele- 379 vantnodes. Computationsonthenetwork(e.g.,networkﬂowsolutions,criticality 380 measures) are performed for each time C82). Values from previous time steps 381 C:may be used as input for computing values in the current time step C8(where 382 C:C8). This permits the method to represent delaysin resource utilization. 383 12 3.2.3. Criticality Ratings 384 Acriticalityfunction 2A:+!Rmapsdemandnodes 32+toacriticality 385 rating2A¹3º. Although it is possible to use binary(e.g., critical, non-critical) or 386 categorical (e.g., low, medium, high) representations, this paper focuses on the 387 continuous variant in which criticality ratings take on values between 0and1. 388 3.2.4. Interdependencies 389 Asystem-of-systems (“SoS”)modelconsistsofasetof :infrastructuresystems 390 S=f(1(2(:g. As shown in Figure 1, two types of dependencies are 391 permitted between pairs of elements from S: 392 1.geospatial dependencies , which arise when elements from network are 393 co-located with those from network . 394 2.physical dependencies , wherein elements in network require resources 395 ﬂowing through network . 396 d3 n6 n5 n2A: S1 n3 S2 n1 d1 n4 n7 d2 n8 d7 n4 n7 n2B: S1 n6 n1 d3 n5 n9 d6 n8 d1 n3 d2 d5 d4colocation colocation colocation colocation physical physical colocation colocation Fig. 1.Infrastructure systems ( ) and (), linked by geospatial and physical dependencies. Dependenciesarerepresentedas interlinks betweennetworks[35]. Incontrast 397 toApostolakisandLemon[45],nodesfromonenetwork (8donotappeardirectly 398 in another network (9. This design choice makes it easier to integrate disparate 399 modeling methods for each individual infrastructure system (see [67]). 400 Interlinksrepresentingphysicaldependenciesareimplementedwiththeuseof 401 interconnectionrecords . Referringto Figure2,let(1representawaterdistribution 402 13 system, and let (2represent an electricity system. A dependency between water 403 nodeE12(1and electricity node E22(2is represented by an interconnection 404 record '¹E1E2º. The amount of resource 'demanded of (2byE1(e.g., the 405 amountofelectricityrequiredtooperateagivenwaterpump)isgivenbyafunction 406 5':(1+!R. For instance, a pump at E1might demand a constant amount of 407 electricity per unit time, or it may require power proportional to the ﬂow 5¹E1º 408 throughE1(e.g.,5'¹Eº=25¹Eº).Delayscanbeaccommodatedbydeferring this 409 demand to later time steps. 410 Block/Lot Water Electricity Source Source Pump Pump Fig. 2.Infrastructure model with two layers, showing resource ﬂow between water pumps and electricity nodes. Dependencies between network elements imply dependencies between sys- 411 tems. If an interconnection record exists that maps elements of (1to elements 412 of(2, we say that (1isphysically dependent on(2, represented as (1!(2. 413 Mutual dependency between systems makes the computational task more diﬃ- 414 cult. The methods of Svendsen and Wolthusen (e.g., [52]) accommodate mutual 415 dependencies using multi-commodity ﬂows, but this approach does not allow for 416 infrastructure-speciﬁc network representations and solution methods. 417 In this paper, the set of physical (resource) dependencies between systems in 418 Sistakentoforma directed,acyclicgraph (“DAG”)Gthatcanbeorderedwitha 419 topological sort (see [68]). In contrast, geospatial dependencies are not restricted 420 in such a fashion. 421 14 3.3. Critical Flow Centrality 422 TheCriticalFlowCentrality (“CFC”)ofacomponentisameasureofitsrole insupplyingresourcestocriticallocations[57]. Recallthatthe ﬂowinnetwork (given assignment ) is the aggregate of all ﬂows reaching the demand nodes: ¹ º=Õ 325¹3º (2) Thecritical ﬂow in network (given assignment ) is the set of ﬂows reaching the demand nodes, weighted by criticality:  ¹ º=Õ 325¹3º2A¹3º (3) A component 2(i.e., node or edge) is deemed to be important to the extent that it carries critical ﬂow. Let 5¹23ºbe the ﬂow that reaches 32from2given assignment, and let»5¹23º¼be its expectation. Then the critical ﬂow centrality (“CFC") of component 2under assignment is: ¹2º=Õ 32+2A¹3º»5¹23º¼ This quantity may be normalized by the critical ﬂow  ¹ º: 0¹2º=¹2º  ¹ º=Í 322A¹3º»5¹23º¼Í 322A¹3º5¹3º(4) Computing the CFC thus reduces to computing the probability ?¹3j2ºthat 423 a unit of commodity passing through component 2ends up in demand node 3. 424 Whiletherearenumerouswaystoaccomplishthistask(e.g.,Markovchains),this 425 paper uses a discrete, search-based approach. 426 For each time step C, a ﬂow solution ¹Cºis generated and represented in a 427 secondary graph 0that gives an adjacency-list representation of the stochastic 428 transition matrix induced by ¹Cº. All nodes of are present in 0, but edges of 429 are only present if they have non-zero ﬂow under ¹Cº. Each node in 0has 430 anoutgoing edge list that lists the probability that a unit of ﬂow travels down an 431 outgoing edge. Finally, each component (i.e., edge or non-demand node) 2in 0 432 hasamapthatliststhedemandnodesreachablefrom 2in 0. Anentryinthismap 433 contains a tupleh3%¹3j2ºithat gives the probability that a unit of ﬂow passing 434 through2reachesdemandnode 3.Equation4 canbecomputedfromthesemaps. 435 15 The CFC calculation proceeds by performing a reverse depth-ﬁrst search 436 (“DFS”) on 0for each demand node 3. When invoked, the DFS computes 437 theprobabilitythateachedgeornon-demandnodein 0sendsﬂowto 3. Agiven 438 node or edge (identiﬁed by an IDnumber) may be visited multiple times in the 439 courseofthesearch,requiringcaretoavoidpushingsuperﬂuousprobability. (This 440 method does not, however, work for graphs 0that contain cycles). 441 Function ComputeProbabilities(G’) Data:G’, a graph with components ¹+ºand absorbing nodes +. foreach32do ReverseSearch (G’, d) end Function ReverseSearch(G’, d) Data:G’, as above. Data:d, an absorbing node. Var excess[] // array of numbers 2»01¼of sizej+j Var stack excess[d.ID] = 1 stack.push(d) whilestack not empty do Var curNode = stack.pop() Var amt = excess[curNode.ID] // probability to push foreachincoming edge curEdge of curNode do curEdge.map.IncrementOrAddProbability(d.ID, amt) excess[curEdge.src.ID] = amt * curEdge.probability stack.push(curEdge.src) end curNode.map.IncrementOrAddProbability(d.ID, amt) excess[curNode.ID] = 0 end Algorithm 1: Probability Calculation. TheIncrementOrAddProbability () function adds an amount of probability to 442 the estimate of %¹3j2ºstored in the map of component 2, creating the map entry 443 if it does not exist. The variable 4G24BBis a lookup table (indexed by node ID) 444 that contains probability values for each node. (The variables excessandamtare 445 designed to deal with overlapping paths). 446 16 In the worstcase, the map at each nodeor link stores jjentries–one for each 447 demand node in 0. Also in the worst case, every edge 42in has ﬂow and 448 is in 0, yielding$¹¹j+j¸jjºjjºin storage space. Since the ReverseSearch 449 procedureisamodiﬁedDFS,thetimerequiredtocomputeprobabilitiesforagiven 450 demandnodeis $¹j+j¸jjº. Thus,calling ReverseSearch onalljjdemandnodes 451 in 0yields$¹¹j+j¸jjºjjº. Infrastructure networks typically have j+jjj 452 andjj/1 2j+j, givingAlgorithm 1 time and space complexity of $¹j+j2º. 453 Therunningtimeoftheentiremethodisthusdominatedbytheﬂowgeneration 454 step, which is typically more expensive than $¹j+j2º. The current paper uses the 455 Edmonds-Karpalgorithm(see[68])forsimplicity,whichis $¹j+j2jjºongeneral 456 graphs and$¹j+j3ºon infrastructure networks. Although ﬂows can be generated 457 using a variety of techniques (e.g., simulation), the method in Algorithm 1 only 458 appliesifthetransitiongraph 0isacyclic. Alternativemethods(e.g.,simulation, 459 Markov chains) can be used if cycles are present. 460 3.4. An Algorithm for Interdependent Critical Flow Centrality 461 Given a modelSwith interdependent sub-systems (1(2(3(=,Algo- 462 rithm 1can be used to compute CFC values for all components in each (8at 463 eachtimestep C. Thisisnotsuﬃcient,however,asphysicaldependenciesmustbe 464 accounted for. Resource demands and criticality values must be propagated from 465 one sub-system to the other. 466 Computation of the CFC for the entire model Sproceeds by computing the 467 CFCforeachindividualinfrastructuresystemintopologicalorder. Dependencies 468 are processed from one system to the next in each iteration, passing demands 469 from higher-level layers to lower-level ones. Algorithm 2 provides a high level 470 overview: 471 17 Function ComputeInterdependentCFC( G) Data:G, a graph with nodes +G=(=f(1(2(:grepresenting individual infrastructure systems, and edges Gformed from physical dependencies between elements of +G. ConvertNetworkRepresentation (G) Var list TopologicalSort (+G) Var t 0 foreach(82listdo ComputeSingleSystemCFC ((8) end Algorithm 2: Computing CFC for a set of interdependent infrastructures. 3.4.1. Converting Network Representations 472 As a pre-processing step, conversion of network representations is performed 473 to transform each individual network (8into a format compatible with maximum 474 ﬂow algorithms. 475 1. Nodes with demands are connected to a supersink node (see [68, 69]). 476 2. Source nodes are connected to a supersource node. 477 3. Nodesinnetwork (1thatrequireresourcesfromnetwork (2arerepresented 478 in(2by corresponding demand nodes. 479 Forstep(3),thecriticalityforthenodesin (1isonlyavailableaftertheCFCforall 480 non-demandnodeshasbeencomputed. Thus,thefullcomputationfor (1mustbe 481 performed before any computations can be performed for (2.Figure 3 provides 482 an illustration of network conversion. 483 3.4.2. Computing CFC Values for a Sub-system 484 Computation of the CFC for sub-system (8proceeds in two stages: (1) ﬂow 485 valuesandcriticalityvaluesarepropagatedfromotherlayers ((8)according 486 todependencies,and;(2)theCFCfor (8iscomputedusingthetechniquediscussed 487 inSection3.3 . Iflayer(8supplieslayer (withresources(e.g.,itisanelectricity 488 networkthatsuppliespowertowaterpumps),thenresourcedemandsfor (appear 489 in(8’s network as sinks with appropriate demands. Topological ordering ensures 490 that(’scriticalityandﬂowvalueshavebeencomputedbefore (8’s.Algorithm3 491 provides an overview of single layer CFC computation. 492 Propagation of criticality and ﬂow values proceeds by examining the set of 493 relevant interconnection records: 494 18 d2 n2 S1 n3 S2 n1 d1 d3 n5 n1 S1 n4 n3 n2 d1 d2pump G/one.denominator: G/two.denominator: d4 d2 n2 S1 n3 S2 n1 d1 d3 n5 n1 S1 n4 n3 n2 d1 d2pump SS ss n6 ss SS Fig. 3.Transformation of ﬂow networks 1and 2with multiple sources ( (8) and sinks ( 38). Supersource (‘SS’) and supersink nodes (‘ss’) are added in the usual manner. An interconnection record '¹E1E2º(whereE12(E22(8) indicates a 495 physical (resource) dependency between systems (and(8(seeFigure 4). De- 496 mand and criticality values for E1must be propagated to E2before the maximum 497 ﬂow and CFC can be computed for (8. 498 Algorithm 4 gives an overview of this process. Criticality values are copied 499 directly,buttheamountofresourcethatmustbeprovidedby E2toE1isdetermined 500 by a function (e.g., the demand induced at E2is half of the ﬂow at E1). 501 19 Function ComputeSingleSystemCFC( (8) PropagateValues ((8) ComputeMaxFlow ((8) ComputeCFC ((8) Algorithm 3: Computing the CFC for a set of interdependent systems. Fig. 4.An Interconnection Record used for Interdependencies. Function PropagateValues( (8) foreachinterconnection record '¹E1E2ºdo ifE22(8+then E2.demand CalculateResultingDemand (¹E1E2º) E2.criticality E1.criticality end end Algorithm 4: Propagation of resource demands. Figure 5 shows a water system and electricity system that are interlinked in 502 two locations: pumps near the source of the water system are fed by electricity 503 nodes labelled and. A ﬂow solution was ﬁrst computed for the water system, 504 yielding ﬂows of 6063litres and 5973liters at the pumps. The induced demand 505 atnodesandoftheelectricitysystemarehalfoftheﬂow—namely, 3031and 506 2986units. 507 Notealsothatedgesandverticeswithnoﬂowareshowninblack. Theexistence 508 of such elements is an artifact of the Edmonds-Karp algorithm [69, 68] used in 509 thissimpleinstantiation,andonethatwouldbecorrectedbyusingdomain-speciﬁc 510 methods (e.g., hydraulic simulation [58]). 511 20 910 11 12 13 910 11 12 13 9 10 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 9 10 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 9 10 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 910 11 12 13 9 10 11 2300 200 1500 4000 1900 2100 180 2100 2100 1200 300 2000 900 2986 3031 4055 11061 91 225 55 75 2174 45 12380 3322 3159 15447 3127 8915 3441 4341 4321 7536 5536 66 66 175 55 11071 11116 3297 119 12370 4266 3444 12359 8961 6815 SOURCE 678910 678910 6 78910 678910 67 8910 678910 678910 6 78910 67 8 910 678910 678 9 10 6 78910 67 8910 678910 678910 678910 678910 67 8910 678910 678 9 10 6 7 8120 900 3000600 600110 600600 562150 2000 2000 1200 240 2380 115 6063 160 6063 2235 6163 2380 6100 6000 1000 5973 1000 42 1042 120 2145 1042 2045 6063 593 5985 5733 1000 635 1000 3140 100 115 1000 2045 1000 730 115 730 900 760 3110 1315 1600 562 1750 5140 3750 240 13263 6063 6000 90 40 PUMP PUMP SOURCEWATER ELECTRICITY ABFig. 5.Interdependent ﬂows. Demand values are in italics while ﬂow values are in regular font. Pumps in the water network are supplied with electricity by nodes A and B. Pumps require electricity proportional to half of their water ﬂow. Black edges/vertices have zero ﬂow. Figure 6 shows the CFC values for the same interdependent infrastructure 512 system under the same ﬂow solution. Criticality levels (ranging from 0to1) are 513 shown in white font for the buildings. (Lot criticality is ﬁxed at 002, and elided 514 for brevity). 515 Thankstothepropagationofbothﬂowandcriticalityvaluesfromonenetwork 516 to the next, the criticality of the water pumps is appropriately represented in the 517 criticality ratings of the electricity system. The electricity nodes andhave 518 inheritedcriticalityvaluesof 032and061fromthecorrespondingpumpvertices 519 in the water system; they require ﬂow of 3031and2986units, which the reader 520 can verify by inspection are half of the ﬂow values at the water pump. 521 21 0.36 0.34 4e-05 0.0032 2.4e-053.3e-050.019 2e-050.3 0.34 0.0870.34 0.160.7 0.160.13 0.087 0.3 0.0890.36 0.32 0.28 2.9e-052.9e-05 2.4e-05 4.9e-060.34 4e-06 0.0870.3 0.36 0.170.3 0.360.12 0.32 0.120.11 0.0230.45 0.320.0047 0.320.23 0.320.45 0.610.61 0.0720.61 0.0720.003 0.0750.23 0.320.043 0.32 0.58 0.120.07 0.469.6e-05 0.00011 0.0720.0720.140.140.069 0.0160.580.0930.081 0.082 0.54 0.089 0.023 1.00.32 0.613.9e-05 6 78910 678910 6 78910 678910 67 8910 678910 678910 6 78910 678910 678910 678 9 10 6 78910 67 8910 678910 678910 678910 678910 67 8910 678910 678 9 10 6 7 80.2 0.4 10.7 0.10.20.10.1 0.020.02 0.2 0.02 0.60.2 0.4 10.7 0.1 0.2 0.1 0.1 0.02 0.02 0.2 0.02 0.6 0.32WATER ELECTRICITY 12002000562120900 3000600 600110 600600 150 2000 230020001200200 1500 40001900 2100180 21002100 300 900SOURCE SOURCE0.230.230.58 6 78910 678910 6 78910 678910 67 8910 678910 678910 6 78910 678910 678910 678 9 10 6 78910 67 8910 678910 678910 678910 678910 67 8910 678910 678 9 10 6 7 8 0.615.3e-05 7.8e-05 303129861.0Fig.6.(Normalized)criticalﬂowcentrality,computedfromtheﬂowsinFigure5. Demandvalues are in italics, CFC values are in regular font, and criticality ratings for buildings are in white. The electricity nodes that supply the water pumps are given criticality ratings of 032and061and demands of 3031and2986via Algorithm 4. While most of the critical demand in the model is for the hospital (criti- 522 cality=1.0) and secondary school (criticality=0.6), the pumps create signiﬁcant 523 critical demand in otherwise non-critical regions of the model. Figure 6 show 524 that the electricity nodes supplying the pumps carry 16.7% and 8.7% of the total 525 critical ﬂow in the electricity network. It would be a poor decision to co-locate 526 electrical assets with water assets when both are carrying highly critical ﬂow. 527 22 4. Evaluation 528 This section demonstrates the method by means of a district-level model of 529 a city containing electricity and water systems. The simplicity of the model is 530 for explanatory purposes; it is possible to use the method on models of greater 531 complexity,providedthatthephysicalinterdependenciescreateadirected,acyclic 532 graph. 533 Each building/lot in the model is given: (1) a type(e.g., hotel); (2) a time 534 seriesrepresenting hourlydemandforwater ;(3)atimeseriesrepresenting hourly 535 demand for electricity , and; (4) a criticality rating in the interval»01¼. Time 536 series are assumed to give average hourly demands over a 24-hour day. However, 537 the method is general, and other scenarios could be supported, such as long-term 538 (i.e., decadal) investigation of urban growth and its eﬀect on capacity. 539 Time series data is assigned to buildings according to type (e.g., secondary 540 school, restaurant), while lots are assigned time series randomly drawn from a 541 library of typical residential demand curves. For simplicity, criticality ratings 542 and vertex/edge capacities are assumed to be static, although they could easily be 543 represented with their own time series. 544 Empirical data for diﬀerent types of buildings in summer was obtained from 545 several sources (e.g., water consumption data was sourced from the California 546 PublicUtilitiesCommission[65],electricitydatafromOntarioPowerGeneration). 547 Examples of water demand curves appear in Figure 7 below: 0.000.020.040.060.08 00:00 04:00 08:00 12:00 16:00 20:00 24:00 TimeValue Laundromat 0.000.020.040.060.08 00:00 04:00 08:00 12:00 16:00 20:00 24:00 TimeValue Hospital Fig.7.Hourlytimeseriesshowingwaterdemandsfromalaundromatandhospitaloveranaverage day. Thetimeserieshavebeennormalizedtocreateaprobabilitydistribution. ForuseintheCFC method, these distributions are scaled by average water usage per day. 548 CFCvaluesarecomputedforeachtimestep C2»1)¼byloadingtherelevant 549 time series data for Cand executing Algorithm 2 . An overview of the process is 550 23 providedin Algorithm5 . Uponterminationofthisprocedure,eachnodeandedge 551 intheinterdependentsystemhasasetofCFCvalues—oneforeachtimestep— 552 that can be used in statistical analysis. 553 Function ComputeCriticaltyForTimeSeries( G) foreachC2»1)¼do LoadDemands (G,C) ComputeInterdependentCFC (G) end Algorithm 5: Computing CFC on a system-of-systems with time-varying demands. Figure 8 shows a graph of CFC values for the water network’s edges over the 554 full 24-hour cycle: 555 0:00 1:00 2:00 3:00 4:00 5:00 6:00 7:00 8:00 9:00 10:00 11:00 12:00 13:00 14:00 15:00 16:00 17:00 18:00 19:00 20:00 21:00 22:00 23:00 TIME Fig. 8.CFC values for each edge in the water network. The edge with a constant criticality rating of 1is the lone edge incident to the 556 source/reservoir. In general, the edges with signiﬁcant criticality values tend to 557 remaincriticalthroughoutthe24-hourcycle,withinterestingbehaviourhappening 558 during the middle of the day. Low criticality nodes become more critical during 559 mid-day, when signiﬁcant water demand begins to push capacity constraints. 560 Incontrast,theedgesoftheelectricitynetworkdisplayamorestabledistribu- 561 tion. InFigure 9, one can clearly see that there are fewer intersections between 562 24 linesintheplotofelectricityedgecriticalityvalues. Theedgetothesinglesource 563 node again has a constant criticality rating of 1, and the ﬂuctuation in criticality 564 valuesofothermajoredgesismuchlesspronounced. Thisislikelyaconsequence 565 ofthefactthatthedemandontheelectricitynetworkdoesnottendtopushcapacity 566 constraints as much as the demand on the water network. 567 0:00 1:00 2:00 3:00 4:00 5:00 6:00 7:00 8:00 9:00 10:00 11:00 12:00 13:00 14:00 15:00 16:00 17:00 18:00 19:00 20:00 21:00 22:00 23:00 TIME Fig. 9.CFC values for each edge in the electricity network. To recap, Algorithm 5 results in a set of CFC values C¹2º, whereCis a timestepand 2isacomponent. Forinstance,theoutputforthewatersystemedges canberepresentedasamatrix 4 F0C4Ainwhichrowsaretimestepsandcolumns are edges: ILI4 F0C4A=  1¹41º 1¹42º 1¹43º  1¹4jjº  2¹41º 2¹42º 2¹43º  2¹4jjº  :¹41º:¹42º:¹43º  :¹4jjº 568 Onemajor issuenot addressedbyclassical workson networkcentrality(e.g., [9]) 569 is the choice of ranking method for component measures taken at diﬀerent times. 570 Themostintuitiveapproachtorankingthecomponentsistotakethe samplemean 571 ofeachcolumnandtosubsequentlyrankcolumnsindescendingorder. Thiswould 572 beanappropriatestrategyifeachrowofthematrixwasasamplefromthespaceof 573 assignments(i.e.,inaMonteCarloapproach)atagiventime C. However,therows 574 25 in the matrix are assessments of the system at diﬀerent points in time. The use of 575 descriptive statistical measures (e.g., average, variance) elides system dynamics. 576 The same is true of various other methods (e.g., spectral analysis, information 577 theory) that might be employed to analyze the matrix. 578 The choice of ranking approach is dependent upon the purpose of analysis. 579 Consider a long-term (e.g., multi-year) analysis that attempts to study the distri- 580 bution of critical ﬂow patterns in response to changing population densities and 581 land-use patterns. In such a setting, the long-term behaviour of the system is of 582 interest. 583 Figure 10 displays a situation in which criticality curves for two diﬀerent 584 components have the same integral but completely diﬀerent trends over time. For 585 along-term(decadal)analysisofinfrastructurecriticality,thecomponentwiththe 586 orange criticality curve is clearly the more important of the two. In this setting, 587 some form of trend-based, multi-variable time series analysis is required. 588 2019 2020 2021 2022 2023 2024 2025 2026 2027 2028 2029 2030 2031 2032 2033 2044 2045 2046 2047 2048 2049 2050 2051 2052 TIME Fig. 10.Two components with similar integrals but diﬀerent long term behavior. Since the time series in this paper represent average demands in a daily cycle , 589 components are ranked according to the integralof their CFC curve. Taking the 590 water network edges as an example, a cubic spline is deﬁned on the set of sample 591 pointsf¹ 1¹4º 2¹4º :¹4ºgcorresponding to edge 42. 592 An integral is calculated from the cubic spline (as shown in Figure 11 ) and 593 normalized by the maximum possible area "-_¹:1º=123=23 594 (recall that CFC values are already normalized, so that the maximum CFC at any 595 time step is 1). The result is then assigned to the edge 4as itsglobal CFC value 596 26  ¹4ºfor the entire time series. The set of all water network edges is then 597 ranked by sorting the edges according to their  values. 598 0:00 1:00 2:00 3:00 4:00 5:00 6:00 7:00 8:00 9:00 10:00 11:00 12:00 13:00 14:00 15:00 16:00 17:00 18:00 19:00 20:00 21:00 22:00 23:00 TIME Fig. 11.Computing the global CFC value for a given component 4. A cubic spline is overlaid on the CFC values for 4. The integral of the spline (blue) is computed and normalized by the total area. The same process is repeated for vertices, and for the other networks in the 599 system-of-systems. Becauseofthewayinwhichthecriticalﬂowcentralitymetric 600 is deﬁned, values for edges and vertices are commensurate, allowing a global 601 ranking of all components in the interdependent system. 602 5. Reliability 603 This section demonstrates that the CFC measure may be combined with stan- 604 dardapproachestonetworkreliability—namely,(1)edgereliabilitymeasures,and; 605 (2) “leave one out” failure analysis. 606 5.1. Edge Reliability 607 Anarbitrarynetworkmodelcanbeaugmentedbyaddinga reliabilityfunction A:!»01¼thatassignsedges 42areliabilityrating A¹4º2» 01¼[70]. One can combine this approach with CFC measures by creating a composite measure thatestimatesthejointreliabilityandcriticalityofacomponent. Forinstance,the (normalized) Unreliable Critical Flow (“UCF”) for an edge 42+is: 0*¹4º=0¹4º¹1A¹4ºº 27 where0¹4ºisthenormalizedCFCforedge 4. (TheUCFis‘normalized’since 608 values lie in the range »01¼, since0¹4ºandA¹4ºare both in»01¼.) Under 609 thismeasure,componentsareimportanttothedegreethattheyare: (1)unreliable, 610 and; (2) instrumental for the delivery of resources to critical locations. 611 The computation of the UCF measure can be accomplished with a slight 612 modiﬁcation to the algorithm for the CFC. Instead of a static value A¹4º, the 613 reliability rating for a network component 4can also be represented as a time 614 series'4=fA41A42A4:g. This allows the modeler to represent diﬀerent 615 processes (e.g., decreasing reliability of components over long time periods). 616 The UCF measures are computed for each timestep Cusing the CFC values 617 and reliability ratings at C. The end result is a matrix in which entry ¹89ºgives 618 the UCF values for each edge 49at timestep8. As in the case of the CFC, a cubic 619 splineisoverlaidonthevaluesforeachedge,creatinganunreliabilitycurve. After 620 computing the integral and dividing it by the maximum possible area, the global 621 UCF value for edge4is computed. 622 Geospatialdependencies betweeninfrastructurecomponentscanbeintroduced 623 into edge reliability analysis in a number of ways. For example, edges that are 624 co-located(e.g.,awaterpipeandelectricitypipesharingthesameservicetunnel) 625 couldbeforcedtosharethesamereliabilityrating. Co-locatedcomponentscould 626 also be assigned a reliability penalty that reﬂects the fact that component failures 627 are no longer completely independent. 628 5.2. “Leave One Out” Failure Analysis 629 CFC measures can also be used with a common form of reliability analysis 630 in which components are deliberately failed or degraded (e.g., by reducing their 631 capacity) in order to assess the eﬀects on the system. A component 4may have a 632 high CFC value under a given assignment, but it may be the case that if 4suﬀers 633 a (partial) failure there are other routes (i.e., fallbacks) through which ﬂow may 634 travel in order to satisfy critical demand. Figure 12 illustrates this situation. 635 Thisformoffailureanalysisprovidesanindicationofwhethertherearefallback 636 routes that can supply critical ﬂow in the event that a component 4fails. If the 637 failure of4consistently results in reduced critical ﬂow across the entire network, 638 one can assume that 4is even more critical than the CFC measure alone might 639 suggest. 640 Algorithm 6 shows a high-level view of a procedure in which capacities of edgesinaninfrastructuresystemaredegradedone-at-a-time. Foreachtime C ), appropriate demands and criticality values are loaded into the graph. Then each edge42isconsideredinorder,degradingitscapacityandperformingtheCFC 28 S n3 d2 d1 n2 n1 e1e2e3cr(d 2)=0.8 cr(d 1)=0.1A: S n3 d2 d1 n2 n1 e1e2e3cr(d 2)=0.8 cr(d 1)=0.1B: n1Fig.12.Twonetworkswithdiﬀerentbehaviourinedgefailurescenarios. NetworkAcarriesmost of its critical ﬂow through the path f414243g. In case of edge failure, no alternative paths are available. Network B has a fallback route in case edge 43fails. computation on the altered network. The critical ﬂow is then used to create a loss measure that indicates the amount of critical ﬂow that is lost when edge 4is degraded. The failure loss!C¹4ºfor edge42at timeCis: !C¹4º=1Í 32+5¹3Cº2A¹3CºÍ 32+X¹3Cº2A¹3Cº where (recalling Section 3.2.2) +is the set of demand nodes in ,X¹3Cºis the 641 demandattime Cfromdemandnode 3,5¹3Cºistheactualﬂowto 3attimeC,and 642 2A¹3Cº2» 01¼isthecriticalityratingfor 3atC. Failurelossvaluesrangefrom 0 643 (no eﬀect on resource delivery) to 1(absolute disruption of resource delivery). 644 Function PerformEdgeFailureAnalysis( () foreachC2»1)¼do LoadDemands ((,C) foreach42do Var originalCapacity e.capacity e.capacity Degrade(e.capacity) ComputeSingleSystemCFC (S) ComputeFailureLoss (S) e.capacity originalCapacity end end Algorithm6: Edgefailureanalysisonnetwork (withtime-varyingdemands. The edge failure mechanism was tested on the network from Figure 5 by 645 degrading the capacity of each edge 4to0. (Demands and criticality ratings were 646 29 thesameasinprevioussections.) Thefailureloss !C¹4ºwascomputedforeach 647 edge4at each time C2 »023¼and averaged over the 24 hour cycle to create 648 an aggregate failure loss metric. The CFC values for each edge 4were likewise 649 averaged over the same time frame. 650 Figure 13 shows both the averaged CFC and averaged FL metrics for the 651 edges of the water network. Two facts are immediately obvious. First, the vast 652 majorityofedgeshavenegligibleaverageCFCandFLvalues. Thesearetypically 653 low-capacityfeedsfromaresidentialstreet’swaterpipetoanindividuallot/parcel. 654 CFCFL Fig.13.Averaged failureloss (FL)andaveraged criticalﬂowcentrality (CFC)ontheedgesofthe water network in Figure 5, computed over a 24-hour period. The majority of edges (e.g., those that feed individual lots) have negligible FL and CFC values. Second, a signiﬁcant percentage of of those edges with high CFC ratings also 655 have low FL values. Although these pipe segments carry a sizable amount of 656 critical ﬂow, alternative routes are available in case they should suﬀer individual 657 failures. Examplesincludethepipesthatdeﬁnetheloopsaroundresidentialblocks; 658 these loops are resistant to individual failure, since there are two paths from the 659 entry point of the loop to any lot/parcel. 660 Ofcourse,themainpipesfromthereservoirhavenobackups,asdemonstrated 661 bytheoverlapofFLandCFCvaluesforedge320. Ingeneral,thecorrelationofFL 662 30 andCFCismildlysigniﬁcantbutalsosomewhatmisleadingasasummarystatistic. 663 With a diﬀerent network topology that included multiple sources and alternative 664 paths, onewouldexpect lesscorrelationbetween theFL andCFCvalues, making 665 the easily computable CFC a poor predictor of the consequences of edge failures. 666 5.3. Reliability Integration Limitations and Assumptions 667 The edge failure analysis presented above was subject to several simplifying 668 assumptions. First, geospatial dependencies were not included in the analysis 669 for reasons of brevity. Second, the failure loss analysis is performed on each 670 subsystem independently, given a ﬂow solution for the entire system-of-systems; 671 handlinginterdependencies requiresiterativemethods. Third,the reliabilitymea- 672 sures could also incorporate component capacity, in order to capture the intuition 673 that a component nearing its maximum load is likely to be less reliable. 674 Fourth, the deﬁnition of the FL metric uses the aggregate of all demands at 675 the network’s demand nodes as the normalizing factor. This is appropriate for a 676 networkwherealldemandsaresatisﬁedinthebaselinestate,butitwilloverestimate 677 losses in networks which exhibit unsatisﬁed demand. For the scenario utilized in 678 this paper, however, this assumption is reasonable. 679 This method outlined in this work does not assume that the methods used to 680 modeleachlayerarecommensurate. Thatis,theelectricitylayermaybemodeled 681 withonesetofdomain-speciﬁctechniques,whilethewaterlayermaybemodeled 682 with another. All that is required is for each layer to provide a means of ﬂow 683 computation and a basic network topology. This design decision, while useful 684 fromasoftwareengineeringperspective,precludestheuseofstandardapproaches 685 to modeling cascading failures. 686 Tomodelphysicaldependenciesandcascadingfailureinsuchasettingrequires 687 theuseofadditionalmachinery. Componentfailureintheelectricitysystemcould 688 result in reduced power levels at the water pumps; this, in turn, could alter water 689 distribution ﬂows, result in reduced electricity demand from other components of 690 the water system—thereby changing demand patterns for the electricity system. 691 Thus,theresultisan equilibriumproblem inwhichchangesinonelayerpercolate 692 through other layers, and then back again. Solving such a problem is well beyond 693 the scope of this paper. 694 6. Conclusion 695 This paper demonstrated how component importance measures based on the 696 notion of critical ﬂow may be applied to interdependent, urban infrastructure 697 31 systems. Themotivationfortheworkwastoprovideurbanplannersandmunicipal 698 engineers with a method of reasoning about the impacts of interventions on the 699 ﬂow of resources to critical locations. The main theme was that network analysis 700 techniques could be combined with criticality and reliability metrics in order to 701 produce composite methods that provide useful information to stakeholders. 702 The perspective of the method was resource-based, focusing on the ways in 703 which system components participate in the delivery of resources. Each individ- 704 ual infrastructure system (8of a composite system Swas represented as a ﬂow 705 networkwithdemands,capacities,supplylimits,andcriticalityratings. Thepaper 706 considered physical dependencies in which one subsystem (8requires resources 707 from another subsystem (9. 708 In the simple variant described in the paper, network ﬂows and the critical 709 ﬂow centrality (“CFC”) measure were computed using a discrete approach. More 710 sophisticatedvariantsarepossible,includingtheuseofdomain-speciﬁcsimulation 711 techniques. For simplicity, the paper assumed that the subsystem dependencies 712 form a directed acyclic graph. 713 The method was demonstrated by use of a simple, district-scale model of a 714 city that contained electricity and water networks. Empirical data was used to 715 estimate resource consumption for diﬀerent types of buildings, yielding a set of 716 demand curves that represent consumption in a 24-hour cycle. This decision 717 simpliﬁed the analysis, and allowed the use of integrals to compute a global CFC 718 value for the entire cycle. For the study of trends in infrastructure systems over 719 time, the integral-based aggregation would need to be supplanted by trend-based, 720 multi-variable time series analysis. 721 Despitethesimplifyingassumption,thesimplemethodpresentedinthepaper 722 satisﬁed the goals outlined in Section 3.1 . First, the computation of CFC metrics 723 for an interdependent system can be computed eﬃciently. For a model S= 724 f(1(2(:gconsistingof :subsystems,computationofCFCmetricsfor Son 725 typicalinfrastructurenetworksis $¹+2º,where+istheaveragenumberofnodes 726 inthesubsystems. Thiscomparesfavorablywithothercentralitymeasures,which 727 can be$¹+3ºor greater. 728 Second, the demonstration showed that the basic method correctly propagates 729 resource demand, criticality ratings and CFC values between systems. Not only 730 areCFCvaluescomparableacrosscomponentswithinagivensystem,buttheyare 731 commensurableacrosssystems—evenincaseswheredisparatemodelingmethod- 732 ologies have been used. 733 Third, the paper showed how common network reliability approaches can be 734 combinedwithCFCmeasurestoyieldcompositemetrics. Edgereliabilitycanbe 735 32 directlyintegratedintotheCFCframeworkbyaddinganotherattributetotheedges 736 andtweakingtheCFCcomputationslightly. Thepaperalsodiscussededgefailure 737 analysis, showing that a composite failure loss metric can be deﬁned that gives 738 an indication of the availability of fallback routes for the delivery of resources to 739 critical locations. 740 Manyavenuesoffutureworkremain,themostimportantofwhichisremoving 741 the restriction ofGto directed, acyclic graphs. To do so invites consideration of 742 equilibrium concerns—that is, changes in one network cause changes in others, 743 altering ﬂow distributions and demand patterns in complex ways. Providing 744 solutions for this type of problem is well outside the scope of the present paper. 745 The version of the CFC computation presented in this paper is suitable for 746 medium/long time horizons only, as it uses integer-valued representations for 747 both demands and capacity constraints. This decision, which was made in order 748 to simplify the problem and avoid numerical instability, means that short-term 749 dynamics are diﬃcult to represent. As a consequence, rates of change (e.g., of 750 ﬂow)onsystemcomponentscannotbeanalyzedprecisely. Theuseofﬂoatingpoint 751 representations and domain-speciﬁc ﬂow computation methods (e.g., simulation) 752 will avoid this limitation. 753 Evenwiththisrestrictioninplace,therearestilladditionalissuestoberesolved. 754 First, a more realistic ﬂow mechanism (e.g., domain-speciﬁc methods) should 755 replace the generic Edmonds-Karp algorithm that favors shortest paths (thereby 756 introducing artifacts into the ﬂow solution). Second, geospatial dependencies 757 shouldbeintroducedintoboththeedgereliabilityandcomponentfailureanalyses. 758 Additional avenues of future research were hinted at throughout the paper. 759 7. Acknowledgements 760 Funding for this work was provided by an Ontario Research Fund — Re- 761 searchExcellenceRound7grantforthe“iCity: UrbanInformaticsforSustainable 762 Metropolitan Growth” project. The author wishes to thank Eric Miller for the 763 opportunity to work on this project. 764