Designing low -cost networks is an essential step in planning linked infrastructure. For the case of capacitated trees, such as oil or gas pipeline networks, the cost is usually a function of both pipeline thickness (i.e. capacity) and pipeline length. Minimizing cost be comes particularly difficult as network topology itself dictates local flow material balances , rendering the optimization space non -linear. The combinatorial nature of potential trees requires the use of graph optimization heuristics to achieve good soluti ons in reasonable time. In this work we perform a comparison of known literature network optimization heuristics and metaheuristics , and propose novel algorithms, including a metaheuristic based on transferring edges of high valency nodes . Our metaheuristi c achieves performance above similar algorithms studied , especially for larger graphs , usually producing a significantly higher proportion of optimal solutions , while remaining in line with time -complexity of algorithms found in the literature. Data points for graph node positions and capacities are first randomly generated , and secondly obtained from the German emissions trading CO 2 source registry. Driven by the increasing necessity to find applications and storage for industry CO 2 emissions, finding mini mum -cost networks increases the business case for large -scale CO 2 transportation pipeline infrastructure. Keywords: Network Design, Heuristics , Cost Minimization, Pipeline networks, CO 2 transport Declarations Funding : The Helmholtz Climate Initiative (HI -CAM) is funded by the Helmholtz Association's Initiative and Networking Fund. The authors are responsible for the content of this publication. Conflicts of interest/Competing interests : There are no conflicts to declare. Availability of data and material : The data used to test the codes here, as well as the results, are provided at : https://git ext.gfz-potsdam.de/yeates/network -design/ - /tree/master/paper_materials_valency_shuffle/ Code availability : The code for all the algorithms used here, as well all data analysis and plotting codes are available at : https://gitext.gfz -potsdam.de/yeates/network -design/ 3 1. Intro duction Network design problems arise in a multitude of scenarios. For fluid transport problems, such as oil or gas transport, the aim is to design networks linking all fixed location points of the network via pipelines in the way that meets system requirements while minimizing a defined cost -function. System requirement s may take different forms such as maximal pressure limitations for safe use, mini mal pressure requirements for ease of transport ( e.g. supercritical CO 2), prohibited transport through certain areas (e.g. natural parks) , or preferential transport through other areas ( e.g. existing pipeline corridors) . While some networks are created ite ratively as new nodes are added (i.e. historical domestic gas distribution networks) , potentially resulting in suboptimal network designs, networks can be designed from scratch as new incentives or legislative restrictions arise. Creating novel networks th at re -purpose existing liquid commodity sources, such as CO 2 storage or utilisation from industry CO 2 emissions (van den Broek et al. 2009 ; Alhajaj and Shah 2020) , slurry collection from dairy plants (Bietresato et al. 2013) , or otherwise networks that dis tribute a new commodity to a series of known recipients such as Hydrogen (André et al. 2013) or biogas networks (Heijnen et al. 2020) , all undergo a step of network design . This involves linking a series of fixed sources/sinks points via pipelines with associated capacity and cost. Finding optimal networks connecting sources to sinks is an example of finding Minimum -cost Capacitated Spanning Trees (MCST) , i.e. the minimum cost tree linking all elements of the network and respecting flow requirements . In this work, t he cost is given by a combined function of the pipeline length and a positive capacity value . Furthermore, as the relative pipeline cost is expected to flatten with increasing capacity , i.e. considering the economy of scale, a concave cost function of the capacity is used. While finding Minimum cost spanning trees without capacity is known to be solvable in polynomial time with algorithms such as Kruskal’s (1956) algorithm , the addition of the capacity const raint adds computational complexity. Such problems are difficult to optimize as the attribution of capacity (and therefore cost) to each pipeline will depend on the exact network topology chosen, creating a highly non-linear optimization space , where switc hing a single edge may alter all the attributed capacities of the network , significantly alter ing the cost. As a result , local minima are easy to find and discovering the optimal MSCT requires iterating over all trees , attributing capacity, then calculating associated cost. The combinatorial nature of trees renders this problem NP -hard. The further addition of intermediate relay nodes, Steiner points, to the network is expected to reduce overall costs, creating Minimu m cost Gilbert trees (MCGT)(Gilbert 1967). As such , each MCST serves as starting point for the introduction of Steiner points. Steiner points do not require any capacity addition to the pipelines and simply provide a way of reducing the overall pipeline le ngths . It has been shown (Heijnen et al. 2020) that the cost of the Minimum cost Gilbert tree (MCGT) is highly dependent on the MCST from which it originates and is expected to reach lower final costs for lower cost input MCSTs . In this regard, finding M CSTs which are as close as possible to the optimal MCST also becomes a valuable undertaking the search for the more general MCGT . Due to the large number of possible trees ( Cayley’s (1857) formula shows that 𝑛𝑛−2 spanning t rees exist for 𝑛 network nodes) , iteration over all possible spanning trees in hope of finding the optimal MCST rapidly becomes impractical. Instead, graph optimization schemes are used to achieve lower cost trees rapidly, to the detriment of not insuring true optimal soluitions . Multiple approaches exist, including multi -Mixed Integer Linear Programming , that make use of traditional solver tools from mathematical programming and optimization (van den Broek et al. 2009 ; Brimberg et al. 2003 ; Sun & Chen 2016 ), agent -based methods such as Ant Swarm Optimization (Maier et al. 2003) , Genetic Algorithms (Weihs et al. 2011) , 4 or graph -based heuristic algorithms (Kazmierczak et al. 2009 ; Bietresato et al. 2013 ; André et al. 2013 ; Heijnen et al. 2020) . A formal mathematical description of the problem can be found in Heijnen et al. (2020) or Xue et al. (1999) . Heuristic algorithms offer the advantages of providing quick, good results while allowing decision -makers and other stakeholders to understand them easily , due to their conceptual simplicity. Finally, they are easily modifiable to integrate extra system requirements due to transparency and code shortness . This work will therefore focus on e stablishing a baseline of MCST heuristic algorithms in combination with our proposed algorithms with associated performance statistics . Comparisons to optimal MCST solutions are given where possible (i.e. when not overly computationally deman ding). 2. Materials and methods 2.1. General procedure In a first step we give a non -exhaustive review of some commonly used applicable heuristics and metaheuristics and introduce some of our own . Heijnen et al. (2020) do a more comprehensive comparison of applicable algorithms, but here we only pick the best performing, with some additions from our own literature review. To do so, we make use of discrete optimization terms and concepts to describe some of the literature algorithms. At the same time, we describe the concepts and algorithms in an accessible way that can be unde rstood by scientists and stakeholders in network design alike. In a second step , we apply the algorithms directly to a series of examples with differing numbers of graph nodes , comparing the minimum -cost network found for each algorithm with an optimal solution, whe n possible. After providing some performance statistics based on a random generation of data points, some one-to-one algorithm comparison results are given using existing cluster s of CO 2 industry emission sources in Germany provided by the European Emissions trading scheme data (German Environment Federal Office 2018) . The general procedure we use to obtain a minimum -cost capacitated spanning tree is given in Fig. 1 . This procedure was also used by Yeates et al. ( 2020) and Heijnen et al. (2020) . A series of sources are initially connected by the unique minimum -spanning tree to a chosen sink , capacity is attributed to the edges, and cost is calculated as a starting point before entering the network optimization algo rithms that we explore in this paper. The flow requirement of the single sink is set by the opposite of the sum of the flow requirements of the sources . The location of the sink is chosen either randomly , along with the positions and flow requirements of the sources (Section 3.1) , or as a set of multiple positions in a regular grid surrounding the sources (Section 3.2) . The capacity attribution is unique for a given network topology and matching total sink -source capacity. This is calculated using the capacity allocation procedure proposed by Heijnen et al. (2020 ). The cost function for a pipeline used throughout this work is a simple continuous concave function of flow capacity 𝐶 as determined by the capa city attribution algorithm, and is proportional to the pipeline length 𝐿 given by: 𝐶𝑜𝑠𝑡 (𝑝𝑖)=𝐿𝐶0.6. The network cost is given by the sum of all pipeline costs. The exponent used here 0.6, is typical of CO 2 networks (Kazmierczak et al. 2009) and integrates notions of economy of scale as it is lower than one. The methods in this paper are applicable to any exponent between 0 and 1. The cheapest networks will represent the best trade -off between overall pipeline length reduction and joining of pipelines into higher capacity , comparatively cheaper pipes. 5 Fig. 1 Procedure for finding a minimum -cost network. Sources are shown in red, while the sink is shown in purple. Inversing requirements and thus sources and sinks yield the equivalent problem 2.2. Algorithms In this section we provi de a series of historical and recent example s of heuristics and metaheuristics to compare against our algorithm proposal s. As seen in Fig. 1, the network optimization algorithms (second last to last tiles in Fig. 1) take a tree with associated cost as input and attempt to find lower - cost capacitated trees through modifications of the tree topology, adaptation of the pipeline capacity requirements where needed, and recalculation of cost. Topology modi fications typically involve breaking the tree at an existing edge and recombining the two disconnected parts of the initial tree via a previously unconsidered edge. Such a transformation (i.e. breaking and recombining once) is known a s 1-neighbourhood tran sformation . These transformations can be described by a formal metric, a distance between solutions, most simply given by the symmetric difference between the sets of edges . As such, a transformation of distance 1 occurs when a n edge is removed, and the distinct parts of the tree are combined in a new manner. A transformation of distance 2 occurs when two such transformations are done simultaneously, providing a solution within the 2 -neighbourhood of the initial solution. Local search heuristics involve only 1 -transformations in a single algorithm step. As soon as a lower - cost solution is found, a first -descent local heuristic would then use the lower -cost solution as the new incumbent solution and repeat the process until no further single step improvement is made, i.e. a local minimum is found. A highest -descent local heuristic scans the entirety of its available 1- neighbourhood , then choosing the cheapest lower -cost solution and restarting the process until no further single step improvement is made. While a heuristic provides a way of making a transformation between two solutions, a metaheuristic considers more general strategies for reaching optimal solutions within the broader optimization space, by attempting carefully chosen jumps into non -immediate neighbourhoods by either invoking memor y, allowing exploratory moves towards higher -cost solutions or doing searches at different starting points in the optimization space and converging towards a single solution. We now provide three exa mples of local heuristics that serve as basis for exploring metaheuristics later, also serving as baselines for comparison to our contributions . Python source code is given for each algorithm in the provided Github repository . 2.2.1. Heuristics a. Delta Change The Delta Change (DC) algorithm was first introduced by Rothfarb et al. (1970) focusing on the optimal design of offshore gas pipelines. The procedure consists of taking an initial spanning tree solution, for example the MST, creat ing new candidate solution s by choosing two disconnected nodes, adding an edge between them to create a cycle in the structure, then breaking the cycle 6 elsewhere by removing one of the pre -existing edges, before recalculating required pipeline capacities and cost . As soon as a better s olution is found, including within the iterations over a given cycle, the algorithm is restarted, categorizing the procedure as a first-descent heuristic . A similar algorithm is described by Andre et al. (2003) although the authors include a randomization of the list of nodes to initiate a cycle from, and limit the number of explored node pairs by only selecting the closest ones to the initial node choice. b. Local Search The Delta Change cycle -based heuristic is a lso used in an algorithm described by Brimberg et a l. (2003) referred to as “Local Search” (LS). However, in Local Search, the steepest -descent form is used , in which all the possible cycles of the 1-neighbourhood are explored before the minimum -cost solut ion is chosen, if an improvement is found. c. Edge Turn Recently, another 1 -neighbourhood search heuristic was proposed by Heijnen et al . (2020) that we label the Edge Turn (ET) algorithm . The author s describe edge turns , which involve choos ing an edge to re move, and reconnect ing the disconnected two parts of the tree via a new edge connecting to one of the initial two nodes. As the breaking and recombin ing is done in the same step , the elementary edge -turn process can be compared to the one used in Delta Change, as each edge turn can be retrieved through an associated corresponding cycle. However, here the number of allowed transformations is reduced as the recombination of the two parts of the tree must necessarily be done at one of the nodes of the remove d edge . Heijnen et al . scan over all possible edge turns before choosing the lowest -cost solution, in a steepest -descent fashion . 2.2.2. Extensions into 2 -neighbourhood - Nested heuristics. Extensions of the above heuristics can be conceived to permit jumps into the 2 -neighbourhood of a graph solution in a single algorithmic step. This can allow discovery of lower cost solution s, inaccessible via a 1 -neighbourhood local search due to existence of local minima . Put simply , such algorithms enable discovery of bette r solutions that involve breaking and recombining the trees twice simultaneously, rather than twice in succession, in which case the first iteration would not lead to a lower -cost solution and th e second iteration would hence not be explored. These algorit hms can be described as nested versions of the 1 -neighbourhood heuristics. The computational complexity of such algorithms increases vastly as another local heuristic (i.e. into the 2 -neighbourhood) is performed for each tentative initial local heuristic (i.e. into the 1 -neighbourhood). We shall include 2- neighbourhood extension s of the Delta Change and Edge Turn algorithms named Nested Delta Change (NDC) and Nested Edge Turn (NET) . A nested version of the Local Search algorithm was attempted but rapidly disregarded as the long computational times were considered impractical. In the case of the Nested Delta Change algorithm, t he first-descent characteristic is maintained , such that as soon as a better solution is f ound, either in the initial 1 -neighbourhood move or any of the associated 2-neighbourhood moves , the algorithm is stopped and restarted at this incumbent solution. For the Nested E dge Turn algorithm , the steepest -descent characteristic is maintained over the full 2-neighbourhood, such that all 2 -neighbourhood solutions are tried, and only the best new solution (if any) is chosen as the incumbent . 2.2.3. Metaheuristics Metaheu ristics are more general search strategies that orient the lower -level heuristics in the context of the larger optimization space, mainly by providing mechanism s to overcome local minima. We give 2 literature examples of metaheuristics, Variable Neighbourhood Search and Tabu Search, and finally provide our contribution algorithm, the High Valency Shuffle. 7 a. Tabu Search Tabu search is a general method proposed by Glover et a l. (1989) and adapted to an oil pipeline design problem by Brimberg et a l. (2003) . It is based on exploring solution s around a local minimum that ma y not immediately be lowe r-cost, but add the least amount of cost to the local minimum . The lowest cost solution is still stored and ultimately chosen if no better solution is found, but exploration to other zones is permitted . In the case of a local minimum, the 1 -transformation that increases the cost the least is chosen as the solution from which the next local search shall be initiated. Th is transformation (or exchange of graph edges) is then added to the Tabu l ist, and the reverse transformation is prohibited by the algorithm . The process is then repeated until a user -defined stopping criterion is reached, and the overall minimum cost solution is returned. The stopping criterion used here is chosen as 10 times the duration of the previously defined Local Search, and the Tabu l ist has a maximal length of 7, identical choices to Brimberg et al. (2003) . b. VNS Variable neighbourhood search was initially proposed by Mladenović & Hansen (1997) and was used by Brimberg et al. (2003 ). Starting from an initial local minimum solution, it mak es tentative jumps into solutions in further away neighbourhoods in a single step , reaching exploratory solutions that can have increased cost versus the current best solution , before attempting a local search in those new location s to find lower -cost s olution. In the simple form proposed by Brimberg el al. implemented in this work , random jumps into new neighbourhoods are performed, successively reaching points at larger and larger radii (or distances) from the local minimum until a better solution is found. N ote that radius is defined by the specific metric used, here the metric described in the introduction. If no better solutions are found up to a chosen neighbourhood radius (here up to 5 consecutive 1 -neighbourhood transformations , as per Brimberg’s implementation ), the algorithm restarts at lower distance jumps again. The algorithm stops when a user -defined stopping criterion is reached. The stopping criterion used here is chosen as 10 times the duration of a Local Search, as used by the authors. c. High Valency Shuffle The High Valency Shuffle (VS) proposed here is based on the observation that many local minima solutions share the characteristic of having at least one high va lency node (>2 edges) that is distinct to a high valency node found in the optimal solution. For a given high valency node i n the optimal solution, a node situated in proximity of it may play a similar role in distributing flow but for a higher cost. The l atter solution represents a local minimum , and metaheuristics and heuristics alike can often be trapped in it. Indeed, none of the local transformations can transfer all the edges in a single step to another node , and incremental transfers over multiple st eps do not necessarily provide lower -cost solutions. In some cases, only a full s ingle -step s witch of high valency nodes can reach a lower cost solution. While this transformation may be accessible via Variable Neighbourhood Search or Tabu search as they explore the wider optimization space outside of the local minima , they are sometimes too general to create such a node switch in a re asonable number of steps. The High Valency Shuffle metaheuristic is used in combination with a lower -level local heuristic such a Local Search, Delta Change or Edge Turn. The High Valency Shuffle is described in Table 1. 8 Table 1 High Valency Shuffle Metaheuristic Step 0: Start from an initial good solution obtained with a local heuristic such as Local Search, Delta Change or Edge Turn. This solution is the starting incumbent solution Step 1: Initiate empty solution list Step 2: Identify all the nodes 𝑛𝐻𝑉 with high valency (3 and above edges) Step 3a: For each node 𝑛𝐻𝑉𝑖 of these 𝑛𝐻𝑉: - Identify closest (in Euclidian distance) nodes 𝑛𝐶 (e.g. 4 closest ones), connected or not to 𝑛𝐻𝑉𝑖 Step 3b: For each node 𝑛𝐶𝑖 of these 𝑛𝐶: - Transfer all the edges from 𝑛𝐻𝑉𝑖 to 𝑛𝐶𝑖 - Connect 𝑛𝐻𝑉𝑖 to 𝑛𝐶𝑖 if not already done Step 3c: If a cycle is detected in the graph: For each of the edges in the cycle: - Tentatively remove edge from cycle - Run a lower -level local heuristic on resulting tree - Add the obtained solution to the solution list Else: - Run a lower -level local heuristic on the tree - Add the obtained solution to the solution list Step 4: Find minimal -cost solution from solution list Step 5: If this solution is better than current incumbent solution: - Set this new solution as incumbent solution - Restart algorithm from Step 1 Else: - End algorithm and return incumbent solution In this study, we use three lower -level local heuristics (seen in Step 3 c) in combination with the High Valency Shuffle: Edge Turn (VSET) , Delta Change (VSDC) , and Local Search (VSLS) , all described previously. An example of a transformation obtained in a single step with a High Valency Shuffle is shown in Fig. 6a . 3. Results 3.1. Random graphs To evaluate the p roposed heuristics and metaheuristics and compare with literature options, we initially perform tests on generated data, with randomly placed source and sink locations as well as source capacities to gauge the algorithm’s effectiveness on a variety of input data. We use the aforementioned algorithms on an increasing number of sources and compare calculation t imes to the fraction of optimal networks found. For the smaller network comparisons, as displayed in Fig. 2, optimal networks were discovered through a brute force method , obtained by comparing all the possible trees, obtained from Prüfer sequences (1918) . For larger network comparisons, shown in Fig. 3, due to impractical computational times for the exhaustive method , the best -found tree for all the algorithms was used for the opti mal comparison . In each random graph the sources and sink were placed randomly on a square and source capacities C were attribute d randomly according to the following law : 𝐶=𝑋3 where X is a random uniformly distributed variable between 0 -100. The sink cap acity was then set as the opposite of the sum of the source capacities to ensure material balance. The minimum spanning tree was used as a starting point for each algorithm , as shown in Fig. 1. 9 Fig. 2 Algorithm performance comparison for small network s wit h a direct comparison to guaranteed optimal network We make a few observation s regarding the performance of the various algorithms on the small networks shown in Fig. 2. There exists a general trend in which increased calculation time leads to increased rates of optima l solutions. All the new algorithms proposed here ( Nested Edge Turn and Nested Delta Change , and all High Valency Shuffle Algorithms ) lead to better optimal ity at the expense of longer calculation times. The High Valency Shuffle Algorithms obtain a 100% optimalit y rate in most of the cases while being faster than the Nested Algorithms . As the number of sources increase s, calculation time for all algorithms in crease s, and the optimality of all algorithms decreases except for the High Valency Shuffle Algorithms and the Nested Algorithms . 10 Fig. 3 Algorithm performance comparison for larger networks. The asterisk refers to the fact that the optimal network comparis ons are not guaranteed to be globally minimal As the number of sources increase, we continue to observe a loss of optimality for all algorithms, except the High Valency Shuffle Algorithms , which have optimali ty rates upwards of 90% in most cases. The Neste d Algorithms , despite their large calculation times , fail to capture the best network solutions in the larger networks. These will therefore be excluded in some of the case study examples as they are impractical. 3.2. Case study: grid search procedure over CO 2 sources As a case study of the algorithms, we use clusters of German industry CO 2 source locations and emission quantities as initial input data. As well as using these sources, we create a grid of potential sink points surrounding the clusters and iterat e the network optimization algorithms over each potential sink location. As such, we both test the speed and efficiency of the algorithms in finding the optimal graph solution and get an indication as to where the optimal sink location is placed. We presen t the CO 2 source clusters used as input data for the sink grid search. The CO 2 source clusters considered are shown by different colours in Fig. 3. They are numbered from 1 to 5 in order of increasing number of points. Source clustering was done using a DB SCAN (Ester et al. 1996) algorithm. We limit ourselves to only a few emission clusters to demonstrate the algorithm effectiveness. The sources in grey without a cluster number are therefore not considered here initially . Emission quantities were taken as a 3-year average of years 2015 -2017 for emitters listed in the EU Emissions Trading Scheme database and georeferenced to provide coordinates. Some further alterations to the dataset were made . Emissions whose primary activity was energy production were remo ved, as these sources (mainly coal plants) have an uncertain future in the German emission landscape with the planned phasing out of large coal -fired plants before 203 8 at the latest (German Federal Ministry for Economic Affairs and Energy 2020 ). Through this omission, 73% of the total emission quantity is not considered. Finally, individual sources emitting annually less than 50000 tons of CO 2 per year were excluded to reduce the numbe r of points in the potential networks . These 11 represent a further reduction of 3% in total emission quantity . In some cases, different emission sources occur in the same location, and are listed as different contributors. For the network optimization proced ure, sources with the same coordinates are combined as a single node with emissions summed together. Fig. 4 Left: CO 2 source clusters considered i n this work. Right: lowest -cost networks obtained for each cluster, with optimal sink location shown as a red cross Cluster details can be found in Table 2. Table 2 Details of CO 2 source cluster s Cluster number Geographical region Number of sources Number of nodes (with sink) 1 Berlin -Brandenburg 7 5 2 Bremen 8 8 3 Hamburg 12 9 4 Southern Germany 28 26 5 Central Germany 50 38 The potential sink locations were chosen within a grid surrounding the overall shape of the cluster up to 50 km around the convex hull of the shape described by the sources in the cluster. The number of sink locations was adapted each time to the network complexity to account for large calculation times . The algorithms and number of sink points tried for each cluster are given in Table 3. The algorithms considered are Edge Turn (ET), Delta Change (DC) , Local Search (LS), Nested Edge Turn (NET), Nested Delta Change (NDC), Tabu Search, Variable Neighbourho od Search (VNS) , High Valency Shuffle with Edge Turn (VSET) , High Valency Shuffle with Delta Change (VSDC) , High Valency Shuffle with Local Search (VSLS) . For the three large st clusters, some of the algorithms required extremely large calculation times. For the clusters 3, 4 and 5 the nested algorithms were therefore omitted , also owing to their lower performance demonstrated in the previous section . Finally, a comparison to a k nown optimal solution was only possible for the smallest cluster . Iterating over all possible trees rapidly becomes computationally too impractical for other clusters . For the 12 other clusters, the optimal tree considered was simply taken as the minimum cost tree for all the algorithms but is not guaranteed to be globally minimal. Table 3 Numerical details for the grid search around the source cluster s Cluster number Algorithms omitted Number of sink locations considered Comparison to a guaranteed global minimum? 1 ∅ 1666 Yes 2 ∅ 1193 No 3 ∅ 1076 No 4 NET,NDC 125 No 5 NET,NDC 65 No Grid search results are summarized in Fig. 5 with colours consistent with clusters s hown in Fig. 4. Fig. 5 Grid search results comparison. Only for Cluster 1 is the comparison to a guaranteed global minimum possibl e Through these grid searches we see that the conclusions obtained with the random graphs are maintained, notably the graphs solutions obtained from the Valency Shuffle algorithms are near - universal ly optima l. Their performance is particularly superior for t he largest clusters , such as cluster 5, and well -above all other literature algorithms considered with similar calculation times such as Tabu and VNS. One valuable comment concerns the choice of lower -level heuristic to use in combination with the high valency shuffle. A fast but mid -fidelity lower -level heuristic (such as the Edge Turn algorithm) achieves almost equivalent optimality rates when combined with the High Valency Shuffle metaheuristic (VSET) as better performing , but slower, lower -level heur istics such as Local Search (combined as VSLS) . In other words, differences in performance of the lower -level heuristic do not seem significantly affect the result obtained when used in combination with the High Valency Shuffle metaheuristic. In this regar d, it becomes useful to employ a fast lower -level heuristic (such as the Edge Turn algorithm) in combination with the High Valency Shuffle for the best combination of 13 optimality and speed. This leads us to speculate that an even faster lower -level heuristic exploring a reduced fraction of the 1 -neighborhood might be sufficient to provide similarly high optimality rates in combination with the high valency shuffle and further decrease the overall calculation time. In the Appe ndices we display the detailed results of the grid searches for all sink locations , presented on top of the original sources . The results are compared through maps showing the location of non - optimal solutions and deviation from optimal cost. 4. Discussion 4.1. High v alency metric – another route for improvement. As described in the Materials and Methods section , the High Valency Shuffle heuristic switches the edges of high valency nodes to their closest neighbours and eliminat es (if any) created cycles before applying a local heuristic and then choosing the lowest cost solution. This metaheuristic is an interplay between two distinct metrics when used with the lower -level heuristics shown here . The first is a classical metric based on the symmetric difference betw een sets of edges. This metric is used by the lower leve l heuristic that breaks single edges and reassembles trees once at each step. The higher -level metaheuristic transfers multiple edges in one step. It can be nonetheless described as making 1 -transformations with in the metric space described by the symmetric difference between the sets of nodes with more than 2 edges . A typical transformation allowed by the High Valency Shuffle metaheuristic with Local Search is shown in Fig. 6a. In this example, taken from the grid search of Cluster 2, the potential sink is represented by a black circle and the sources are coloured as displayed in Fig. 4. The sink location remains the same within both panel s. The graph to the left of panel a) is the result obtained from the Local Search . The graph to the right of panel a) is the result obtained from the High Valency Shuffle with Local Search , which is found to have a cost 5% lower. The edges of a high valency node are transferred to a different node , as expected. To obtain this graph from the first, the number of 1 - transformations in the symmetric difference of edges metric is thre e. Within the proposed high valency metric, this is equivalent to a local 1-transformation . Indeed, the set of nodes possessing more than two edges has changed by one element . However, i n panel b), the graph to the left is a suboptimal result obtained from the High Valency Shuffle with Local Search . The graph on the right is a lower cost solution found with the Nested Delta Change algorithm. We can therefore see that the High Valency Shuffle was unable to add an extra element in the set of high valency node s (i.e. two high valency nodes rather than one) as it instead relies on switching each high valency node to another node. The transformation seen in panel b) is nonetheless a 1-transformation in the high valency metric space but is not captured by the High Valency Shuffle metaheuristic . A number of suboptimal solutions provided by the High Valency Shuffle metaheuristic can be included in this category of transformations . The inclusion of a heuristic that c ould test selected candidate nodes to become high v alency nodes without transfer from a prior high valency node would be beneficial for achieving higher optimality rates. 14 Fig. 6 Final graphs obtained using different algorithms. In both cases the transformation to the lower cost graph involves a local 1 -transformation in the proposed high valency metric space 4.2. Time complexity analysis and application to large graphs The calculation times from randomly generated graphs in the previous section are displayed as a function of the number of sources in the network. The straight lines seen in the log -log plot are indicative of polynomial time algorithms . The curves are the refore fitted to a power function of the form: 𝑡(𝑁)=𝐴𝑁𝐵, where 𝑡(𝑁) is the calculation time in seconds, given in terms of the number of sources 𝑁. 𝐴 and 𝐵 are the variables to fit . To achieve this, the curves are fitted by minimum least - squares to a linear function in log -log space, with a weigh t function given by 𝑤=√𝑡. The weighting is included to avoid the intrinsic bias of fitting a transformed function in log -space in which a disproportional weight is given to the points with smaller x -axis values. The results from the curve - fitting are shown in Table 4. Table 4 Power law fit parameters for the time complexity analysis: The algorithms’ performance is fitted with the function 𝑡(𝑁)=𝐴𝑁𝐵 with t in seconds LS DC ET Tabu VNS DCN ETN VSLS VSDC VSET A 2.28 ×10−4 2.09 ×10−4 0.64 ×10−4 4.80 ×10−4 8.82 ×10−4 2.05 ×10−4 0.51 ×10−4 25.82 ×10−4 0.76 ×10−4 2.34 ×10−4 B 3.80 3.69 3.78 4.39 4.19 5.75 5.70 4.07 5.15 4.42 15 Fig. 7 Time complexity comparison of the algorithms described in this stud y Finally, a test is done on a large graph composed of CO 2 sources previously not included in the clusters from the previous section . Within this CO 2 source cluster, there are 143 emitters at 94 distinct locations. This cluster is seen at the top left of Fig. 7. Combining the sources from the same location gives a network of 95 nodes , including the sink. Due to the large number of nodes in the cluster, only a single sink location was considered and placed at the centre of gravity of all the sources, weighted by emission quantity. All algorithms were tested but stopped if not completed before a reasonable period of time (48h). The lowest -cost graph was only obtained using the High Valency metaheuristic in combination with the Edge Turn heuristic. Predicted calculation times, measure d calculation times obtained relative graph costs, when applicable, are shown in Table 5. The cost difference of the network created from the Minimum Spanning Tree (MST), which simply minimizes network length, is also given for comparison. Table 5 Details of the network design modelling of the large final CO 2 source clu ster Algorithm MST LS ET DC Tabu VNS VSET Predicted time (h) 0 2.03 0.51 1.14 61.22 45.36 33.51 Measured time (h) 0 2.09 0.31 0.93 25.89 24.89 25.24 Cost difference +33.6% +0.19% +2.6% +0.94% +0.03% +0.19% 0 To better visualise the results from this final example, the spatial networks are shown in Fig. 8. The compared algorithms are Local Search (LS), Variable Neighbourhood Search (VNS), Delta Change (DC), Edge Turn (ET), Tabu Search (Tabu) and High Valency Shuffle with Edge Turn (VSET). A comparison with the solution obtained from the Minimum Spanning Tree (MST) is also given. The lowest -cost network obtained with the VSET algorithm is shown with pipelines in grey. Network differences with the lowest -cost network for other solutions are shown in red, whereas sections of the networks shared with the lowest -cost solution are also left in grey. For clarity, the sources are 16 displayed with a small, constant size, despite contributing differently to the net work flow. However, pipeline thickness remains representative of flow quantity. The sink node is shown as a black square. Fig. 8 Final comparison on algorithms with reasonable (< 48h) calculation times on a graph with 95 nodes 5. Conclusion In this paper we h ave conducted performance comparisons on graph optimization heuristics and metaheuristics aimed at achieving minimum -cost capacitated trees within the context of pipeline network design. In the search for lower -cost solutions, we attempt nested solutions o f literature algorithms that enable unoriented explorations into wider optimization spaces in a single algorithmic step. These solutions provided higher rates of optimal networks than literature options, at the expense of large r calculation times. A new me taheuristic and associated distance metric were developed that captured useful network transformations easily . Such transformations had required a large amount of algorithm steps previously and were often not accessed at all . The High Valency 17 Shuffle , switching the edges of nodes with more than two edges to neighbours, proved to achieve almost 100% optimal networks in most examples used. The calculation time s for the new metaheuristic , while large, remain comparable to current literature option s such as Tabu or Variable Neighbourhood search while achieving greater proportions of optimal solutions . The new metaheuristic is used in combination of a lower -level heuristic , the choice of which seems to influence only the resulting solution cost but h as a large effect on calculation time. It seems therefore beneficial to combine the High Valency Shuffle with a fast low -level metaheuristic which may by itself only achieve low rates of optimal solutions .