A Large and Complex Wind Farm¶

This notebook demonstrates how to apply OptiWindNet to optimize electrical network (cable routing) for a large-scale wind farm with 122 turbines and 2 substations, set in a concave area containing an internal obstacle. We will use EWRouter for warm-starting, then run MILPRouter with the COIN-OR CBC solver to find a high-quality branched cable network. The process covers loading real-world location data, generating and visualizing warm-start solutions, configuring solver and model options, and plotting the final optimized network.

Note: Any additional geometric complexity will have negligible impact on solving time, as they do not increase the number of variables in the MILP optimization model.

This notebook uses COIN-OR CBC as solver, but examples for other solvers (gurobi, cplex, ortools, scip, highs) can be found in the MILPRouter notebooks. Gurobi and IBM ILOG CPLEX are comercial solvers (academic license available). Google’s OR-Tools, COIN-OR CBC, SCIP and HiGHS are open-source software.

[1]:
from optiwindnet.api import WindFarmNetwork, EWRouter, MILPRouter, ModelOptions, load_repository
[2]:
# Display figures as SVG in Jupyter notebooks
%config InlineBackend.figure_formats = ['svg']

Load input data¶

Note: the load_repository() functionality of OptiWindNet is used to load a prebuilt Networkx.Graph of the avaible locations. For more details on this functionality look into the notebook about Load repositories containing location data.

[3]:
locations = load_repository()

The example location from https://onlinelibrary.wiley.com/doi/abs/10.1049/rpg2.12593 (Taylor_2023) provides suitable specifications for our purpose; a large wind farm with multiple substations and complex geometry.

[4]:
wfn = WindFarmNetwork(L=locations.taylor_2023, cables=[(5, 1800.0), (8, 2200.0)])

Note: To view output plots, trust this notebook first.

[5]:
wfn
[5]:
../_images/notebooks_a09_large_windfarm_complex_geometry_9_0.svg

There are two options for generating an initial solution to warm-start the MILPRouter:

  • heuristic: EWRouter()

  • meta-heuristic: HGSRouter()

EWRouter has the downside of not limiting the number of feeders, thus the model to be warm-started must also have an unlimited number of feeders. This is usually not a problem, as most cable route sets produced by the MILPRouter uses the minimum number of feeders or only one more than the minimum.

Generate the warm-start solution¶

To generate a warm-start solution, simply run wfn.optimize() with the desired solver e.g.:

wfn.optimize(router=warmstart_router)

In this example, we use EWRouter to warm-start the optimization model. Alternatively, HGSRouter() can be used for the same purpose.

The resulting solution is stored within the wfn object. The next time wfn.optimize() is called, the stored solution (S which contains selected-links) will be reused as a warm start; provided that the solver in MILPRouter supports warm-starting.

[6]:
warmstart_router = EWRouter()
res_warmstart= wfn.optimize(router=warmstart_router)
[7]:
wfn.length()
[7]:
113112.15297348917

Visualizing the warm-start solution¶

S can be visualized using wfn.plot_selected_links().

Check the plotting tutorial for details on the available plots.

[8]:
wfn
[8]:
../_images/notebooks_a09_large_windfarm_complex_geometry_17_0.svg

MILP Router¶

Create an instance of MILPRouter and choose COIN-OR Branch and Cut (CBC) as solver.

When initializing a MILPRouter three arguments are required to be given:

  • solver name

  • time limit

  • gap

Solver options, including those set by OptiWindNet as well as additional configurable parameters, can be modified by creating a dictionary and passing it to the router. The same approach applies to model options. For a detailed explanation about model_options and solver_options of MILP routers see the notebook about ModelOptions vs SolverOptions.

[9]:
solver_options=dict( # these are in addition to the default solver.options
        # if repeatable results are desired, set the seed
        RandomCbcSeed=1234567,
        # CBC works better if the number threads is set to the number of physical cores
        threads=6,
    )

model_options = ModelOptions(
                        topology='branched',
                        feeder_limit='unlimited',
                        feeder_route='segmented',
                    )

milp_router = MILPRouter(solver_name='cbc', time_limit=30, mip_gap=0.005, solver_options=solver_options, model_options=model_options)

Run optimizatin with MILP router¶

[10]:
# Set logging level to INFO to display messages about solver options
import logging
logging.basicConfig(level=logging.INFO)

# run optimization with the milp router (Set verbose=True to display messages about solver progress)
res = wfn.optimize(router=milp_router, verbose=True)
Using warm start: the model is initialized with the provided solution S.

