3.2 Comparison with MILP model from literature¶

This is used in the paper Flexible cable routing framework for wind farm collection system optimization.

Preamble¶

[1]:
import urllib.request
import lzma
import shutil
import os
from pathlib import Path
import numpy as np
import xarray as xr
[2]:
from optiwindnet.db import database_connection, G_from_routeset, NodeSet, RouteSet
from optiwindnet.plotting import gplot
[3]:
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import make_axes_locatable

%config InlineBackend.figure_formats = ['svg']
plt.rcParams['svg.fonttype'] = 'none'

style = Path('mplstyles/jupyter_dark.mplstyle')
style.is_file() and plt.style.use(style)

Setup¶

[4]:
def fetch_zenodo_file(record, file_name, target_path: str = '.'):
    url = f'https://zenodo.org/records/{record}/files/{file_name}.xz?download=1'
    try:
        with open(os.path.join(target_path, file_name), 'xb') as f_out:
            with urllib.request.urlopen(url) as response:
                with lzma.open(response) as decompressor:
                    shutil.copyfileobj(decompressor, f_out)
    except FileExistsError:
        return
[5]:
record = '19798734'
DB_FILE = 'MILP_model_comparison.v4.sqlite'
[6]:
fetch_zenodo_file(record, DB_FILE)
[7]:
handle_from_name = {
    'BIG Ronne Bank North': 'rbn',
    'Ormonde': 'ormonde',
    'BIG Ronne Bank South': 'rbs',
    'DanTysk': 'dantysk',
    'Horns Rev 1': 'horns',
    'Thanet': 'thanet',
    'West of Duddon Sands': 'sands',
    'Anholt': 'anholt',
    'SynthTess': 'tess',
    'London Array': 'london',
    'SynthTess (3 OSS)': 'tess3',
}

Axes definition¶

[8]:
Ns_q = NodeSet.select()
[9]:
with database_connection(DB_FILE):
    name2handle = dict((ns.name, next(iter(ns.routesets)).handle) for ns in Ns_q)
name2handle
[9]:
{'BIG Ronne Bank North': 'rbn',
 'Ormonde': 'ormonde',
 'BIG Ronne Bank South': 'rbs',
 'DanTysk': 'dantysk',
 'Horns Rev 1 legacy': 'horns',
 'Thanet': 'thanet',
 'West of Duddon Sands': 'sands',
 'Anholt legacy': 'anholt',
 'SynthTess': 'tess',
 'London Array legacy': 'london',
 'SynthTess (3 OSS)': 'tess3'}
[10]:
node_count = {handle: NodeSet.get(NodeSet.name == name).T
              for name, handle in name2handle.items()}
handleAx = list(node_count.keys())
handleAx.sort(key=lambda k: node_count[k])
handleAx
[10]:
['rbn',
 'ormonde',
 'rbs',
 'dantysk',
 'horns',
 'thanet',
 'sands',
 'anholt',
 'tess',
 'london',
 'tess3']
[11]:
capacityLim = (3, 16)
capacityAx = list(range(*capacityLim))
capacityAx
[11]:
[3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15]
[12]:
solverAx = ['global_optimizer', 'MILP']
[13]:
data = 'length', 'graph', 'runtime'
dims= 'handle', 'capacity', 'solver'
shape = len(handleAx), len(capacityAx), len(solverAx)
optima = xr.Dataset(data_vars=dict(length=(dims, np.full(shape, np.nan)),
                                   graph=(dims, np.empty(shape, dtype=object)),
                                   runtime=(dims, np.full(shape, np.nan))),
                    coords=dict(zip(dims, (handleAx, capacityAx, solverAx))))

Fetch data from database¶

[14]:
with database_connection(DB_FILE):
    for method_idx in solverAx:
        for rs in RouteSet.select().where(
            RouteSet.creator.startswith(method_idx)
        ):
            idx = (rs.handle, rs.capacity, method_idx)
            print(idx, rs.length)
            optima.length.loc[idx] = rs.length
            G = G_from_routeset(rs)
            optima.graph.loc[idx].values.__setitem__((), G)
            optima.runtime.loc[idx] = rs.runtime
