experiments.py

from os import path
import string
import networkx as nx
import matplotlib as mpl
import matplotlib.pyplot as plt
mpl.use('Agg') # disable pylab's attempt to open a window
from matplotlib.gridspec import GridSpec
from matplotlib.pylab import figure, savefig

examples = [ (nx.generators.classic.complete_graph(12), 'nxcomplete'),
             (nx.generators.lattice.triangular_lattice_graph(5, 4), 'nxtri'),
             (nx.generators.lattice.grid_2d_graph(5, 4), 'nxsquare'),
             (nx.generators.lattice.hexagonal_lattice_graph(3, 3), 'nxhex'),
             (nx.generators.random_graphs.connected_watts_strogatz_graph(20, 4, 0.15), 'nxws'),
             (nx.generators.random_graphs.random_regular_graph(3, 16), 'nxregular'),
             (nx.generators.random_graphs.gnm_random_graph(20, 35), 'nxerdos'),             
             (nx.generators.classic.circular_ladder_graph(12), 'nxcirc'),
             (nx.generators.classic.ladder_graph(12), 'nxladder'),
             (nx.generators.random_graphs.barabasi_albert_graph(15, 3), 'nxba'),
             (nx.generators.random_graphs.powerlaw_cluster_graph(20, 4, 0.15), 'nxtree'),
             (nx.generators.classic.star_graph(20), 'nxstar'),             
             (nx.generators.classic.barbell_graph(6, 2), 'nxbarbell'),
             (nx.generators.classic.wheel_graph(12), 'nxwheel'),
             (nx.generators.classic.path_graph(12), 'nxpath')]

redraw_examples = True
EPS = True
labels = True
single = True
if redraw_examples:
    pos = 0
    if single: 
        fig = plt.figure(constrained_layout = False, figsize = (50, 30))
        gs = GridSpec(3, 5, figure = fig, wspace = 0.2, hspace = 0.2)
        plt.gca().set_axis_off()
        plt.subplots_adjust(top = 0.95, bottom = 0, right = 1, left = 0.02, 
                            hspace = 0, wspace = 0)
        plt.margins(0,0)
        plt.gca().xaxis.set_major_locator(plt.NullLocator())
        plt.gca().yaxis.set_major_locator(plt.NullLocator())
    for (G, label) in examples:
        if not single:
            plt.figure(figsize = (10, 10))
            ax = plt.gca()
        else:
            row = pos // 5
            col = pos % 5
            ax = fig.add_subplot(gs[row, col])
            
        if labels:
            l = string.ascii_uppercase[pos]
            pos += 1
            plt.title(l, loc = 'left', fontsize = 50)
            # ax.text(0.9, 0.9, l, fontsize = 50)
        nx.draw(G, ax = ax)
        _ = ax.axis('off')
        if not single:
            if EPS:
                savefig(label + '.eps', bbox_inches = 'tight', pad_inches = 0)
            else:
                savefig(label + '.png', width = 1000, bbox_inches = 'tight', pad_inches = 0)
    if single:
        if EPS:
            savefig('examples.eps')                
        else:
            savefig('examples.eps')                
    quit()

models = [ (nx.generators.random_graphs.connected_watts_strogatz_graph, 3),
           (nx.generators.random_graphs.barabasi_albert_graph, -2),
           (nx.generators.classic.barbell_graph, 2),
           (nx.generators.classic.circular_ladder_graph, 1),
           (nx.generators.classic.complete_graph, 1),
           (nx.generators.random_graphs.gnm_random_graph, -3),
           (nx.generators.lattice.hexagonal_lattice_graph, 2),
           (nx.generators.classic.ladder_graph, 1),
           (nx.generators.classic.path_graph, 1),
           (nx.generators.random_graphs.random_regular_graph, -1),
           (nx.generators.classic.star_graph, 1),
           (nx.generators.lattice.grid_2d_graph, 2),
           (nx.generators.lattice.triangular_lattice_graph, 2),
           (nx.generators.random_graphs.powerlaw_cluster_graph, 3),
           (nx.generators.classic.wheel_graph, 1) ]

from math import floor, ceil, sqrt, log