INFO:optiwindnet.MILP.pyomo:>>> cbc solver options <<<
Dins: 'on'
RandomCbcSeed: 1234567
Rens: 'on'
Rins: 'on'
VndVariableNeighborhoodSearch: 'on'
cliqueCuts: 'off'
flowCoverCuts: 'on'
gomoryCuts: 'on'
knapsackCuts: 'off'
liftAndProjectCuts: 'off'
mixedIntegerRoundingCuts: 'on'
nodeStrategy: 'downFewest'
pivotAndComplement: 'off'
probingCuts: 'off'
proximitySearch: 'off'
ratioGap: 0.005
residualCapacityCuts: 'off'
seconds: 30
threads: 6
timeMode: 'elapsed'
twoMirCuts: 'off'
zeroHalfCuts: 'off'

Welcome to the CBC MILP Solver
Version: 2.10.8
Build Date: Jan  1 1970

command line - C:\my_program_files\cbc\bin\cbc.exe -threads 6 -timeMode elapsed -nodeStrategy downFewest -Dins on -VndVariableNeighborhoodSearch on -Rens on -Rins on -pivotAndComplement off -proximitySearch off -gomoryCuts on -mixedIntegerRoundingCuts on -flowCoverCuts on -cliqueCuts off -twoMirCuts off -knapsackCuts off -probingCuts off -zeroHalfCuts off -liftAndProjectCuts off -residualCapacityCuts off -RandomCbcSeed 1234567 -seconds 30 -ratioGap 0.005 -printingOptions all -import C:\Users\amia\AppData\Local\Temp\tmpallszsg2.pyomo.lp -mipstart \Users\amia\AppData\Local\Temp\tmpkwuu819e.cbc.soln -stat=1 -solve -solu C:\Users\amia\AppData\Local\Temp\tmpallszsg2.pyomo.soln (default strategy 1)
threads was changed from 0 to 6
Option for timeMode changed from cpu to elapsed
Option for nodeStrategy changed from fewest to downfewest
Option for Dins changed from off to on
Option for VndVariableNeighborhoodSearch changed from off to on
Option for Rens changed from off to on
Option for gomoryCuts changed from ifmove to on
Option for mixedIntegerRoundingCuts changed from ifmove to on
Option for flowCoverCuts changed from ifmove to on
Option for cliqueCuts changed from ifmove to off
Option for twoMirCuts changed from root to off
Option for knapsackCuts changed from ifmove to off
Option for probingCuts changed from on to off
Option for zeroHalfCuts changed from ifmove to off
randomCbcSeed was changed from -1 to 1234567
seconds was changed from 1e+100 to 30
ratioGap was changed from 0 to 0.005
Option for printingOptions changed from normal to all
opening mipstart file \Users\amia\AppData\Local\Temp\tmpkwuu819e.cbc.soln.
MIPStart values read for 244 variables.
Presolve 4314 (-2) rows, 2891 (-1) columns and 14858 (-1442) elements
Statistics for presolved model
Original problem has 2892 integers (1446 of which binary)
Presolved problem has 2891 integers (1445 of which binary)
==== 1446 zero objective 840 different
==== absolute objective values 840 different
==== for integers 1446 zero objective 840 different
==== for integers absolute objective values 840 different
===== end objective counts