('rbn', 3, 'global_optimizer') 33298.17836033638
('rbn', 4, 'global_optimizer') 29860.204861551734
('rbn', 5, 'global_optimizer') 28423.402888795787
('rbn', 6, 'global_optimizer') 27716.325905266214
('rbn', 7, 'global_optimizer') 27168.285951533515
('rbn', 8, 'global_optimizer') 26552.257544399687
('rbn', 9, 'global_optimizer') 26552.244636810836
('rbn', 10, 'global_optimizer') 26201.089671194622
('rbn', 11, 'global_optimizer') 26201.086106804836
('rbn', 12, 'global_optimizer') 26156.787801059225
('rbn', 13, 'global_optimizer') 25805.619229896412
('ormonde', 3, 'global_optimizer') 29324.48987050354
('ormonde', 4, 'global_optimizer') 24688.407433332213
('ormonde', 5, 'global_optimizer') 21429.106016125974
('ormonde', 6, 'global_optimizer') 19474.12162862974
('ormonde', 7, 'global_optimizer') 18039.395753991095
('ormonde', 8, 'global_optimizer') 16921.93611519869
('ormonde', 9, 'global_optimizer') 16921.93611519869
('ormonde', 10, 'global_optimizer') 16921.93611519869
('ormonde', 11, 'global_optimizer') 16739.07924499595
('ormonde', 12, 'global_optimizer') 16760.418413486368
('ormonde', 13, 'global_optimizer') 16703.761325471824
('ormonde', 14, 'global_optimizer') 16733.895536985965
('ormonde', 15, 'global_optimizer') 16408.416545612094
('rbs', 3, 'global_optimizer') 112191.83135823737
('rbs', 4, 'global_optimizer') 96763.06977124394
('rbs', 5, 'global_optimizer') 88555.6831332727
('rbs', 6, 'global_optimizer') 83608.47225858844
('rbs', 7, 'global_optimizer') 81373.62178541177
('rbs', 8, 'global_optimizer') 79779.67989422882
('rbs', 9, 'global_optimizer') 78342.94151532424
('rbs', 10, 'global_optimizer') 77816.19093515755
('rbs', 11, 'global_optimizer') 76847.45689962857
('rbs', 12, 'global_optimizer') 76461.16977503322
('rbs', 13, 'global_optimizer') 75934.43060698418
('rbs', 14, 'global_optimizer') 75684.76644574714
('rbs', 15, 'global_optimizer') 75684.78493691
('dantysk', 3, 'global_optimizer') 186965.4607179916
('dantysk', 4, 'global_optimizer') 150107.1555436056
('dantysk', 5, 'global_optimizer') 128883.88192879168
('dantysk', 6, 'global_optimizer') 114424.59372828214
('dantysk', 7, 'global_optimizer') 105055.65476872603
('dantysk', 8, 'global_optimizer') 97072.5662850161
('dantysk', 9, 'global_optimizer') 92475.3988276232
('dantysk', 10, 'global_optimizer') 88952.72298602328
('dantysk', 11, 'global_optimizer') 84861.28342989853
('dantysk', 12, 'global_optimizer') 81738.84210431519
('dantysk', 13, 'global_optimizer') 79813.65633309666
('dantysk', 14, 'global_optimizer') 76909.60581391084
('dantysk', 15, 'global_optimizer') 76516.82497873653
('horns', 3, 'global_optimizer') 114842.47426147909
('horns', 4, 'global_optimizer') 92186.60718333787
('horns', 5, 'global_optimizer') 79828.4366963156
('horns', 6, 'global_optimizer') 70608.30322850394
('horns', 7, 'global_optimizer') 64103.554765415436
('horns', 8, 'global_optimizer') 59003.343585534516
('horns', 9, 'global_optimizer') 55976.13421886173
('horns', 10, 'global_optimizer') 53410.53765921339
('horns', 11, 'global_optimizer') 53048.65219744707
('horns', 12, 'global_optimizer') 50830.62670976542
('horns', 13, 'global_optimizer') 50281.22439480547
('horns', 14, 'global_optimizer') 49151.833684431
('horns', 15, 'global_optimizer') 49028.77324477393
('thanet', 3, 'global_optimizer') 102067.59439020033
('thanet', 4, 'global_optimizer') 83379.53051809293
('thanet', 5, 'global_optimizer') 72952.83343753527
('thanet', 6, 'global_optimizer') 66434.23934955901
('thanet', 7, 'global_optimizer') 61454.790693540715
('thanet', 8, 'global_optimizer') 58072.80219849351
('thanet', 9, 'global_optimizer') 55632.41817377201
('thanet', 10, 'global_optimizer') 53604.9204210903
('thanet', 11, 'global_optimizer') 52935.0469647424
('thanet', 12, 'global_optimizer') 52519.5723244418
('thanet', 13, 'global_optimizer') 52282.26732291174
('thanet', 14, 'global_optimizer') 51626.11768783746
('thanet', 15, 'global_optimizer') 51172.45741670388
('sands', 3, 'global_optimizer') 172952.6188869118
('sands', 4, 'global_optimizer') 141460.61072459244
('sands', 5, 'global_optimizer') 122476.62244061175
('sands', 6, 'global_optimizer') 110983.14147307574
('sands', 7, 'global_optimizer') 103171.00610692591
('sands', 8, 'global_optimizer') 98437.45024262191
('sands', 9, 'global_optimizer') 94451.19707146697
('sands', 10, 'global_optimizer') 91457.25135532257
('sands', 11, 'global_optimizer') 88469.51895972031
('sands', 12, 'global_optimizer') 86681.12424763002
('sands', 13, 'global_optimizer') 85176.76772892375
('sands', 14, 'global_optimizer') 83223.59901516212
('sands', 15, 'global_optimizer') 82324.85351264001
('anholt', 3, 'global_optimizer') 307601.80897840013
('anholt', 4, 'global_optimizer') 243877.9157801405
('anholt', 5, 'global_optimizer') 206265.11081878003
('anholt', 6, 'global_optimizer') 179735.23789478844
('anholt', 7, 'global_optimizer') 162568.48904861844
('anholt', 8, 'global_optimizer') 147645.0805647098
('anholt', 9, 'global_optimizer') 138261.72441922815
('anholt', 10, 'global_optimizer') 130097.51022406039
('anholt', 11, 'global_optimizer') 122180.8378072613
('anholt', 12, 'global_optimizer') 118563.72956328407
('anholt', 13, 'global_optimizer') 113833.91835451077
('anholt', 14, 'global_optimizer') 109626.84326015905
('anholt', 15, 'global_optimizer') 108125.04204062361
('tess', 3, 'global_optimizer') 186886.02338171666
('tess', 4, 'global_optimizer') 156969.81104506835
('tess', 5, 'global_optimizer') 139734.5976706505
('tess', 6, 'global_optimizer') 129458.15141623422
('tess', 7, 'global_optimizer') 125185.3078497028
('tess', 8, 'global_optimizer') 121802.42817771433
('tess', 9, 'global_optimizer') 119283.33979264004
('tess', 10, 'global_optimizer') 118023.72038947423
('tess', 11, 'global_optimizer') 116942.07457303673
('tess', 12, 'global_optimizer') 115860.40631070288
('tess', 13, 'global_optimizer') 115778.73745447413
('tess', 14, 'global_optimizer') 114778.75672959436
('tess', 15, 'global_optimizer') 114251.24384624511
('london', 3, 'global_optimizer') 273692.07420971623
('london', 4, 'global_optimizer') 225450.88061668695
('london', 5, 'global_optimizer') 195552.13374809423
('london', 6, 'global_optimizer') 178137.77587870634
('london', 7, 'global_optimizer') 166836.08976842053
('london', 8, 'global_optimizer') 156600.64395518746
('london', 9, 'global_optimizer') 150561.5961252862
('london', 10, 'global_optimizer') 143869.12719558683
('london', 11, 'global_optimizer') 140641.09868662452
('london', 12, 'global_optimizer') 137570.67701081812
('london', 13, 'global_optimizer') 137730.35823912258
('london', 14, 'global_optimizer') 133542.84353291648
('london', 15, 'global_optimizer') 135257.49448724478
('tess3', 3, 'global_optimizer') 362392.4229015522
('tess3', 4, 'global_optimizer') 306837.1517530203
('tess3', 5, 'global_optimizer') 275707.70952610986
('tess3', 6, 'global_optimizer') 260145.915885964
('tess3', 7, 'global_optimizer') 254498.86815371612
('tess3', 8, 'global_optimizer') 247119.33918442467
('tess3', 9, 'global_optimizer') 244211.01081018447
('tess3', 10, 'global_optimizer') 242302.63827966966
('tess3', 11, 'global_optimizer') 241070.63252210253
('tess3', 12, 'global_optimizer') 240102.00313237103
('tess3', 13, 'global_optimizer') 241101.99428584546
('tess3', 14, 'global_optimizer') 239865.4727318498
('tess3', 15, 'global_optimizer') 239865.3461584033
('rbn', 3, 'MILP') 33266.36371548993
('rbn', 4, 'MILP') 29860.204861551734
('rbn', 5, 'MILP') 28391.588592928623
('rbn', 6, 'MILP') 27551.248664223054
('rbn', 7, 'MILP') 27168.2535875646
('rbn', 8, 'MILP') 26552.257544399687
('rbn', 9, 'MILP') 26552.23142749388
('rbn', 10, 'MILP') 26201.089322215314
('rbn', 11, 'MILP') 26156.76489956474
('rbn', 12, 'MILP') 26156.76420160615
('rbn', 13, 'MILP') 25805.619927854998
('ormonde', 3, 'MILP') 28468.75482702624
('ormonde', 4, 'MILP') 23857.474053100257
('ormonde', 5, 'MILP') 21276.36757136436
('ormonde', 6, 'MILP') 19474.12162862974
('ormonde', 7, 'MILP') 18039.395753991095
('ormonde', 8, 'MILP') 16921.93611519869
('ormonde', 9, 'MILP') 16921.93611519869
('ormonde', 10, 'MILP') 16921.93611519869
('ormonde', 11, 'MILP') 16730.024237975187
('ormonde', 12, 'MILP') 16694.362130485166
('ormonde', 13, 'MILP') 16694.362130485166
('ormonde', 14, 'MILP') 16694.362130485166
('ormonde', 15, 'MILP') 16408.416545612094
('rbs', 3, 'MILP') 110902.33921821232
('rbs', 4, 'MILP') 96489.50490846757
('rbs', 5, 'MILP') 88462.09862949734
('rbs', 6, 'MILP') 83608.44233166556
('rbs', 7, 'MILP') 81132.15038929948
('rbs', 8, 'MILP') 79779.66991858807
('rbs', 9, 'MILP') 78342.92156404235
('rbs', 10, 'MILP') 77576.12125079533
('rbs', 11, 'MILP') 76752.00687240942
('rbs', 12, 'MILP') 76235.9774658932
('rbs', 13, 'MILP') 75934.45757235313
('rbs', 14, 'MILP') 75684.76718464145
('rbs', 15, 'MILP') 75684.76574816529
('dantysk', 3, 'MILP') 183828.57515589587
('dantysk', 4, 'MILP') 148767.30533416467
('dantysk', 5, 'MILP') 127925.05054530728
('dantysk', 6, 'MILP') 113450.59849386492
('dantysk', 7, 'MILP') 104050.23450286816
('dantysk', 8, 'MILP') 96997.71300844806
('dantysk', 9, 'MILP') 92433.19918314437
('dantysk', 10, 'MILP') 88245.91837868161
('dantysk', 11, 'MILP') 84804.05437119018
('dantysk', 12, 'MILP') 81545.88644903792
('dantysk', 13, 'MILP') 79248.97831538455
('dantysk', 14, 'MILP') 76875.24994486028
('dantysk', 15, 'MILP') 76425.01648507122
('horns', 3, 'MILP') 114135.23903799473
('horns', 4, 'MILP') 91430.7470537994
('horns', 5, 'MILP') 79298.7577318621
('horns', 6, 'MILP') 70483.52870416982
('horns', 7, 'MILP') 64067.0520781751
('horns', 8, 'MILP') 58884.66672358922
('horns', 9, 'MILP') 55928.08103681947
('horns', 10, 'MILP') 53322.23018440208
('horns', 11, 'MILP') 52862.54053053543
('horns', 12, 'MILP') 50866.643435902406
('horns', 13, 'MILP') 50283.19009713246
('horns', 14, 'MILP') 49027.6426352129
('horns', 15, 'MILP') 49028.06565600351
('thanet', 3, 'MILP') 99743.37885680527
('thanet', 4, 'MILP') 82626.23761831594
('thanet', 5, 'MILP') 72359.70992142372
('thanet', 6, 'MILP') 65930.14790575013
('thanet', 7, 'MILP') 61504.0323170166
('thanet', 8, 'MILP') 58024.89002746611
('thanet', 9, 'MILP') 55454.90538854051
('thanet', 10, 'MILP') 53493.64385257402
('thanet', 11, 'MILP') 52854.72835292261
('thanet', 12, 'MILP') 52515.405574957396
('thanet', 13, 'MILP') 51862.64712074536
('thanet', 14, 'MILP') 51624.265800852416
('thanet', 15, 'MILP') 51159.83126677292
('sands', 3, 'MILP') 171522.54983954018
('sands', 4, 'MILP') 140796.20242514837
('sands', 5, 'MILP') 121889.14290338813
('sands', 6, 'MILP') 110405.63024608347
('sands', 7, 'MILP') 103085.02144206477
('sands', 8, 'MILP') 97875.32129300473
('sands', 9, 'MILP') 93822.85401656787
('sands', 10, 'MILP') 90722.38586933207
('sands', 11, 'MILP') 88657.34055081209
('sands', 12, 'MILP') 86277.54368427079
('sands', 13, 'MILP') 85003.89230321214
('sands', 14, 'MILP') 83211.52829881277
('sands', 15, 'MILP') 82149.87665908411
('anholt', 3, 'MILP') 300639.1427423917
('anholt', 4, 'MILP') 240130.73344383217
('anholt', 5, 'MILP') 202978.28824163205
('anholt', 6, 'MILP') 178293.27077721423
('anholt', 7, 'MILP') 160444.70000917363
('anholt', 8, 'MILP') 146734.81319304628
('anholt', 9, 'MILP') 139011.23328317754
('anholt', 10, 'MILP') 130019.40149400773
('anholt', 11, 'MILP') 122371.00453172733
('anholt', 12, 'MILP') 117808.84297279986
('anholt', 13, 'MILP') 113332.92557026623
('anholt', 14, 'MILP') 109665.57303570553
('anholt', 15, 'MILP') 108111.31851484787
('tess', 3, 'MILP') 186198.53958616164
('tess', 4, 'MILP') 155809.8463116733
('tess', 5, 'MILP') 139394.3765547251
('tess', 6, 'MILP') 129458.16599251606
('tess', 7, 'MILP') 124770.88397191273
('tess', 8, 'MILP') 121802.47776368656
('tess', 9, 'MILP') 119283.31769569566
('tess', 10, 'MILP') 118023.72807192127
('tess', 11, 'MILP') 116942.0583395688
('tess', 12, 'MILP') 115860.41130221337
('tess', 13, 'MILP') 114778.73263117071
('tess', 14, 'MILP') 114778.7359845226
('tess', 15, 'MILP') 114251.2244999111
('london', 3, 'MILP') 269312.4251166492
('london', 4, 'MILP') 222515.34310332066
('london', 5, 'MILP') 194018.68429151856
('london', 6, 'MILP') 177414.05766684632
('london', 7, 'MILP') 165754.2212189447
('london', 8, 'MILP') 156194.69277115766
('london', 9, 'MILP') 149734.5493588601
('london', 10, 'MILP') 143835.04216111975
('london', 11, 'MILP') 140271.64771892317
('london', 12, 'MILP') 137387.45485098675
('london', 13, 'MILP') 136027.17574108153
('london', 14, 'MILP') 132234.81972404307
('london', 15, 'MILP') 130009.69988158008
('tess3', 3, 'MILP') 358098.38768541877
('tess3', 4, 'MILP') 304387.5213575615
('tess3', 5, 'MILP') 273742.07777123415
('tess3', 6, 'MILP') 258234.37143086415
('tess3', 7, 'MILP') 251784.56995000975
('tess3', 8, 'MILP') 246387.25213710437
('tess3', 9, 'MILP') 243478.94111512793
('tess3', 10, 'MILP') 241570.61751547427
('tess3', 11, 'MILP') 240424.89920448486
('tess3', 12, 'MILP') 240102.00215962718
('tess3', 13, 'MILP') 239456.19478320668
('tess3', 14, 'MILP') 239456.2253206955
('tess3', 15, 'MILP') 239133.26242696456