def create(m, n, pr, k = None, h = None):
    G = None
    try:
        if pr == 1: # models that take only the graph order
            G = m(n)
        elif pr == 2: # models that are k by k 
            if k * k % 2 == 0:
                G = m(k, k)
            else:
                G = m(k, k + 1)                    
        elif pr == 3: #  WS
            p = log(1 + (11 - h) / 50)
            G =  m(n, k // 2, p)
        elif pr == -1: # regular random 
            G =  m(h, n)
        elif pr == -2: # BA
            G = m(n, h)
        elif pr == -3: # ER
            G = m(n, 2 * n)
    except Exception as e:
        print('# ERROR', m.__name__, 'failed to create a graph:', e)        
        pass
    if G is None:
        print('# OMITTING', m.__name__, 'as it produced a nil output with', n, pr, k, h)
        return None
    CC = sorted(nx.connected_components(G), key=len, reverse=True)
    if len(G) > len(CC[0]):
        print('# WARNING', m.__name__, 'resulted in a disconnected graph with', n, pr, k, h)        
        return G.subgraph(CC[0])
    else:
        return G

# <insert here the contents of import.lst as generated by list.sh>
# <import.lst>
from componentBasedCharacteristics import normalizedGCC
from componentBasedCharacteristics import splittingNumber
from componentBasedCharacteristics import randomRobustnessIndex
from componentBasedCharacteristics import robustnessMeasure53
from componentBasedCharacteristics import connectivityRobustnessFunction
from componentBasedCharacteristics import kResilienceFactor
from componentBasedCharacteristics import resilienceFactor
from componentBasedCharacteristics import perturbationScore
from componentBasedCharacteristics import robustnessIndex
from componentBasedCharacteristics import robustnessMeasureR
from degreeBasedCharacteristics import degreeEntropy
from degreeBasedCharacteristics import relativeEntropy
from densityBasedCharacteristics import hubDensity
from densityBasedCharacteristics import definition523
from distanceBasedCharacteristics import geographicalDiversity
from distanceBasedCharacteristics import effectiveGeographicalPathDiversity
from distanceBasedCharacteristics import totalGraphGeographicalDiversity
from distanceBasedCharacteristics import compensatedTotalGeographicalGraphDiversity
from distanceBasedCharacteristics import globalFunctionalityLoss
from distanceBasedCharacteristics import temporalEfficiency
from distanceBasedCharacteristics import deltaEfficiency
from distanceBasedCharacteristics import dynamicFragility
from distanceBasedCharacteristics import vulnerability
from flowBasedCharacteristics import electricalNodalRobustness
from flowBasedCharacteristics import relativeAreaIndex
from otherCharacteristics import entropy
from otherCharacteristics import effectiveGraphResistance
from otherCharacteristics import viralConductance
from otherCharacteristics import RCB
from pathBasedCharacteristics import localConnectivity
from pathBasedCharacteristics import globalConnectivity
from pathBasedCharacteristics import pairwiseDisconnectivityIndex
from pathBasedCharacteristics import fragmentation
from pathBasedCharacteristics import kVertexFailureResilience
from pathBasedCharacteristics import vertexResilience
from pathBasedCharacteristics import kEdgeFailureResilience
from pathBasedCharacteristics import edgeResilience
from pathBasedCharacteristics import pathDiversity
from pathBasedCharacteristics import percolatedPath
from pathBasedCharacteristics import percolationCentrality
from pathBasedCharacteristics import treeness
from randomWalkBasedCharacteristics import networkCriticality
from spectralCharacteristics import reconstructabilityCoefficient
from spectralCharacteristics import normalizedSubgraphCentrality
from spectralCharacteristics import generalizedRobustnessIndex
from spectralCharacteristics import redundancyOfAlternativePaths
from spectralCharacteristics import naturalConnectivity
from spectralCharacteristics import subgraphCentrality
from spectralCharacteristics import normalizedLocalNaturalConnectivity
# </import.lst>

# <import2.lst>
from componentBasedCharacteristics import maximumPerturbationScore
# </import2.lst>

import igraph as ig
from time import time
from sys import argv, stderr
from multiprocessing import Process, freeze_support
from math import floor, ceil, sqrt, log

def measure(replica, c, descr, g, g2 = None):
    before, after, value = None, None, None
    if g2 == None:
        try:
            before, value, after = time(), c(g), time()
        except Exception as e:
            print('# ERROR measuring', descr, e, replica)
            return
    else:
        try:
            before, value, after = time(), c(g, g2), time()
        except Exception as e:
            print('# ERROR measuring', descr, e, replica)
            return            
    if value is not None:
        try:        
            if hasattr(value, "__iter__"):
                value = '{:f} (avg)'.format(sum(value) / len(value))
            print('{:s} {:s} {:f}'.format(descr, str(value), after - before), replica)
        except Exception as e:
            print('# ERROR reporting', descr, e, replica)
            return                        

def run(permitted = 1):
    freeze_support()
    # <insert in the following list the contents of char.lst>
    characteristics = [
        normalizedGCC,
        splittingNumber,
        randomRobustnessIndex,
        robustnessMeasure53,
        connectivityRobustnessFunction,
        kResilienceFactor,
        resilienceFactor,
        perturbationScore,
        robustnessIndex,
        robustnessMeasureR,
        degreeEntropy,
        relativeEntropy,
        hubDensity,
        definition523,
        geographicalDiversity,
        effectiveGeographicalPathDiversity,
        totalGraphGeographicalDiversity,
        compensatedTotalGeographicalGraphDiversity,
        globalFunctionalityLoss,
        temporalEfficiency,
        deltaEfficiency,
        dynamicFragility,
        vulnerability,
        electricalNodalRobustness,
        relativeAreaIndex,
        entropy,
        effectiveGraphResistance,
        viralConductance,
        RCB,
        localConnectivity,
        globalConnectivity,
        pairwiseDisconnectivityIndex,
        fragmentation,
        kVertexFailureResilience,
        vertexResilience,
        kEdgeFailureResilience,
        edgeResilience,
        pathDiversity,
        percolatedPath,
        percolationCentrality,
        treeness,
        networkCriticality,
        reconstructabilityCoefficient,
        normalizedSubgraphCentrality,
        generalizedRobustnessIndex,
        redundancyOfAlternativePaths,
        naturalConnectivity,
        subgraphCentrality,
        normalizedLocalNaturalConnectivity
    ] 
    # </char.lst> do not forget to remove the last comma
    # <insert in this list the content of char2.lst>
    comparative = [
        maximumPerturbationScore
    ]
    # </char2.lst> do not forget to remove the last comma

    executed = set()
    filename = f'results_{permitted}sec.txt'
    if path.isfile(filename):
        with open(filename) as data:
            for line in data:
                fields = line.split() # order and model and replica
                keep = fields[:2] + [fields[-1]]
                case = ' '.join(keep)
                executed.add(case)
    for power in range(5, 10):
        n = 2**power
        k = int(floor(sqrt(n)))
        h = int(ceil(sqrt(k)))
        for (m, p) in models:
            for r in range(10 - power // 2):
                d = '{:s} {:d}'.format(m.__name__, n)
                case = f'{d} {r}'
                if case not in executed:
                    G = create(m, n, p, k, h)
                    if G is not None:
                        G = nx.convert_node_labels_to_integers(G)
                        for vertex in G:
                            if 'pos' in G.nodes[vertex]:
                                del G.nodes[vertex]['pos']
                        nx.write_graphml(G, 'G.graphml')
                        g = ig.read('G.graphml', format="graphml")
                        for c in characteristics:
                            proc = Process(target=measure, args=(r, c, d + ' ' + c.__name__, g))
                            proc.start()
                            proc.join(permitted)
                            if proc.is_alive():
                                proc.terminate()
                                proc.join()
    for power in range(4, 9, 2):
        n = 2**power
        k = int(floor(sqrt(n)))
        h = int(ceil(sqrt(k)))
        for (m1, p1) in models:
            for (m2, p2) in models:
                for r in range(10 - power // 2): # less replicas for larger graphs
                    d = '{:s} {:s} {:d}'.format(m1.__name__ , m2.__name__, n)
                    case = f'{d} {r}'
                    if case not in executed:
                        G1 = create(m1, n, p1, k, h)
                        if G1 is not None:
                            G1 = nx.convert_node_labels_to_integers(G1)
                            for vertex in G1:
                                if 'pos' in G1.nodes[vertex]:
                                    del G1.nodes[vertex]['pos']
                            nx.write_graphml(G1, 'G1.graphml')
                            g1 = ig.read('G1.graphml', format="graphml")
                            G2 = create(m2, n, p2, k, h)
                            if G2 is not None:
                                G2 = nx.convert_node_labels_to_integers(G2)
                                for vertex in G2:
                                    if 'pos' in G2.nodes[vertex]:
                                        del G2.nodes[vertex]['pos']
                                nx.write_graphml(G2, 'G2.graphml')
                                g2 = ig.read('G2.graphml', format="graphml")
                                for c in comparative:
                                    p = Process(target=measure, args=(r, c, d + ' ' + c.__name__, g1, g2))
                                    p.start()
                                    p.join(permitted)
                                    if p.is_alive():
                                        p.terminate()
                                        p.join()
                            
if __name__ == '__main__':
    permitted = None
    try: # check for how many seconds until the execution of an individual measurement is terminated
        permitted = int(argv[1])
    except:
        permitted = 1 # default
        stderr.write('using the default one-second timeout')
    run(permitted)