Problem has 4314 rows, 2891 columns (1445 with objective) and 14858 elements
Column breakdown:
0 of type 0.0->inf, 1446 of type 0.0->up, 0 of type lo->inf,
0 of type lo->up, 0 of type free, 0 of type fixed,
0 of type -inf->0.0, 0 of type -inf->up, 1445 of type 0.0->1.0
Row breakdown:
0 of type E 0.0, 243 of type E 1.0, 0 of type E -1.0,
1 of type E other, 0 of type G 0.0, 0 of type G 1.0,
1 of type G other, 2891 of type L 0.0, 1054 of type L 1.0,
124 of type L other, 0 of type Range 0.0->1.0, 0 of type Range other,
0 of type Free
Continuous objective value is 104281 - 0.06 seconds
Cgl0003I 0 fixed, 0 tightened bounds, 528 strengthened rows, 0 substitutions
Cgl0003I 0 fixed, 0 tightened bounds, 528 strengthened rows, 0 substitutions
Cgl0003I 0 fixed, 0 tightened bounds, 202 strengthened rows, 0 substitutions
Cgl0003I 0 fixed, 0 tightened bounds, 121 strengthened rows, 0 substitutions
Cgl0004I processed model has 3881 rows, 2892 columns (2892 integer (1446 of which binary)) and 15680 elements
Cbc0045I MIPStart provided solution with cost 113017
Cbc0012I Integer solution of 113017.26 found by Reduced search after 0 iterations and 0 nodes (0.28 seconds)
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 255 rows 151 columns
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 255 rows 151 columns
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 255 rows 151 columns
Cbc0038I Full problem 3881 rows 2892 columns, reduced to 405 rows 226 columns
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 292 rows 181 columns
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 292 rows 181 columns
Cbc0038I Full problem 3882 rows 2892 columns, reduced to 292 rows 181 columns
Cbc0031I 88 added rows had average density of 88.568182
Cbc0013I At root node, 88 cuts changed objective from 104281.43 to 104895.91 in 19 passes
Cbc0014I Cut generator 0 (Gomory) - 230 row cuts average 346.8 elements, 0 column cuts (0 active)  in 0.486 seconds - new frequency is 1
Cbc0014I Cut generator 1 (MixedIntegerRounding2) - 2 row cuts average 356.0 elements, 0 column cuts (0 active)  in 0.070 seconds - new frequency is 3
Cbc0014I Cut generator 2 (FlowCover) - 163 row cuts average 16.7 elements, 0 column cuts (0 active)  in 0.030 seconds - new frequency is 1
Cbc0010I After 0 nodes, 1 on tree, 113017.26 best solution, best possible 104895.91 (1.99 seconds)
Cbc0010I After 100 nodes, 60 on tree, 113017.26 best solution, best possible 104952.14 (7.87 seconds)
Cbc0010I After 200 nodes, 109 on tree, 113017.26 best solution, best possible 104952.14 (10.44 seconds)
Cbc0010I After 300 nodes, 162 on tree, 113017.26 best solution, best possible 104952.14 (12.46 seconds)
Cbc0010I After 400 nodes, 213 on tree, 113017.26 best solution, best possible 104952.14 (14.34 seconds)
Cbc0010I After 500 nodes, 264 on tree, 113017.26 best solution, best possible 104952.14 (15.93 seconds)
Cbc0010I After 600 nodes, 313 on tree, 113017.26 best solution, best possible 104952.14 (17.23 seconds)
Cbc0010I After 700 nodes, 360 on tree, 113017.26 best solution, best possible 104952.14 (18.35 seconds)
Cbc0010I After 800 nodes, 375 on tree, 113017.26 best solution, best possible 104952.14 (19.10 seconds)
Cbc0010I After 900 nodes, 400 on tree, 113017.26 best solution, best possible 104952.14 (19.70 seconds)
Cbc0010I After 1000 nodes, 423 on tree, 113017.26 best solution, best possible 104952.14 (20.19 seconds)
Cbc0010I After 1100 nodes, 416 on tree, 113017.26 best solution, best possible 104952.14 (20.50 seconds)
Cbc0010I After 1200 nodes, 410 on tree, 113017.26 best solution, best possible 104952.14 (20.87 seconds)
Cbc0010I After 1300 nodes, 406 on tree, 113017.26 best solution, best possible 104952.14 (21.23 seconds)
Cbc0010I After 1400 nodes, 407 on tree, 113017.26 best solution, best possible 104952.14 (21.71 seconds)
Cbc0010I After 1500 nodes, 391 on tree, 113017.26 best solution, best possible 104952.14 (22.33 seconds)
Cbc0010I After 1600 nodes, 413 on tree, 113017.26 best solution, best possible 104952.14 (22.98 seconds)
Cbc0010I After 1700 nodes, 417 on tree, 113017.26 best solution, best possible 104952.14 (23.43 seconds)
Cbc0010I After 1800 nodes, 428 on tree, 113017.26 best solution, best possible 104952.14 (23.91 seconds)
Cbc0010I After 1900 nodes, 461 on tree, 113017.26 best solution, best possible 104952.14 (24.36 seconds)
Cbc0010I After 2000 nodes, 481 on tree, 113017.26 best solution, best possible 104952.14 (24.80 seconds)
Cbc0010I After 2100 nodes, 484 on tree, 113017.26 best solution, best possible 104952.14 (25.23 seconds)
Cbc0010I After 2200 nodes, 480 on tree, 113017.26 best solution, best possible 104952.14 (25.71 seconds)
Cbc0010I After 2300 nodes, 482 on tree, 113017.26 best solution, best possible 104952.14 (26.12 seconds)
Cbc0010I After 2400 nodes, 476 on tree, 113017.26 best solution, best possible 104952.14 (26.48 seconds)
Cbc0010I After 2500 nodes, 496 on tree, 113017.26 best solution, best possible 104952.14 (26.95 seconds)
Cbc0010I After 2600 nodes, 498 on tree, 113017.26 best solution, best possible 104952.14 (27.27 seconds)
Cbc0010I After 2700 nodes, 517 on tree, 113017.26 best solution, best possible 104952.14 (27.72 seconds)
Cbc0010I After 2800 nodes, 537 on tree, 113017.26 best solution, best possible 104952.14 (28.22 seconds)
Cbc0010I After 2900 nodes, 571 on tree, 113017.26 best solution, best possible 104952.14 (28.77 seconds)
Cbc0010I After 3000 nodes, 590 on tree, 113017.26 best solution, best possible 104952.14 (29.30 seconds)
Cbc0030I Thread 0 used 515 times,  waiting to start 0.6256721,  2947 locks, 0.21976614 locked, 0.018067122 waiting for locks
Cbc0030I Thread 1 used 517 times,  waiting to start 0.42010474,  2975 locks, 0.20623255 locked, 0.027341366 waiting for locks
Cbc0030I Thread 2 used 495 times,  waiting to start 1.1103261,  2838 locks, 0.20609736 locked, 0.03138876 waiting for locks
Cbc0030I Thread 3 used 511 times,  waiting to start 1.4357114,  2921 locks, 0.20849109 locked, 0.020550966 waiting for locks
Cbc0030I Thread 4 used 523 times,  waiting to start 1.6045356,  2968 locks, 0.22195411 locked, 0.031655073 waiting for locks
Cbc0030I Thread 5 used 536 times,  waiting to start 1.752903,  3048 locks, 0.20871186 locked, 0.032920122 waiting for locks
Cbc0030I Main thread 26.558228 waiting for threads,  6230 locks, 0.003030777 locked, 0.025787354 waiting for locks
Cbc0020I Exiting on maximum time
Cbc0005I Partial search - best objective 113017.26 (best possible 104952.14), took 188757 iterations and 3092 nodes (30.16 seconds)
Cbc0032I Strong branching done 17988 times (489750 iterations), fathomed 292 nodes and fixed 891 variables
Cbc0035I Maximum depth 99, 88349 variables fixed on reduced cost
Cuts at root node changed objective from 104281 to 104896
Gomory was tried 1161 times and created 2367 cuts of which 0 were active after adding rounds of cuts (6.767 seconds)
MixedIntegerRounding2 was tried 1581 times and created 922 cuts of which 0 were active after adding rounds of cuts (3.641 seconds)
FlowCover was tried 1161 times and created 5082 cuts of which 0 were active after adding rounds of cuts (1.684 seconds)