runtime¶

[15]:
fig, axs = plt.subplots(1, len(solverAx))
fig.suptitle('Time (h) to reach gap (clipped at 3 h)')
for solver, ax in zip(solverAx, axs):
    img = ax.pcolormesh(
        np.clip(optima.runtime.sel(solver=solver).values.T, 0, 10800)/3600,
        )
    ax.tick_params(top=False, labeltop=False, bottom=True, labelbottom=True)
    ax.set_aspect('equal')
    ax.set_xlabel('site')
    ax.set_xticks(range(len(optima.handle)), optima.handle.data, rotation='vertical')
    ax.set_ylabel(r'$\kappa$')
    ax.set_yticks(np.arange(len(optima.capacity)) + 0.5, range(*capacityLim))
    divider = make_axes_locatable(ax)
    # cax = divider.append_axes("top", size="5%", pad=0.05)
    # plt.colorbar(img, cax=cax, orientation='horizontal', location='top')
    # cax.set_title(Title[i])
    cax = divider.append_axes("right", size="5%", pad=0.05)
    plt.colorbar(img, cax=cax, orientation='vertical')
    ax.set_title(solver)
    ax.grid(False)
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_21_0.svg

total length¶

[16]:
reference = optima.length.sel(solver='MILP').values.T
fig, axs = plt.subplots(1, len(solverAx) - 1)
if len(solverAx) == 2:
    axs = [axs]
