重要代码

"""
主要变量说明：
_C  : 道路设计承载量
_D  : 道路长度
_OD : O/D矩阵
_DEG: 结点度数列表
_B  : 结点介数列表
_Q  : 当前道路承载量
_L  : 道路时间开销
_SL : 两点之间最短时间开销
"""
def dijkstra(self, G):
    """
    dijkstra算法
    """
    # 邻接表
    adj = {}
    # 距离
    dists = {}
    # 前置点
    prevs = {}
    for i in xrange(self.size):
        adj[i] = []
        for j in xrange(self.size):
            if G[i][j] != float('inf'):
                adj[i].append(j) 
    for o in xrange(self.size):
        vis = {x: False for x in xrange(self.size)}
        min_heap = []
        heapq.heappush(min_heap, (0, o))
        dist = {x: float('inf') for x in xrange(self.size)}
        dist[o] = 0
        prev = {x: x for x in xrange(self.size)}
        while len(min_heap) > 0:
            d, i = heapq.heappop(min_heap)
            if vis[i]: continue
            vis[i] = True
            for j in adj[i]:
                if not vis[j] and dist[j] > dist[i] + G[i][j]:
                    # 松弛操作
                    dist[j] = dist[i] + G[i][j]
                    prev[j] = i
                    heapq.heappush(min_heap, (dist[j], j))
        dists[o] = dist
        prevs[o] = prev
    return dists, prevs
def gen_D(self, is_city):
    """
    随机生成路网
    """
    self._D = np.full((self.base_size, self.base_size), np.inf)
    self._DEG = np.zeros(self.base_size)
    for i in xrange(self.base_size):
        self._D[i][i] = 0
        deg = random.randint(0, min(4, self.base_size-1))
        connection_list = random.sample(range(self.base_size), deg)
        while i in connection_list:
            connection_list = random.sample(range(self.base_size), deg)
        for j in connection_list:
            if self._DEG[i] >= 4 or self._DEG[j] >= 4: continue
            if is_city:
                # 100 ~ 500 m
                self._D[i][j] = self._D[j][i] = random.random() * 0.4 + 0.1
            else:
                # 30 ~ 300 m
                self._D[i][j] = self._D[j][i] = random.random() * 0.27 + 0.03 
            self._DEG[i] += 1
            self._DEG[j] += 1
def cal_B(self):
    """
    计算介数
    """
    self._B = np.zeros(self.size)
    dists, prevs = self.dijkstra(self._D)
    for o in xrange(self.size):
        for d in xrange(self.size):
            if o == d or dists[o][d] == float('inf'): continue
            u = d
            while u != o:
                self._B[u] += 1
                u = prevs[o][u]
            self._B[o] += 1
    total = sum(self._B)
    self._B /= total
def cal_sp_fast(self):
    """
    计算任意两点最短路以及相应路径矩阵
    """
    self._PATH = np.full(self._D.shape, -1, 'int32')
    self._SL = np.array(self._L)
    adj = {}
    dists, prevs = self.dijkstra(self._L)
    for o in xrange(self.size):
        for d in xrange(self.size):
            self._SL[o][d] = dists[o][d]
            if o == d or self._SL[o][d] == float('inf'): continue
            u = d
            while u != o:
                self._PATH[prevs[o][u]][d] = u
                u = prevs[o][u]
def iteration(self, eta=0.2, it_no=1, sp_func=None):
    """
    Frank-Wolfe算法迭代
    """
    _NEW_Q = np.zeros(self._Q.shape)
    # 更新最短时间路径的函数(即cal_sp_fast)
    sp_func()
    for i in xrange(self.base_size):
        for j in xrange(self.base_size):
            if i == j: continue
            if not self.is_reachalbe(i, j): continue                
            s = i
            while s != j:
                t = self._PATH[s][j]
                _NEW_Q[s][t] = self._OD[i][j]
                s = t
    _DETAL_Q = _NEW_Q - self._Q
    for i in xrange(self.base_size):
        for j in xrange(self.base_size):
            if i == j or self._D[i][j] == float('inf'): continue
            self._Q[i][j] += eta * _DETAL_Q[i][j] / it_no ** 0.3
            self._L[i][j] = bpr(self._D[i][j] / SPEED, self._Q[i][j], self._C[i][j])
    ev = self.evaluate()
    return ev
def converge_fast(self):
    """
    尝试收敛
    """
    self._L = self._D / SPEED
    pre = self.first_iteration(sp_func=self.cal_sp_fast)
    now = self.iteration(sp_func=self.cal_sp_fast)
    cnt = 1
    while now < pre:
        pre = now
        now = self.iteration(it_no=cnt, sp_func=self.cal_sp_fast)
        cnt += 1
    return cnt
def get_hub(self, num, mode=1):
    """
    获取介数最小，适中，最大的num个结点
    """
    self.cal_B()
    b = zip(self._B, range(len(self._B)))
    sorted_b = [y for x, y in sorted(b)]
    if num > len(sorted_b):
        num = len(sorted_b)
    if mode == 0:
        ret = sorted_b[:num]
    elif mode == 1:
        mid = len(sorted_b) / 2
        l = mid - num/2
        r = mid + num/2
        while r - l > num: r -= 1
        while r - l < num: r += 1
        ret = sorted_b[l:r]
    elif mode == 2:
        ret = sorted_b[-num:]
    return np.array(ret, dtype=int)
def merge(self, other, conn_num=3, mode=0):
    """
    按照不同mode合并城市与小区
    """
    origin_size = self.size
    self._D = self.merge_matrix(self._D, other._D, np.inf)
    self._C = self.merge_matrix(self._C, other._C, 0.0)
    self.size = self._D.shape[0]
    city_v = self.get_hub(conn_num, mode/3)
    area_v = other.get_hub(conn_num, mode%3) + origin_size
    conn_num = min(len(city_v), len(area_v))
    for i in xrange(conn_num):
        self._D[city_v[i]][area_v[i]] = self._D[area_v[i]][city_v[i]] = 0
        self._C[city_v[i]][area_v[i]] = self._C[area_v[i]][city_v[i]] = np.inf