Result - Stopped on time limit

Objective value:                113017.25676634
Lower bound:                    104952.135
Gap:                            0.08
Enumerated nodes:               3092
Total iterations:               188757
Time (CPU seconds):             30.25
Time (Wallclock seconds):       30.25

Total time (CPU seconds):       30.69   (Wallclock seconds):       30.69

INFO:optiwindnet.MILP.pyomo:>>> Solution <<<
SolutionInfo(runtime=30.69, bound=104952.14, objective=113017.25676634, relgap=0.07136181674462683, termination='maxTimeLimit')

[11]:
wfn.solution_info()
[11]:
{'runtime': 30.69,
 'bound': 104952.14,
 'objective': 113017.25676634,
 'relgap': 0.07136181674462683,
 'termination': 'maxTimeLimit'}
[12]:
wfn.G.size(weight='length')
[12]:
113112.15297348917
[13]:
wfn.length()
[13]:
113112.15297348917

Terse link of the optimized network is returned as output of wfn.optimize(). If not assigned to a variable (e.g. res), it can be easily accessed via wfn.terse_links().

[14]:
res
[14]:
array([  6,   0,   1,  10,   3,  -2,   5,  14,   7,  10,  11,  20,  -2,
        12,  13,  16,  17,  18,  19,  29,  -2,  14,  31,  15,  25,  26,
        27,  28,  -2,  -2,  -2,  41,  42,  34,  35,  36,  37,  38,  -2,
        -2,  -2,  42,  53,  44,  45,  55,  47,  48,  49,  39,  40,  50,
        64,  54,  66,  67,  68,  69,  59,  49,  -2,  60,  51,  62,  65,
        77,  -1,  68,  81,  70,  83,  59,  60,  61,  73,  63,  75,  78,
        79,  80,  -1,  -1,  -1,  82,  83,  84,  85,  74,  87,  90,  91,
        -1,  -1,  -1,  93,  94,  95,  96,  87, 100, 101,  92,  -1,  94,
       103, 104, 114, 106, 100,  -1, 102, 110, 111, 112, 113, 109, 115,
       121, 117, 118, 116, 120])

Plot the optimized network graph¶

[15]:
wfn
[15]:
../_images/notebooks_a09_large_windfarm_complex_geometry_31_0.svg