divnorm=plt.cm.colors.TwoSlopeNorm(vmin=0.94, vcenter=1., vmax=1.06)
for solver, ax in zip(['global_optimizer'], axs):
    img = ax.pcolormesh(
        optima.length.sel(solver=solver).values.T/reference,
        cmap=plt.cm.seismic,
        norm=divnorm,
    )
    ax.tick_params(top=False, labeltop=False, bottom=True, labelbottom=True)
    ax.set_aspect('equal')
    ax.set_xlabel('site')
    ax.set_xticks(np.arange(len(optima.handle)) + 0.5, optima.handle.data, rotation='vertical')
    ax.set_ylabel(r'$\kappa$')
    ax.set_yticks(np.arange(len(optima.capacity)) + 0.5, range(*capacityLim))
    divider = make_axes_locatable(ax)
    # cax = divider.append_axes("top", size="5%", pad=0.05)
    # plt.colorbar(img, cax=cax, orientation='horizontal', location='top')
    # cax.set_title(Title[i])
    cax = divider.append_axes("right", size="5%", pad=0.05)
    plt.colorbar(img, cax=cax, orientation='vertical',
                 boundaries=np.linspace(0.97, 1.03, 13),
                )
    ax.set_title(solver)
    ax.grid(False)
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_23_0.svg

