-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathbridge_ranking_preprocessing2.py
123 lines (103 loc) · 3.88 KB
/
bridge_ranking_preprocessing2.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
"""
Created on Thu Jun 04 12:19:41 2015
@author: cedavidyang
"""
__author__ = 'cedavidyang'
import os
import sys
import psycopg2
import numpy as np
import scipy.stats as stats
import pyNBI.bridge as pybridge
import pyNBI.traffic as pytraffic
import pyDUE.generate_graph as g
import pyDUE.ue_solver as ue
from pyNataf.nataf import natafcurve
from pyNataf.robust import semidefinitive
from pyNBI.risk import social_cost, bridge_cost
from cvxopt import matrix, mul
from multiprocessing import Pool, Manager, freeze_support
import itertools
import time
import datetime
import shelve
# load existing metadata
filename = os.path.join(os.path.abspath('./'), 'Data', 'Python', 'metadata.out')
my_shelf = shelve.open(filename, 'r')
for key in my_shelf:
globals()[key]=my_shelf[key]
my_shelf.close()
# update according to new data
# global variables for parallel computing... stupid multiprocessing in Python
# open databases
conn_gis = psycopg2.connect("dbname='gisdatabase' user='postgres' host='localhost' password=''")
cur_gis = conn_gis.cursor()
conn_nbi = psycopg2.connect("dbname='nbi' user='postgres' host='localhost' password=''")
cur_nbi = conn_nbi.cursor()
# retrieve initial condition states of bridges
bridge_db = pytraffic.retrieve_bridge_db(cur_gis, cur_nbi)
onlinks = [[(15,16,1), (16,15,1)], [(24,25,1), (25,24,1)], [(24,25,1), (25,24,1)], [(24,27,1),(27,24,1)],
[(8,22,1), (22,8,1)], [(37,38,1), (38,37,1)], [(30,37,1), (37,30,1)], [(30,37,1), (37,30,1)],
[(30,37,1), (37,30,1)], [(22,30,1), (30,22,1)], [(22,23,1), (23,22,1)], [(22,23,1), (23,22,1)],
[(22,23,1), (23,22,1)], [(33,36,1), (36,33,1)], [(29,33,1), (33,29,1)], [(29,26,1), (26,29,1)],
[(29,26,1), (26,29,1)], [(13,14,1), (14,13,1)], [(12,13,1), (13,12,1)], [(13,14,1), (14,13,1)]]
for bridge,onlink in zip(bridge_db,onlinks):
bridge[-1] = onlink
# get transition matrix
if os.path.isfile('./pmatrix.npy'):
pmatrix = np.load('pmatrix.npy')
else:
pmatrix = pybridge.transition_matrix()
#create graph
theta = matrix([0.0,0.0,0.0,0.15])
delaytype = 'Polynomial'
graph0 = g.LA_county(parameters=theta,delaytype=delaytype, cur_gis = cur_gis)
nlink = len(graph0.links)
conn_gis.close()
conn_nbi.close()
# capacity drop
cap_drop_array = np.ones(np.asarray(bridge_db, dtype=object).shape[0])*0.5
# initial capacity without failed bridges
all_capacity = np.zeros(nlink)
for link, link_indx in graph0.indlinks.iteritems():
all_capacity[link_indx] = graph0.links[link].capacity
# initial delay
res0 = ue.solver_fw(graph0, full=True, x0=res0[0])
delay0 = res0[1][0,0]
length_vector = np.zeros(len(graph0.links.keys()))
for link_key, link_indx in graph0.indlinks.iteritems():
length_vector[link_indx] = graph0.links[link_key].length
distance0 = (res0[0].T * matrix(length_vector))[0,0]
#res_bench = ue.solver(graph0)
## create bookkeeping dict
#bookkeeping = {}
## correlation
corr_length = 8.73
correlation = pybridge.bridge_correlation(bridge_db, corr_length)
norm_cov = semidefinitive(correlation, tol=1e-14, deftol=1e-12)
#norm_cov = None
# with nataf
popt = np.load('nataf_popt.npy')
def nataf(x):
return natafcurve(x,*popt)
nataf=None
if nataf is not None:
norm_cov = nataf(correlation)
try:
tmp = np.linalg.cholesky(norm_cov)
except:
norm_cov = semidefinitive(norm_cov, tol=1e-14, deftol=1e-12)
filename = os.path.join(os.path.abspath('./'), 'Data', 'Python', 'metadata.out')
my_shelf = shelve.open(filename,'n') # 'n' for new
keys = ['bridge_db', 'all_capacity', 'cap_drop_array', 'res0', 'length_vector', 'pmatrix', 'theta',
'delaytype', 'graph0', 'popt', 'nataf', 'norm_cov', 'corr_length', 'delay0', 'distance0']
for key in keys:
try:
my_shelf[key] = globals()[key]
except:
#
# __builtins__, my_shelf, and imported modules can not be shelved.
#
if not key.startswith("_"):
print('ERROR shelving: {0}'.format(key))
my_shelf.close()