Engineering Negative Cycle Canceling for Wind Farm Cabling
Sascha Gritzbach, Torsten Ueckerdt, Dorothea Wagner, Franziska Wegner, Matthias Wolf
EEngineering Negative Cycle Canceling for WindFarm Cabling
Sascha Gritzbach
Karlsruhe Institute of Technology, Germany
Torsten Ueckerdt
Karlsruhe Institute of Technology, Germany
Dorothea Wagner
Karlsruhe Institute of Technology, Germany
Franziska Wegner
Karlsruhe Institute of Technology, Germany
Matthias Wolf
Karlsruhe Institute of Technology, Germanyall authors: fi[email protected]
Abstract
In a wind farm turbines convert wind energy into electrical energy. The generation of eachturbine is transmitted, possibly via other turbines, to a substation that is connected to the powergrid. On every possible interconnection there can be at most one of various different cable types.Each type comes with a cost per unit length and with a capacity. Designing a cost-minimalcable layout for a wind farm to feed all turbine production into the power grid is called the
Wind Farm Cabling Problem ( WCP ).We consider a formulation of
WCP as a flow problem on a graph where the cost of a flow onan edge is modeled by a step function originating from the cable types. Recently, we presented aproof-of-concept for a negative cycle canceling-based algorithm for
WCP [14]. We extend key stepsof that heuristic and build a theoretical foundation that explains how this heuristic tackles theproblems arising from the special structure of
WCP .A thorough experimental evaluation identifies the best setup of the algorithm and compares itto existing methods from the literature such as
Mixed-integer Linear Programming ( MILP )and Simulated Annealing ( SA ). The heuristic runs in a range of half a millisecond to approximatelyone and a half minutes on instances with up to 500 turbines. It provides solutions of similar qualitycompared to both competitors with running times of one hour and one day. When comparing thesolution quality after a running time of two seconds, our algorithm outperforms the MILP - and SA -approaches, which allows it to be applied in interactive wind farm planning. Mathematics of computing → Network flows; Mathematics ofcomputing → Graph algorithms; Mathematics of computing → Network optimization
Keywords and phrases
Negative Cycle Canceling, Step Cost Function, Wind Farm Planning
Related Version
A preliminary version of this work was published in the Proceedings of the27th Annual European Symposium on Algorithms (ESA’19).
Funding
This work was funded (in part) by the Helmholtz Program Storage and Cross-linkedInfrastructures, Topic 6 Superconductivity, Networks and System Integration, by the Helmholtzfuture topic Energy Systems Integration, and by the German Research Foundation (DFG) as partof the Research Training Group GRK 2153: Energy Status Data – Informatics Methods for itsCollection, Analysis and Exploitation.
Acknowledgements
We thank our colleagues Lukas Barth and Marcel Radermacher for valuablediscussions and Lukas Barth for providing a code interface for simultaneous
MILP -solver experiments. a r X i v : . [ c s . D S ] A ug Engineering Negative Cycle Canceling for Wind Farm Cabling
Wind energy becomes increasingly important to help reduce effects of climate change. Asof 2017, 11 . internal cabling ) where multiple turbinesmay be connected in series. The designer of a wind farm has various cable types available,each of which with respective costs and thermal capacities. The latter restricts the amountof energy that can be transmitted through a cable. Planning a wind farm as a wholeconsists of various steps, including determining the locations for turbines and substations,layouting the connections from substations to the grid point, and designing the internalcabling. The planning process comes with a high level of complexity, which automatedapproaches struggle with [24]. Therefore, one might opt for decoupling the planning steps.We call the task of finding a cost-minimal internal cabling of a wind farm with given turbineand substation positions, as well as given turbine production and substation capacities, the Wind Farm Cabling Problem ( WCP ). Since
WCP is a generalization of the
N P -hardproblem
Capacitated Minimum Spanning Tree [21], it is
N P -hard as well.Due to the overall cost of a wind farm, using one day of computation time or morearguably is a reasonable way to approach
WCP . Such computation times, however, are notappropriate for an interactive planning process: Imagine a wind farm planner uses a planningtool which allows altering turbine positions to explore their influence on possible cable layouts.In that case, computation times of at most several seconds are desirable.
We extend our recent proof-of-concept, in which negative cycle canceling is applied to aformulation of
WCP as a network flow problem (cf. Section 3) with a step cost functionrepresenting the cable types [14]. The idea of negative cycle canceling is to iteratively identifycycles in a graph in which the edges are associated with the costs of (or gains from) changingthe flow. Normally, a cycle of negative total cost corresponds to a way to decrease the costof a previously found flow. Due to the step cost function, however, not every negative cyclehelps improve a solution to
WCP . We explore this and other issues for negative cycle cancelingthat arise from the step cost function in the flow problem formulation for
WCP . We presenta modification of the Bellman-Ford algorithm [3, 9] and build a theoretical foundation thatexplains how the modified algorithm addresses the aforementioned issues, e. g., by being ableto identify cycles that actually improve a solution. This modification works on a subgraph ofthe line graph (cf. page 7) of the input graph and can be implemented in the same asymptoticrunning time as the original Bellman-Ford algorithm.We further extend that heuristic by identifying two key abstraction layers and applyingdifferent strategies in those layers. Using different initializations is hinted at in the section onfuture work in [14]. We follow this hint and design eight concrete initialization strategies. Inanother layer, we propose a total of eight so-called “delta strategies” that specify the orderin which different values for flow changes are considered.In [14] we compared the Negative Cycle Canceling (
NCC ) algorithm to a
Mixed-integer . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 3
Linear Program ( MILP ) using the
MILP solver Gurobi with one-hour running times onbenchmark sets from the literature [17]. We extend this evaluation by identifying the best ofour variants and by comparing its results to the results of
MILP experiments after runningtimes of two seconds, one hour, and one day on the same benchmark sets. A running timeof two seconds helps identify the usefulness of the
NCC algorithm to an interactive planningprocess. The other running times stand for non-time-critical planning. We also compare thealgorithm to an approach using Simulated Annealing [17] with different running times. The res-ults show that our heuristic is very fast since it terminates on instances with up to 500 turbinesin under 100 seconds. At two seconds our algorithm outperforms its competitors, making itfeasible for interactive wind farm planning. Even with longer running times for the
MILP - and SA -approaches, our algorithm yields solutions to WCP of similar quality but in tens of seconds.In Section 2 we review existing work on
WCP and negative cycle canceling. In Section 3we define
WCP as a flow problem. We give theoretical insights on the difference to standardflow problems and present and analyze our Negative Cycle Canceling algorithm in Section 4.An extensive experimental evaluation of the algorithm is given in Section 5. We concludewith a short summary of the results and outline possible research directions (see Section 6).
In one of the first works on
WCP , a hierarchical decomposition of the problem was intro-duced [4]. The layers relate to well-known graph problems and heuristics for various settingsare proposed. Since then, considerable effort has been put into solving variants of
WCP . Exactsolutions can be computed using
Mixed-integer Linear Program ( MILP ) formulations in-cluding various degrees of technical constraints, e. g., line losses, component failures, and windstochasticity [18]. However, sizes of wind farms that are solved to optimality in reasonabletime are small. Metaheuristics such as Genetic Algorithms [26, 6] or Simulated Annealing [17]can provide good but not necessarily optimal solutions in relatively short computation times.We applied negative cycle canceling to a suitable flow formulation for
WCP [14], butthere is still an extensive agenda of open questions such as investigating the effect of othercable types, a comparison to existing heuristics, and using the solution as warm start for a
MILP solver. Originally, negative cycle canceling is proposed in the context of minimum costcirculations when linear cost functions are considered [16]. The algorithm for the Minimum-Cost Flow Problem based on cycle canceling with strongly polynomial running time runsin O ( nm (log n ) min { log( nC ) , m log n } ) time on a network with n vertices, m edges, andmaximum absolute value of costs C [12]. The bound for the running time of this algorithmwas later tightened to Θ(min { nm log( nC ) , nm } ) [22]. Negative cycle canceling has also beenused for problems with non-linear cost functions. Among these are multicommodity flowproblems with certain non-linear yet convex cost functions based on a queueing model [20]and the Capacity Expansion Problem for multicommodity flow networks with certain non-convex and non-smooth cost functions [7]. A classic algorithm for finding negative cycles isthe Bellman-Ford algorithm [3, 9] with heuristic improvements [13, 11]. An experimentalevaluation of these heuristics and other negative cycle detection algorithms is given in [5].A step cost function similar to the one in WCP appears in a multicommodity flow problem,for which exact solutions can be obtained by a procedure based on Benders Decomposition [10].However, this procedure is only evaluated on instances with up to 20 vertices and 37 edges andsome running times exceed 13 hours. While our approach does not guarantee to solve
WCP tooptimality, our evaluation shows that the solution quality is very good compared to the
MILP with running times not exceeding 100 seconds on wind farms with up to 500 turbines.
Engineering Negative Cycle Canceling for Wind Farm Cabling
The model presented in this paper is based on an existing flow model for
WCP [14]. Webriefly recall the model. Given a wind farm, let V T and V S be the sets of turbines andsubstations, respectively. We define a vertex set V of a graph by V = V T ∪ V S . For anytwo vertices u and v that can be connected by a cable in the wind farm, we define exactlyone directed edge e = ( u, v ), where the direction is chosen arbitrarily. We obtain a directedgraph G = ( V, E ) with V = V T ∪ V S and E ⊆ ( V × V ) \ ( V S × V S ) such that ( u, v ) ∈ E implies ( v, u ) / ∈ E . There are no edges between any two substations since we consider the windfarm planning step in which all positions of turbines and substations, as well as the cablingfrom substations to the onshore grid point have been fixed. We assume that all turbinesgenerate one unit of electricity. Note that our algorithm can be easily generalized to handlenon-uniform integral generation. Substations have a capacity cap sub : V S → N representingthe maximum amount of turbine production they can handle and each edge has a length givenby len : E → R ≥ representing the geographic distance between the endpoints of the edge.A flow on G is a function f : E → R and for an edge ( u, v ) with f ( u, v ) > < f ( u, v ) units of flow go from u to v (resp. − f ( u, v ) units go from v to u ). For a flow f and a vertex u we define the net flow in u by f net ( u ) = P ( w,u ) ∈ E f ( w, u ) − P ( u,w ) ∈ E f ( u, w ).A flow f is feasible if the conditions on flow conservation for both turbines (Equation (1))and substations (Equation (2)) are satisfied and if there is no outflow from any substation(Equations (3) and (4)). f net ( u ) = − ∀ u ∈ V T , (1) f net ( v ) ≤ cap sub ( v ) ∀ v ∈ V S , (2) f ( u, v ) ≥ ∀ ( u, v ) ∈ E : v ∈ V S , (3) f ( v, u ) ≤ ∀ ( v, u ) ∈ E : v ∈ V S . (4)Let c : R ≥ → R ≥ ∪ {∞} be a non-decreasing, left-continuous step function with c (0) = 0.This function represents the cable costs and sup { x ∈ R ≥ : c ( x ) < ∞} is the maximum cablecapacity, which we assume to be a natural number. Note that such a function is neitherconvex nor concave in general. The cost of a flow on a wind farm graph is then given bycost( f ) = X e ∈ E c ( | f ( e ) | ) · len( e ) . (5)The value of c ( | f ( e ) | ) stands for the cost per unit length of the cheapest cable type withsufficient capacity to transmit | f ( e ) | units of turbine production. With all that, WCP isthe problem of finding a feasible flow f on a given wind farm graph that minimizes thecost. There is an analogon to the linear-cost integer flow theorem (e. g. [2, Thm. 9.10]) thatguarantees an optimal flow with integral values. (cid:73) Lemma 1.
Suppose the cost function is discontinuous only at integers and there is afeasible flow. Then, there is a cost-minimal integral flow.
Proof.
Suppose f is a (possibly non-integral) flow of minimum costs. We define anotherflow network on the same graph by setting the capacity cap( e ) of every edge e to d| f ( e ) |e .Each turbine requires a net flow of −
1. We model the substation capacities by adding anew vertex s and edges from all substations to s with capacities equal to the substationcapacities. The net flow shall be 0 at all substations and | V T | at s . We further define zerocosts for flows on all edges. By the integrality property of min-cost flow problems with linear . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 5 cost functions (e. g., [2, Thm. 9.10]) there is a feasible integral flow f in this network. Dueto the construction of the flow network, f satisfies the constraints in Equations (1) to (4).Since the cost function c is non-decreasing, it holds for all e ∈ E that c ( x ) ≤ c ( | f ( e ) | )for all x ∈ [0 , | f ( e ) | ]. Since c is a left-continuous step function that is discontinuous onlyat integers, we also have c ( x ) ≤ c ( | f ( e ) | ) for all x ∈ [ | f ( e ) | , cap( e )]. It holds in particularthat c ( | f ( e ) | ) ≤ c ( | f ( e ) | ). Thus, cost( f ) ≤ cost( f ) and f is optimal in the original network. (cid:74) Given a wind farm graph G we define the residual graph R of G with vertices V ( R ) andedges E ( R ) by V ( R ) = V ( G ) ∪ { s } and E ( R ) = { e, ¯ e : e ∈ E ( G ) } ∪ { ( v, s ) , ( s, v ) : v ∈ V S } where ¯ e is the reverse of e . The new vertex s , the super substation , is a virtual substationwithout capacity, that is connected to all substations. The edges to and from s are used tomodel the substation capacity constraints and to allow the production of one turbine to bereassigned to another substation.For a given feasible flow f in G of finite cost and ∆ ∈ N we further define residualcosts , which represent by how much the cost for the edge changes if the flow on the edgeis increased by ∆ (cf. Figure 1 (a) – (d) for an example). Note that for negative quant-ities of flow this implies that the absolute value of the flow may be reduced or even thedirection of the flow on an edge may change. More formally, we define γ : E ( R ) → R by γ ( e ) = (cid:0) c ( | f ( e ) + ∆ | ) − c ( | f ( e ) | ) (cid:1) · len( e ) for all e ∈ E ( R ) that are neither incident to s nor lead to a substation where we alias f (¯ e ) = − f ( e ) for all e ∈ E ( G ). By this definition theresidual costs are infinite if c ( | f ( e ) + ∆ | ) = ∞ , i. e., if the maximal capacity on e is exceeded.For u ∈ V S and v ∈ V T , we set γ ( u, v ) = ∞ whenever f ( v, u ) < ∆ because sending f ( u, v )+∆units from u to v would otherwise imply that flow leaves a substation. On edges into s , weset γ ( u, s ) = 0 if and only if f ( u, s ) + ∆ ≤ cap sub ( u ) and γ ( u, s ) = ∞ otherwise. On edgesleaving the super substation, we set γ ( s, u ) = 0 if and only if f ( u, s ) ≥ ∆ and γ ( s, u ) = ∞ otherwise to prevent flow from leaving the substation.In a nutshell, the Negative Cycle Canceling ( NCC ) algorithm (Algorithm 1) starts withan initial feasible flow and some value of ∆, computes the residual costs, and looks for anegative cycle in the residual graph. If the algorithm finds a negative cycle, it cancels thecycle, i. e., it changes the flow by adding ∆ units of flow on all (residual) edges of the cycle.Note that this may decrease the actual amount of flow on edges of G . Then this procedureis repeated with the new flow and some value of ∆ which may but not need to differ fromthe previous one. If no negative cycle is found, a new value of ∆ is chosen and new residualcosts are computed. This loop is repeated until all sensible values of ∆ have been consideredfor a single flow, which is then returned by the algorithm. This flow is of integer value,since the initial flow is designed to only have integer values and we solely consider naturalvalues for ∆. Without loss of generality we can restrict ourselves to integer flows accordingto Lemma 1, even though our algorithm does not necessarily find an optimal solution of WCP .One question we answer is to what extent the algorithm benefits from different initial flowsand different orders in which the values of ∆ are chosen. We present various initializationsand orders for ∆ in Sections 4.3 and 4.4. A cycle is a sequence of consecutive edges such that the first edge starts at the same vertex where thelast edge ends and such that no two edges start at the same vertex. That is, all cycles are simple. Acycle is said to be negative if the sum of residual costs over all edges is negative. Engineering Negative Cycle Canceling for Wind Farm Cabling (d) (f) (h)(a) (c) (e) v v v v v v u uv v v u ∞ -220 ∞ ∞ -2 v v v u v v v u (b) (g) v v v u v v v u ∞ -2 31 11 112 3 11 Figure 1
Examples of flows and corresponding residual graphs. (a) shows a wind farm graph.Edges between turbines are of length 2, edges between the substation u and any turbine are oflength 3. (b) depicts a cost function induced by two cable types. (c) displays a feasible flow. Dashedlines do not carry any flow. The thickness of solid lines represent the necessary cable type to carrythe respective flow. (d) is the residual graph for the flow in (c) and ∆ = 1. The super substationis omitted for ease of presentation. There are three negative cycles: uv u , uv v u , and uv v u .(e) shows the flow obtained by sending one unit of flow along uv v u in (c). (f) is the residual graphfor (e) and ∆ = 1. (g) depicts the flow obtained by sending one unit of flow along uv v u in (c).(h) displays the residual graph for (g) and ∆ = 1. The details of the algorithm (cf. Sections 4.1 and 4.2) address problems that arise from thespecial structure of
WCP , namely the non-linear cost function c . Firstly, in classical min-costflow problems, when c is linear, the cost for changing flow by a certain amount is proportionalto the amount of flow change (and the length of the respective edge) and does not depend onthe current amount of flow on that edge. Hence, there is no need for computing residual costsfor different values of ∆. Secondly, short cycles, i. e., cycles of two edges, may have non-zerototal cost in WCP (cf. cycle uv u in Figure 1 (d)). Canceling such a cycle, however, doesnot change the flow and therefore does not improve the solution. Hence, only cycles of atleast three edges ( long cycles) are interesting to us because they do not contain both an edgeand its reverse. Finding any negative cycle can be done in polynomial time but finding longnegative cycles is N P -hard for general directed graphs [15, Theorem 4 for k = 3]. Thirdly,the order of canceling cycles matters (Figure 1 (c) – (g)). In (d), there are two long negativecycles: uv v u and uv v u . After canceling uv v u (Figure 1 (e)), the other cycle uv v u isnot negative anymore. Ultimately, Figure 1 (e) and (f) show that the non-existence of negativecycles in (all) residual graphs does not imply that the underlying flow is optimal—contraryto min-cost flow problems with linear cost functions. In other words, there are flows thatrepresent local but not global minima. We assume that the reader is familiar with the standard Bellman-Ford algorithm [3, 9], whichis a common approach to finding negative cycles. We observed in preliminary experimentsthat it mostly reports short cycles even if long cycles exist. The reason is that negativeresidual costs on an edge are repeatedly used if the cost of the reverse edge is, say, zero. Inthat case, the negative residual cost strongly influences the distance labels on close verticesand overshadows long cycles (see cycle uv u in comparison to cycle uv v u in Figure 1 (d)). . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 7 One solution is to prohibit propagating the residual cost of an edge over its reverse edge.To this end, we employ the Bellman-Ford algorithm on the subgraph L of the directed linegraph of R which we obtain from the line graph by removing all edges representing U-turns,i. e., edges of the form ( e, ¯ e ) for e ∈ E ( R ). We define the cost of an edge ( e , e ) in L as γ ( e ).At every vertex e of L we maintain a distance label ‘ ( e ) initialized as γ ( e ). Thus, throughoutthe Bellman-Ford algorithm, ‘ ( e ) represents the length of some walk in L starting at anyvertex of L and ending at e . By construction of L , the label ‘ ( e ) also stands for some walk in R which ends at the target vertex of e and which does not traverse an edge of R directly after itsreverse. Consequently, a cycle C in L corresponds to a closed walk W without U-turns of thesame cost in R . In particular, W is not a short cycle, which is what we wanted. It may stilloccur, however, that W includes an edge and its reverse. In that case, W consists of more thanone cycle that may be negative themselves. Therefore, we decompose the closed walk W intocycles, which, in turn, can be canceled one after another. For more details, refer to Section 4.2.A downside of running the Bellman-Ford algorithm on the line graph is that more labelshave to be stored and the running time of the algorithm is in O ( | V ( L ) | · | E ( L ) | ), which isworse than the running time on R . We present how to implement an algorithm that directlyworks on R , that is equivalent to the Bellman-Ford algorithm on L , and that has the sameasymptotic running time as the original Bellman-Ford algorithm on R . To this goal, we usethe special structure of L to analyze what the steps of the Bellman-Ford algorithm on L meanfor R . When running the Bellman-Ford algorithm on L , there is one label for every vertexof L . Each of those labels gives rise to a label on an edge of R . The labels at incoming edgesof v ∈ V ( R ) are used to compute the labels at outgoing edges of v . Let ( v, w ) and ( v, x ) betwo edges leaving v . Let us assume that ( x, v ) has the smallest label of all edges entering v .Then, ( x, v ) is used to relax ( v, w ). But it cannot be used to relax ( v, x ). To do so, we needthe second smallest label of all edges entering v . This yields the following observation. (cid:73) Observation 2.
For each vertex v of R only the two smallest labels of incoming edges of v are required to correctly update the labels on outgoing edges of v . We call these labels relevant . Consequently, throughout our modified version of theBellman-Ford algorithm, we maintain two distance labels ‘ ( v ) and ‘ ( v ), and two parentpointers parent ( v ) and parent ( v ) for every v ∈ V ( R ), respectively. As above, ‘ i ( v ) with i =1 , i ( v ) , v ). That means that the parent pointers hold the edges that have been used tobuild the values of the distance labels. The algorithm ensures that parent ( v ) = parent ( v )and ‘ ( v ) ≤ ‘ ( v ) for every v ∈ V ( R ). In every iteration of the Bellman-Ford algorithm,each edge of R is considered for relaxation: For an edge e = ( u, v ) take ‘ ( u ) = ‘ ( u )if parent ( u ) = v and ‘ ( u ) = ‘ ( u ) otherwise. Then, check if ‘ ( u ) + γ ( e ) yields a new relevantlabel at v . If, during a relaxation step, several incumbent labels and a newly computedcandidate label have the same value, we break ties in favor of the older labels—as in theoriginal algorithm. For each edge, checking if it yields a new relevant label at its end vertexcan be done in constant time. With Observation 2 we show reduced bounds for the number ofiterations and the overall running time compared to a straightforward implementation on L . (cid:73) Theorem 3.
If after · | V ( R ) | iterations there is an edge that allows reducing a label, thenthere is a negative cycle in L . The line graph L ( G ) of a directed graph G shows which edges are incident to each other. It is definedby V ( L ( G )) = E ( G ) and E ( L ( G )) = { (( u, v ) , ( v, w )): ( u, v ) , ( v, w ) ∈ E ( G ) } . A walk is a sequence of consecutive edges. A walk is called closed if the start vertex of the first edgeequals the target vertex of the last edge. In particular, every cycle is a closed walk. Engineering Negative Cycle Canceling for Wind Farm Cabling
Proof.
Let n = | V ( R ) | and suppose e n +1 is an edge that allows reducing a label after2 n iterations. We iteratively construct a walk backwards starting from e n +1 by repeatedlyapplying the following procedure. At an edge e i = ( v, w ) we define e i − as the incoming edgeof v other than ( w, v ) with the smallest label. If there are several possibilities, we pick the edgewith the oldest label among them. The label at e i − is relevant by definition. We stop whenan edge would be repeated. At this point, the walk contains a closed subwalk W = ( e k , . . . , e l )for suitable k, l ∈ Z with k < l ≤ n + 1. By Observation 2 there are at most 2 n edges withrelevant labels. Hence k ≥ e n +1 can be updated after 2 n iterations, the label at e n must havebeen updated in iteration 2 n . Repeating this argument inductively shows that for i ≥ k the label at e i was updated in or after iteration i . Therefore, all labels of edges in W wereupdated after the initialization. By the way the labels are computed, we therefore have ‘ ( e i − ) + γ ( e i ) ≤ ‘ ( e i ) (6)for all edges e i on W where we alias e k − = e l .If one of these inequalities is strict, i. e., ‘ ( e j − ) + γ ( e j ) < ‘ ( e j ) for some j ∈ { k, . . . , l } ,then summing over the inequalities for all edges in W will give X e ∈ W (cid:0) ‘ ( e ) + γ ( e ) (cid:1) < X e ∈ W ‘ ( e ) , (7)which can be simplified to X e ∈ W γ ( e ) < . (8)Hence, the total costs of W will be negative, which will complete the proof.It remains to show that there is some edge e j = ( v, w ) for which the inequality is strict.To this aim let e j be the edge with the oldest label among edges in W . The label ‘ ( e j ) wascomputed from the label ‘ ( e ) of an edge e = ( u, v ) with u = w , which may or may notbe e j − . That means ‘ ( e ) + γ ( e j ) = ‘ ( e j ) (9)where ‘ denotes the labels after the algorithm finishes and ‘ denotes the labels when ‘ ( e j ) iscomputed. Note that the label at e may have been updated afterwards, i. e., ‘ ( e ) ≤ ‘ ( e ).For the sake of contradiction assume ‘ ( e j − ) ≥ ‘ ( e ). Then, ‘ ( e ) ≥ ‘ ( e ) ≥ ‘ ( e j − ) ≥ ‘ ( e )where the first inequality holds since labels at the same edge do not increase during thealgorithm and the second inequality follows from e being an incoming edge of v . Hence, allthese labels are equal. The first equality implies that the label on e was not updated afterthe point in time when ‘ ( e j ) was computed and that ‘ ( e ) is older than ‘ ( e j ). Using thesecond equality, we distinguish two cases: If e = e j − , then ‘ ( e j − ) is older than ‘ ( e j ), whichcontradicts the choice of e j . If e = e j − , then e should have been included in W insteadof e j − . Thus, the assumption of ‘ ( e j − ) ≥ ‘ ( e ) is wrong and it holds that ‘ ( e j − ) < ‘ ( e ).Combining this inequality and Equation (9) completes the proof. (cid:74)(cid:73) Corollary 4.
A negative cycle in L can be computed in O ( | V ( R ) | · | E ( R ) | ) time if one exists. . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 9 The previously described Bellman-Ford algorithm on L is encapsulated in Algorithm 1. Wefirst compute some initial flow (line 1) using one of eight initialization strategies presentedin Section 4.3. In line 3 we compute the residual graph R using a given flow f and a given ∆and run the modified Bellman-Ford algorithm (line 4). In the repeat-loop, we consider oneedge after another and check in line 7 if it can be relaxed (again). In that case, we extract awalk W in R with negative costs leading to that edge by traversing parent pointers. However,canceling W directly may not improve the costs of the flow as W may still contain an edgeand its reverse. We decompose W into a set of simple cycles C in line 8 and cancel each cycleindependently if it is long and has negative costs (lines 9 to 12). Note that even though W has negative costs, it may happen that only short cycles in C have negative costs and all longcycles have non-negative costs. In this case we search for another negative cycle in L (line 13).If no negative cycle in the current graph L is canceled, a new value for ∆ is determinedaccording to the delta strategy (cf. Section 4.4) in line 14 and new residual costs γ arecomputed. Line 14 also checks if every possible value for ∆ has been used after the lastupdate of f without improving the solution. If so, f is returned. Algorithm 1:
Negative Cycle Canceling
Input:
Graph G , costs c , edge lengths len Result:
A feasible flow f in G f := InitializeFlow ( G, len), ∆ := InitialDelta while ∆ = N U LL do ( R, γ ) :=
ComputeResidualGraph ( G, c, f, ∆) RunBellmanFord ( R, γ ) found := false foreach e ∈ E ( R ) do W := FindNegativeClosedWalk ( R, e ) C := DecomposeWalkIntoCycles ( W ) foreach C ∈ C do if | C | ≥ and γ ( C ) < then f := AddFlowOnCycle ( f, C, ∆) found := true if found then break ∆ := NextDelta (∆ , found ) return f We apply two well-known speed-up techniques to the Bellman-Ford algorithm. Firstly, ifone iteration does not yield any update of any label, then the computation is aborted andno negative cycle can be found in the current residual graph. Secondly, after sorting edgesby start vertices, we track whether the labels at a vertex v have been updated since lastconsidering its outgoing edges. If not, then there is no need to relax the outgoing edges. Before we can start searching for and canceling negative cycles, we need some feasible initialflow. To obtain such a flow, we consider eight strategies, which all roughly work as follows.We pick a turbine u whose production has not been routed to a substation yet. We then search for a shortest path P from u to a substation v with free capacity using Dijkstra’salgorithm [8]. The search only considers edges on which the production of the turbine canbe routed, i. e., it ignores congested edges. We then route the production of u along P to v .We consider two metrics to compute shortest paths. Either we use the lengths definedby len (cf. Section 3) or we assume a length of 1 for every edge. Turbine production can eitherbe routed to a nearest or a farthest (in the sense of the respective metric) substation withfree capacity. There are two ways in which the flow is updated: The simpler variant routesonly the production of u along P , i. e., the flow along P is increased by 1. The other variantgreedily collects as much production from u and other turbines on P as possible withoutviolating any capacity constraints. The resulting flows are integral since the substationcapacities and the maximum cable capacity are natural numbers. If no feasible flow of finitecost is found during the initialization, the algorithm returns without a result.This yields eight initialization strategies, which we name as follows. The base part ofeach name is either BFS if unit distances are used or
Dijkstra (abbr.
Dijk ) if the distancesgiven by len are used. This part is followed by a suffix specifying the target substation:
Any (abbr. A ) for the nearest and Last (abbr. L ) for the farthest substation. An optionalprefix of Collecting (abbr. C ) means that the production is greedily collected along shortestpaths. For example, CollectingDijkstraLast (abbr.
C-Dijk-L ) iterates over all turbinesand for each turbine u it finds the substation v such that the shortest path given by lenfrom u to v is longest among all substations. Along a shortest path from u to v , turbineproduction is collected greedily. A delta strategy consists of two parts: an initial value for ∆ and a function that returns thevalue of ∆ for the following iteration. We discuss eight delta strategies. The simplest onestarts with ∆ = 1 and increments ∆ until a negative cycle is canceled. Then, ∆ is reset to 1.We call this strategy
Inc (as in increasing). Similarly,
Dec (as in decreasing) starts with thelargest possible value for ∆, which is twice the largest cable capacity. Then, ∆ is decrementeduntil a cycle is canceled and reset to the largest value. The third strategy
IncDec behaves like
Inc until a negative cycle is canceled. Then, it decrements ∆ until ∆ = 1 and behaves like
Inc again. To improve performance, all ∆ can be skipped during incrementation up to thelast value of ∆ for which a negative cycle was canceled. The fourth strategy
Random returnsrandom natural numbers between one and the maximum possible value for ∆. Between anytwo cycle cancellations, no value is repeated.For each strategy, we consider the following modification: After canceling a negative cycle,we retain the current value of ∆, recompute the residual costs with the new flow, and runthe Bellman-Ford algorithm again. We repeat this, until ∆ does not yield a negative cycle.In that case, ∆ is changed according to the respective delta strategy. We call the strategiesafter the modification
StayInc , StayDec , StayIncDec , and
StayRandom (or
S-Inc , S-Dec , S-IncDec , and
S-Random for short).
In the previous sections, we introduced a heuristic with various strategies for the
WCP . Wefirst use statistical tests to evaluate these strategies and identify the best ones (Section 5.1).Using the result we compare the best variant (i. e., best combination of initialization anddelta strategy) with different base line algorithms for the
WCP namely solving an exact
MILP formulation (Section 5.3) and a Simulated Annealing algorithm [17] (Section 5.4). In . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 11 preliminary experiments (Section 5.2) we determine which of the
MILP solvers Gurobi andCPLEX works better for
WCP to establish which solver we compare the
NCC algorithm to.For our evaluation we use benchmark sets for wind farms from the literature [17] consistingof wind farms of different sizes and characteristics: small wind farms with exactly onesubstation ( N : 10–79 turbines), wind farms with multiple substations ( N : 20–79 turbines, N : 80–180 turbines, N : 200–499 turbines), and complete graphs ( N : 80–180 turbines).Our code is written in C ++
14 and compiled with GCC 7 . . -O3 -march=native flags. All simulations are run on a 64-bit architecture with four 12-core CPUs of AMDclocked at 2 . .
0. All computations arerun in single-threaded mode to ensure comparability of the different algorithms.
In a first step, we want to determine which delta strategy works best. To this end, werandomly select 200 instances per benchmark set. We run our algorithm on each instance withevery pair of delta and initialization strategy. The first eight rows of Table 1 show for everybenchmark set the minimum, average, and maximum running times for each delta strategyacross all initialization strategies. We first observe that all variants are fast, with runningtimes between tenths of milliseconds to 4 . Dec is always the slowest strategy on average, which can be explained by the factthat
Dec often tries large values for ∆, for which negative cycles are found rarely. The otherstrategies all roughly complete in the same time on average. It seems to be slightly faster torepeat the same ∆. However, for our purpose all variants have small enough running times.We therefore base our decision, which variant to choose, solely on their solution qualities.To compare the variants in terms of solution quality, we compute for each delta strategy i and instance m the mean X ( i ) m of the solution values over all eight initialization strategies.This gives us 1000 data points per delta strategy. For delta strategies i, j we perform aBinomial Sign Test counting instances with X ( i ) m < X ( j ) m and X ( j ) m < X ( i ) m (Appendix A),that means for this test we are rather interested in whether strategy i performs better thanstrategy j on instance m and not by how much i is better than j on m . Table 2 summarizesthe results of all tests after Bonferroni-correction by 112 (the number of tests from bothdelta and initialization strategies). The percentage given in an entry in row i and column j states on how many instances i performes strictly better than j after averaging over allinitialization strategies. Note that entries ( i, j ) and ( j, i ) need not represent 1000 instances,as two variants may return equal solution values.In the row IncDec , all values are above 50 %, three of which are significant at the 10 − and another one at the 10 − -level. The smallest value (50 . StayIncDec ) standsfor 460 instances on which
IncDec performs better than
StayIncDec . To the contrary, thereare 446 instances on which
StayIncDec yields better solutions (cf. entry 49 . StayIncDec and column
IncDec ). While the differences between the four delta strategiesinvolving
Inc and
IncDec are not statistically significant,
IncDec does seem to have a slightadvantage over the others. Hence we consider
IncDec as the best delta strategy.In Figure 2 (left), for the dark green curve all instances are ordered by X ( Random ) m /X ( IncDec ) m in ascending order. For a given value α on the abscissa, the curve shows the relative costfactor of the instance at the α -quantile in the computed order. The other curves workaccordingly. We see, for example, that IncDec works strictly better than
StayInc on 49 . . . IncDec outperforms
Inc by at least 0 . Random ) and 0.947(
Inc ) and the maximum ratios are between 1.027 (
Random ) and 1.104 (
StayIncDec ). Table 1
Minimum, average and maximum of running times in milliseconds of different variants.Running time measurement starts before the initial flow is computed and ends with the terminationof the algorithm prior to outputting the solution. The first eight rows represent running times acrossall initialization strategies per delta strategies and benchmark sets. The best delta strategy in termsof solution quality is marked in green; minimal values per column are marked in yellow. The lastrow represents the algorithm variant
IncDec , CollectingDijkstraAny .DeltaStrategy N N N N N min avg max min avg max min avg max min avg max min avg max Dec
Inc
IncDec
Random
S-Dec
S-Inc
S-IncDec
S-Random
BestVar
Table 2
Comparison of delta strategies over all initialization strategies. An entry in row i andcolumn j shows on how many instances strategy i produces better solutions than strategy j . Valuesare marked by a star if they are significant with p < − and by two stars if p < − . The beststrategy is marked in green. Inc Dec IncDec Random S-Inc S-Dec S-IncDec S-RandomInc — 60 . ?? . . ?? . . ?? . . ? Dec . . . . . . . IncDec . . ?? — 59 . ?? . . ?? . . ? Random . . . . . . S-Inc . . ?? . . ? — 58 . ?? . . S-Dec . . . . . . . S-IncDec . . ?? . . ? . . ?? — 55 . S-Random . . ?? . . ? . . ?? . Table 3
Comparison of the initialization strategies when the delta strategy
IncDec is fixed. Anentry in row i and column j shows on how many instances strategy i produces better solutionsthan strategy j . Values are marked by a star if they are significant with p < − and by two starsif p < − . The best strategy is marked in green. Dijk-A BFS-A C-Dijk-A C-BFS-A Dijk-L BFS-L C-Dijk-L C-BFS-LDijk-A — 55 . . . . . . . ? BFS-A . . . . . . . C-Dijk-A . . ? — 55 . . . ? . . C-BFS-A . . . . . . . Dijk-L . . . . . . . BFS-L . . . . . . . C-Dijk-L . . . . . . . C-BFS-L . . . . . . . . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 13 Instances in % R e l a t i v e C o s t s IncS-IncDecS-DecIncDecS-IncDecRandS-Rand
Instances in % R e l a t i v e C o s t s Dijk-ADijk-LBFS-ABFS-LC-Dijk-AC-Dijk-LC-BFS-AC-BFS-L
Figure 2
Evaluation of the
NCC
Algorithm using different strategies. For each strategy and foreach instance, the ratio of the best solution value found by that
NCC variant to the best solution valuefound by the reference variant (marked in red) are computed. They are shown in increasing order.The dashed lines represent the 25% and 75% quantiles of the instances.
Left:
The delta strategies arepresented relative to the
IncDec strategy. Solution values represent the average over all initializationstrategies.
Right:
The initialization strategies are presented relative to the
CollectingDijkstraAny strategy with fixed delta strategy
IncDec . Next, we want to find the best initialization strategy after fixing
IncDec as the deltastrategy. We pair each initialization strategy with
IncDec on the same 1000 instances and sum-marize the results of all pairwise tests after Bonferroni-correction with factor 112 in Table 3.We see that both initialization strategies using Euclidean distances and routing turbineproduction to the nearest free substation, i. e.,
DijkstraAny and
CollectingDijkstraAny ,seem to work best. In particular, these are the only initialization strategies that show somesignificant advantage over other strategies. In Figure 2 (right) we depict ratios of solution val-ues compared to
CollectingDijkstraAny . The minimum ratios are between 0.886 and 0.923for all strategies other than
DijkstraAny (0.974). The maximum ratios range between 1.054and 1.085. For the main part there is hardly any difference between collecting strategiesand their non-collecting counterparts. The figure shows, e. g., that on roughly 22 % of allinstances
CollectingDijkstraAny is better than
BFSAny and
CollectingBFSAny by 0 . CollectingDijkstraAny has a slight but not significant advantage over
DijkstraAny . Wetherefore declare
CollectingDijkstraAny paired with
IncDec as our best variant.The last row in Table 1 shows the running time characteristics of
CollectingDijkstraAny paired with
IncDec . Running times range between tenths of milliseconds and 100 seconds.
We conduct preliminary experiments to determine which
MILP solver we use as a baselinefor our algorithm. To this goal, we randomly choose 35 instances each from benchmarksets N , and N and 70 instances each from benchmark sets N , N , and N . We compareGurobi 8.0.0 and IBM ILOG CPLEX Optimization Studio v12.8 with a running time of oneday per instance and solver using the MILP formulation from Appendix B. Since computingan optimal solution to the
MILP takes too long in almost all instances, we restrict the solvers -0.25-0.20-0.15-0.10-0.050.00 1 hour 12 hours 1 day
Maximum Running Time R e l a t i v e C o s t D i ff e r e n c e N N N N N Figure 3
Comparison solution values found by
MILP solvers CPLEX and Gurobi after differentrunning times. For each instance and running time the abscissa shows a normalized difference insolution values, i. e., (sol
Gurobi − sol CPLEX ) / max(sol Gurobi , sol CPLEX ) . There are twelve instances from N with a value between -0.29 and -0.49 and another three instances from N with a value less than -0.99.Four instances from N are infeasible. to different maximum running times. Each solver uses one thread per instance and nodefiles are written to disk after the solver uses more than 0.5 GB of memory to store node files.Other than that, default values are used.During the experiments, we consider three time stamps: one hour, twelve hours, and oneday. For each solver, instance, and time stamp we record the value of the best incumbent solu-tion and the MIP gap. If a solver terminates with a proven optimal solution after time stamp t ,then the respective values during termination are assigned to all subsequent time stamps.The results of the experiment are depicted in Figure 3 for the quality of the best solutionfound by the respective solver and in Figure 4 for a comparison of MIP gaps. In Figure 3 eachdata point corresponds to an instance and a time stamp. The value on the abscissa stands fora normalized difference in solution values, i. e., (sol Gurobi − sol CPLEX ) / max(sol Gurobi , sol CPLEX ) . Thisyields a value in [ − , MIP gaps computed by CPLEX and Gurobi for each instance andtime stamps. MIP gaps (or relative gaps ) are a standard notion from
Mixed-integer LinearProgramming . The best feasible solution the solver finds yields an upper bound (ub) onthe optimal value. The solver also tries to prove lower bounds (lb). Combining the bestupper and the best lower bound yield the MIP gap ub − lb / ub . This value is in the unit intervaland gives information on how “bad” the solution value can be compared to the (unknown)optimal value. A value of zero shows that the best feasible solution found by the solver isoptimal. Note, however, that a solution might be optimal even though the gap is positive.Evidently, Gurobi performs better across all benchmark sets and time stamps. Whilethere is evidence that the best incumbent solutions computed by Gurobi and CPLEX becomemore similar the longer the experiments run, we also see that Gurobi seems to work betterthan CPLEX the bigger the instances become. We therefore use Gurobi as the MILP solverto compute the baseline to which we compare the negative cycle canceling-based algorithm. . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 15
Gap Gurobi G a p C P L E X Gap Gurobi N N N N N Gap Gurobi
Figure 4
Comparison of gaps between solution values and lower bounds on the optimal value forsolutions computed by CPLEX and Gurobi separated by benchmark sets after different maximumrunning times:
Left: one hour,
Middle: twelve hours,
Right: one day.
We compare our algorithm in its best variant, i. e.,
CollectingDijkstraAny with
IncDec ,with Gurobi on the
MILP formulation in Appendix B. We randomly select 200 instances perbenchmark set from the benchmark sets in [17].In Figure 5 we plot the ratio of the best solution value found by our algorithm to Gurobi’sbest solution at running times of two seconds, one hour, and one day for each benchmark setseparately. These running times represent both interactive and non-time-critical planning.Since our algorithm terminates in under 100 seconds, the comparisons in Figure 5 (middleand right) use the solution our algorithm provides at termination. While discussing theplots, we also discuss an adaptation of the relative gaps ub − lb / ub we introduced in Section 5.2.For each instance, we use the lower bounds from the one-day MILP experiments. For eachinstance, each maximum running time and for both the
MILP and the
NCC algorithm takebest solution value (ub) found at the maximum running time. We refer to the relative gapsas
MILP gap and
NCC gap , respectively, and show them in Figure 6.After two seconds our algorithm outperforms Gurobi on all benchmark sets as it findsbetter solutions on 89 % of all instances with the lowest percentage on benchmark set N .On N the NCC gaps are on average 14 . . MILP gaps of 16 . . N , the NCC gaps are on average 27 . . . MILP gap. The values for N range between those for N and N . The ratiosof solution values range between 0.699 and 1.019 for N , N , and N . On N , which containsthe largest instances, our algorithm computes better solutions on 62 % of the instances. Onsix instances Gurobi does not find a solution. The instances on which Gurobi is better areon average larger than the other instances in N . There are 18 instances on which the ratioof solution values exceeds 1.1 with a maximum of 1.228. On those very large instances,detecting negative cycles takes longer and fewer iterations are performed in two seconds. The NCC gaps spread between 31 . . . MILP gaps areeven worse with a mean value of 48 . . Instances in % R e l a t i v e C o s t s Instances in % N N N N N Instances in %
Figure 5
Comparison of the
NCC algorithm to Gurobi on 200 instances per benchmark set. Theordinate shows the ratio of objective values at various maximum running times of our algorithm toobjective values of Gurobi. Running times:
Left: two seconds,
Middle: one hour,
Right: one day. graphs of N , our algorithm produces solutions that are at least 75 % cheaper than Gurobi’son all but one instance (which has a ratio of 0.411). The gaps are on average at 53 . NCC algorithm and at 92 . N , N , and N . On 25 % ofthe instances from N , on 18 . N , and on one instance from N the solution values are equal. On N and N , our algorithm still yields better solutionson 87 . . . . N ), the ratio exceeds 1.10 with a maximum of 1.165. That means, whilethe NCC algorithm is comparable to Gurobi in solution quality on small instances, it provesbetter on larger wind farms. Furthermore, our algorithm is much faster since it terminatesin under 100 seconds—compared to one hour of maximum running time for Gurobi.After running times of one day (right plot in Figure 5), while our algorithm is at least asgood as Gurobi only on between 25 % ( N ) and 38 . N ) of the instances, it is within 1 %of Gurobi’s solution on 87 . MILP solver is the better choiceif more time is available. Between running times of one hour and one day, the gaps lookvastly the same and there is hardly any difference between
NCC gaps and
MILP gaps. Theyrange between zero and 25 . N , clot around 28 % for N and N and around 34 %for N with seven outliers to the worse by the NCC algorithm.In summary, these experiments show that the
NCC algorithm is a viable option comparedto Gurobi with long running times and that it yields better solutions than the
MILP solver ifonly a short amount of time is given. . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 17
NCC GapMILP Gap t = 2s
NCC GapMILP Gap t = 1h
NCC GapMILP Gap t = 1d
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
NCC GapMILP Gap
Figure 6
Comparison of gaps between solution values and lower bounds on the optimal value forsolutions computed by Gurobi (
MILP ) and the
NCC algorithm separated by benchmark sets ( N inrow 1 through N in row 5) after maximum running times of two seconds, one hour, and one day.Lower bounds are taken from MILP experiments with running times of one day.
Instances in % R e l a t i v e C o s t s N N N N N Instances in % R e l a t i v e C o s t s Figure 7
Comparison of Negative Cycle Canceling algorithm to the Simulated Annealing algorithmon 200 instances per benchmark set. The ordinate represents the ratio of objective values at differentmaximum running times of our algorithm to objective values of the Simulated Annealing algorithm.
Left:
Running time of two seconds.
Right:
Running time of one hour.
We compare our best algorithm variant with the best variant of a Simulated Annealing ( SA )algorithm [17]. We run the SA algorithm on 200 randomly selected instances per benchmarkset (independently selected from other experiments). We compare the best solutions foundafter two seconds and one hour (Figure 7).After two seconds, the NCC algorithm performs at least as good as the SA algorithm onall instances from N and on 74 . . N and N , respectively. The minimumratios are 0.381 for N , 0.911 for N , and 0.875 for N with one instance in N where the SA algorithm does not find a solution. The maximum ratio on those benchmark sets is atmost 1.034. On the larger instances of N and N , our algorithm presumably cannot performsufficient iterations, as the SA algorithm is better on 71 % of those instances. Yet, the SA algorithm does not find feasible solutions on 38 . N . The ratios havea wide spread: from 0.203 to 1.261 for N and from 0.838 to 1.480 for N (save for theinstances without a solution from the SA algorithm).After one hour, the SA algorithm provides better solutions than our algorithm on 67 . N and N , respectively. Our algorithm, however, stays within 1 %in solution quality on 84 . N – N . Again, our algorithm seemsto be stuck in local minima. On N and N , our algorithm performs better than the SA algorithm on 86 % and 74 . SA algorithm needs more timeto explore the solution space. The minimum ratios of solution values are as low as 0.716for N and between 0.905 and 0.995 for the other benchmark sets. The maximum ratios areat most 1.057 for all benchmark sets except N (1.159). This supports our findings fromthe MILP experiments that our algorithm is competitive to other approaches to solving
WCP within very short amounts of time. In view of an interactive planning process, it stands outthat the SA algorithm struggles to find solutions quickly in dense graphs. . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 19 Based on recently presented ideas [14] we propose and compare numerous variants of aNegative Cycle Canceling heuristic for the
Wind Farm Cabling Problem . While allvariants run in the order of milliseconds up to 4.5 minutes, they differ significantly in quality.We identify the best variant and use it to compare our heuristic to the
MILP solver Gurobiand a Simulated Annealing algorithm from the literature. With these comparisons we areable to solve several open questions [14]. While the
MILP solver Gurobi has the potentialto find optimal solutions if it runs long enough, our heuristic is able to find solutions ofcomparable quality in only a fraction of the time. Our algorithm beats Gurobi in findinggood solutions in a matter of seconds. We make similar observations when we compareourselves to a Simulated Annealing approach.Moving forward, one may investigate how to improve the solution quality of our heuristic.Visually comparing flows from our algorithm and other solution methods may help to identifywhat kind of more complex circulations improve the solution. It then remains to investigatehow these circulations can be detected. Also, methods for escaping local minima such astemporarily allowing worse solutions could help to improve our algorithm. It also remainsopen whether one can prove any theoretical guarantees on the solution quality or the numberof iterations. Along the same lines, any theoretical insights on why one delta or initializationstrategy works better than another, or on the order in which cycles should be canceled couldhelp improve the
NCC algorithm.In a broader algorithmic view, the heuristic can be easily generalized to minimum-costflow problems with other types of cost functions provided that one searches for integral flows.It would be interesting to see how well the heuristic performs there.
References
4C Offshore Ltd.
Hornsea Project Three Offshore Wind Farm , 2018. , Accessed: 2018-08-15. Ravindra K. Ahuja, Thomas L. Magnanti, and James B. Orlin.
Network flows: theory,algorithms, and applications . Prentice Hall, Upper Saddle River, NJ [u.a.], 1993. Richard Bellman. On a routing problem.
Quarterly of Applied Mathematics , 16:87–90, 1958. doi:10.1090/qam/102435 . Constantin Berzan, Kalyan Veeramachaneni, James McDermott, and Una-May O’Reilly.Algorithms for cable network design on large-scale wind farms. Technical report, MassachusettsInstitute of Technology, 2011. Boris V. Cherkassky and Andrew V. Goldberg. Negative-cycle detection algorithms.
Mathem-atical Programming , 85(2):277–311, Jun 1999. doi:10.1007/s101070050058 . Ouahid Dahmani, Salvy Bourguet, Mohamed Machmoum, Patrick Guerin, Pauline Rhein,and Lionel Josse. Optimization of the connection topology of an offshore wind farm network.
IEEE Systems Journal , 9(4):1519–1528, 2015. doi:10.1109/JSYST.2014.2330064 . Mauricio C. de Souza, Philippe Mahey, and Bernard Gendron. Cycle-based algorithms formulticommodity network flow problems with separable piecewise convex costs.
Networks ,51(2):133–141, 2008. doi:10.1002/net.20208 . Edsger W. Dijkstra. A note on two problems in connexion with graphs.
Numerische Mathematik ,1(1):269–271, Dec 1959. doi:10.1007/BF01386390 . Lester R. Ford, Jr. and Delbert R. Fulkerson.
Flows in Networks . Princeton University Press,Princeton, NJ, USA, 2010. Virginie Gabrel, Arnaud Knippel, and Michel Minoux. Exact solution of multicommoditynetwork optimization problems with general step cost functions.
Operations Research Letters ,25(1):15 – 23, 1999. doi:10.1016/S0167-6377(99)00020-6 . Andrew V. Goldberg and Tomasz Radzik. A heuristic improvement of the Bellman-Fordalgorithm.
Applied Mathematics Letters , 6(3):3 – 6, 1993. doi:10.1016/0893-9659(93)90022-F . Andrew V. Goldberg and Robert E. Tarjan. Finding minimum-cost circulations by cancelingnegative cycles.
Journal of the ACM , 36(4):873–886, October 1989. doi:10.1145/76359.76368 . Donald Goldfarb, Jianxiu Hao, and Sheng-Roan Kai. Shortest path algorithms using dynamicbreadth-first search.
Networks , 21(1):29–50, 1991. doi:10.1002/net.3230210105 . Sascha Gritzbach, Torsten Ueckerdt, Dorothea Wagner, Franziska Wegner, and Matthias Wolf.Towards negative cycle canceling in wind farm cable layout optimization. In
Proceedingsof the 7th DACH+ Conference on Energy Informatics , volume 1 (Suppl 1). Springer, 2018. doi:10.1186/s42162-018-0030-6 . Longkun Guo and Peng Li. On the complexity of detecting k -length negative cost cycles.In Combinatorial Optimization and Applications, COCOA 2017 , volume 10627 of
LectureNotes in Computer Science , pages 240–250. Springer International Publishing, 2017. doi:10.1007/978-3-319-71150-8_21 . Morton Klein. A primal method for minimal cost flows with applications to the assignmentand transportation problems.
Management Science , 14(3):205–220, 1967. doi:10.1287/mnsc.14.3.205 . Sebastian Lehmann, Ignaz Rutter, Dorothea Wagner, and Franziska Wegner. A simulated-annealing-based approach for wind farm cabling. In
Proceedings of the Eighth InternationalConference on Future Energy Systems , e-Energy ’17, pages 203–215, New York, NY, USA,2017. ACM. doi:10.1145/3077839.3077843 . Sara Lumbreras and Andres Ramos. Optimal design of the electrical layout of an offshore windfarm applying decomposition strategies.
IEEE Transactions on Power Systems , 28(2):1434–1441, 2013. doi:10.1109/TPWRS.2012.2204906 . New York State Energy Research and Development Authority.
New York State OffshoreWind Master Plan , 2017. ,Accessed: 2018-08-15. Adam Ouorou and Philippe Mahey. A minimum mean cycle cancelling method for nonlinearmulticommodity flow problems.
European Journal of Operational Research , 121(3):532 – 548,2000. doi:10.1016/S0377-2217(99)00050-8 . Christos H. Papadimitriou. The complexity of the capacitated tree problem.
Networks ,8(3):217–230, 1978. doi:10.1002/net.3230080306 . Tomasz Radzik and Andrew V. Goldberg. Tight bounds on the number of minimum-meancycle cancellations and related results.
Algorithmica , 11(3):226–242, Mar 1994. doi:10.1007/BF01240734 . David J. Sheskin.
Handbook of parametric and nonparametric statistical procedures . A Chapman& Hall Book. CRC Press, Taylor & Francis, Boca Raton [u.a.], 5. ed. edition, 2011. Pedro Santos Valverde, António J. N. A. Sarmento, and Marco Alves. Offshore wind farmlayout optimization – state of the art.
Journal of Ocean and Wind Energy , 1(1):23–29, 2014. WindEurope asbl/vzw.
Wind in power 2017 , 2018. https://windeurope.org/wp-content/uploads/files/about-wind/statistics/WindEurope-Annual-Statistics-2017.pdf ,Accessed: 2018-08-15. Menghua Zhao, Zhe Chen, and Frede Blaabjerg. Optimization of electrical system for a largeDC offshore wind farm by genetic algorithm. In
Proceedings of NORPIE 2004 , pages 1–8,2004. . Gritzbach, T. Ueckerdt, D. Wagner, F. Wegner, and M. Wolf 21
A Binomial Sign Test for Two Dependent Samples
Statistical tests help to find the best strategy variant for our algorithm with regards toavailable initialization strategies (Section 4.3) and delta strategies (Section 4.4). In our case,we use the
Binomial Sign Test for two dependent samples [23, p. 303].We explain this test in a general setting here and specify how we apply the test in moredetail below. Generally speaking, we compare k variants of an algorithm. In our case theseare the different initialization and delta strategies (Sections 4.3 and 4.4). We apply eachvariant to each instance. For every instance m , we denote the total cost of the resulting flowcomputed by variant i on instance m by X ( i ) m .For any ordered pair of two variants ( i, j ) running on a fixed instance m , we calculateits solution difference D = X ( i ) m − X ( j ) m and increment—depending on the sign of D —either D i
1) tests, one for each ordered pair of variants, and always test thenull hypothesis H : θ = 0 . : θ > . θ is theprobability in the underlying hypothesized distribution D i
Recall that we introduced the notion of cable types in Section 1. Let K denote the set of cabletypes. Each cable type k ∈ K has a capacity on the amount of turbine production that can betransmitted through it, which we denote by cap k , as well as a cost per unit length c k for layinga cable of type k . The MILP formulation for