Instances where Perez-Rua 2020 is shorter than our B&C result¶

[17]:
solverAx
[17]:
['global_optimizer', 'MILP']
[18]:
juru_is_better = ((optima.length[:,:, 0]/optima.length[:,:, 1] - 1) < -5e-4)
juru_is_better.to_pandas()
[18]:
capacity 3 4 5 6 7 8 9 10 11 12 13 14 15
handle
rbn False False False False False False False False False False False False False
ormonde False False False False False False False False False False False False False
rbs False False False False False False False False False False False False False
dantysk False False False False False False False False False False False False False
horns False False False False False False False False False True False False False
thanet False False False False True False False False False False False False False
sands False False False False False False False False True False False False False
anholt False False False False False False True False True False False False False
tess False False False False False False False False False False False False False
london False False False False False False False False False False False False False
tess3 False False False False False False False False False False False False False
[19]:
juru_is_better.sum().item()
[19]:
5
[20]:
merit = optima.length[:,:, 1]/optima.length[:,:, 0] - 1
merit.min().item(), merit.max().item()
[20]:
(-0.038798549578038966, 0.0054209425428310976)
[21]:
juru_idx = np.argwhere(juru_is_better.data)
juru_idx
[21]:
array([[4, 9],
       [5, 4],
       [6, 8],
       [7, 6],
       [7, 8]])
[22]:
for handle_idx, κ_idx in juru_idx:
    length_juru = optima.length[handle_idx, κ_idx, 0].item()
    solver_idx = 1
    length_contestant = optima.length[handle_idx, κ_idx, solver_idx].item()
    G = optima.graph[handle_idx, κ_idx, 0].item()
    for _, _, r in G.edges(data="rogue"):
        if r is not None:
            print(handle_idx, κ_idx)
    num_rogues = sum(1 for _, _, kind in G.edges(data='kind') if kind == 'rogue')
    print(f'{length_juru:.0f} -> {length_contestant:.0f} '
          f'({solverAx[solver_idx]}) '
          f'{100*(length_contestant/length_juru - 1):+.2f}% '
          f'{handleAx[handle_idx]}, κ = {κ_idx + 3}, '
          f'num_rogues = {num_rogues}')
50831 -> 50867 (MILP) +0.07% horns, κ = 12, num_rogues = 0
61455 -> 61504 (MILP) +0.08% thanet, κ = 7, num_rogues = 0
88470 -> 88657 (MILP) +0.21% sands, κ = 11, num_rogues = 0
138262 -> 139011 (MILP) +0.54% anholt, κ = 9, num_rogues = 2
122181 -> 122371 (MILP) +0.16% anholt, κ = 11, num_rogues = 2
[23]:
num_instances = (optima.length[:,:, 0].size - optima.length[:,:, 0].isnull().sum()).item()
num_instances
[23]:
141
[24]:
merit.median().item()
[24]:
-0.0024347664901663846

Create quality comparison plots¶

[25]:
with plt.style.context('pdf_1col'):
    fig, ax = plt.subplots(facecolor='w')
    divnorm=plt.cm.colors.TwoSlopeNorm(vmin=-4, vcenter=0., vmax=4)
    img = ax.pcolormesh(100*merit, cmap=plt.cm.seismic, norm=divnorm)
    ax.tick_params(top=False, labeltop=False, bottom=True, labelbottom=True)
    ax.set_ylabel('Site')
    ax.set_yticks(np.arange(len(optima.handle)) + 0.5, optima.handle.data)
    ax.set_xlabel('Capacity')
    ax.set_xticks(np.arange(len(optima.capacity)) + 0.5, range(*capacityLim))
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("right", size="5%", pad=0.05)
    plt.colorbar(img, cax=cax, orientation='vertical',
                 boundaries=np.linspace(-4, 1, 11),
                )
    cax.set_ylabel('Length change (%)')
    ax.scatter(*(juru_idx.T[::-1] + 0.5), marker='s',
               facecolor='none', color='darkorange', s=140)
    ax.set_aspect('equal')
    ax.grid(False)
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_34_0.svg
[26]:
with plt.style.context('pdf_1col'):
    fig, ax = plt.subplots(facecolor='w')
    hist = ax.hist((100*merit.data).flat, orientation='horizontal', bins=np.linspace(-4, 1, 69),
                   edgecolor='none');
    median = 100*merit.median().item()
    ax.axhline(median, *ax.get_xlim(), color='k', ls=':', lw=0.5,
               label=f'median: {median:.2f}%')
    ax.legend(loc='lower right')
    ax.set_ylabel('Length change (%)')
    ax.set_xlabel('Number of instances')
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_35_0.svg
[27]:
with plt.style.context('pdf_2col'):
    fig, (ax, ax2) = plt.subplots(1, 2, facecolor='w', gridspec_kw={'wspace':  0.1})
    divnorm = plt.cm.colors.TwoSlopeNorm(vmin=-4, vcenter=0., vmax=4)
    img = ax.pcolormesh(100*merit, cmap=plt.cm.seismic, norm=divnorm)
    ax.tick_params(top=False, labeltop=False, bottom=True, labelbottom=True)
    ax.set_ylabel('Site')
    ax.set_yticks(np.arange(len(optima.handle)) + 0.5, optima.handle.data)
    ax.set_xlabel('Capacity')
    ax.set_xticks(np.arange(len(optima.capacity)) + 0.5, range(*capacityLim))
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("right", size="5%", pad=0.05)
    plt.colorbar(img, cax=cax, orientation='vertical',
                 boundaries=np.linspace(-4, 1, 11))
    ax.scatter(*(juru_idx.T[::-1] + 0.5), marker='s',
               facecolor='none', color='darkorange', s=140)
    ax.set_aspect('equal')
    ax.grid(False)

    hist = ax2.hist((100*merit.data).flat, orientation='horizontal', bins=np.linspace(-4, 1, 69),
                   edgecolor='none');
    median = 100*merit.median().item()
    ax2.axhline(median, *ax.get_xlim(), color='k', ls=':', lw=0.5,
                label=f'median: {median:.2f}%')
    ax2.legend(loc='lower right', fontsize='small')
    ax2.set_ylabel('Length difference (%)')
    ax2.set_xlabel('Number of instances')
    ax2.set_ylim(-4, 1)
    fig.savefig('quality_comparison.pdf')
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_36_0.svg

Plot the instance with the highest positive difference¶

[28]:
with plt.style.context('pdf_1col'):
    G = optima.graph.loc['anholt', 9, 'global_optimizer'].item()
    ax = gplot(G, infobox=False, legend=False, dark=False, figsize=(4.2, 50))
    ax.figure.canvas.draw()
    bbox = ax.get_window_extent().transformed(ax.figure.dpi_scale_trans.inverted())
    print(f'{bbox.width:.2f}, {bbox.height:.2f}')
4.18, 1.95
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_38_1.svg
[29]:
ax.figure.savefig('anholt_9_ref.pdf', pad_inches=0., transparent=True)
[30]:
with plt.style.context('pdf_1col'):
    G = optima.graph.loc['anholt', 9, 'MILP'].item()
    ax = gplot(G, infobox=False, legend=False, dark=False, figsize=(4.2, 2.2))
    ax.figure.canvas.draw()
    bbox = ax.get_window_extent().transformed(ax.figure.dpi_scale_trans.inverted())
    print(f'{bbox.width:.2f}, {bbox.height:.2f}')
4.18, 1.95
../_images/notebooks_42-Paper_3.2_MILP_model_comparison_40_1.svg
[31]:
ax.figure.savefig('anholt_9_our.pdf', pad_inches=0., transparent=True)