利用Frank-Wolfe求解UE用户均衡模型,以SiouxFalls网络为例(Python)

本资源中利用Frank-Wolfe算法求解了SiouxFalls网络的交通分配结果(UE用户均衡结果,即没有用户可以通过单方面改变出行路径从而降低出行费用)。网络基本信息如txt文件所示,路阻函数采用了经典美国联邦公路局BPR函数,出行需求如ODPairs.txt所以。最终输出结果为输出每条路段上的流量。
资源下载链接链接: link

Python源代码如下:


```python
#第一步:基本类定义(路段,节点,OD对)
import os
class node():#节点类,包含每个节点的坐标信息
    def __init__(self):
        self.ID=None
        self.pos_x=None
        self.pos_y=None
        self.origin=-1
        self.innode=[]#进入节点路段集合
        self.outnode=[]#离开节点路段集合
class link(): #路段类,包含路段的起讫点、自由流行驶时间、容量、BRP路阻函数的参数
    def __init__(self):
        self.ID=None
        self.O_link=None#起点
        self.D_link=None#终点
        self.FFT=None
        self.Traveltime=None
        self.Capacity=None
        self.B=0.15#SiouxFalls网络中路网B与Power系数值均相等,方便期间,统一设置一个初始值
        self.Power=4
class OD():#OD类,记录每一个OD需求信息,包含起点、重点、流量需求
    def __init__(self):
        self.ID=None
        self.origin_node=None
        self.Destination=[]
        self.odlinks_demand=[]
#第二步:定义路网核心函数(路网信息导入、UE分配计算等)
class network():#网络类,核心函数
    def __init__(self):
        self.Nodes=[]#网络包含节点集合
        self.Links=[]#路段集合
        self.Origins=[]#起点集合,即流量发点
        self.num_Nodes=0#节点个数
        self.num_Links=0#路段个数
        self.num_Origins=0#起点个数
        self.Linkflows=[]#路段流量集合
        self.Linktimes=[]#路段行驶时间集合
        self.max_err=0.001#UE最大误差
        self.err=1#初始UE误差
        self.Djpathcost=[]#最短路阻抗集合(Dijkstra算法求得)
        self.Djpath=[]#最短路集合
    def read_nodes(self,path):#读取节点信息
        f=open(path)
        file=f.readlines()
        for i in file:
            row=i.strip('\n').split('\t')
            newnode=node()
            newnode.ID=int(row[0])
            newnode.pos_x=float(row[1])
            newnode.pos_y=float(row[2])
            self.Nodes.append(newnode)
        f.close()
    def read_link(self,path):#读取路段信息
        f=open(path)
        file=f.readlines() 
        linkid=1
        for i in file:
            row=i.strip('\n').split('\t')
            newlink=link()
            newlink.ID=linkid
            newlink.O_link=self.Nodes[int(row[0])-1]
            newlink.D_link=self.Nodes[int(row[1])-1]
            newlink.FFT=float(row[2])
            newlink.Capacity=float(row[3])
            newlink.O_link.outnode.append(newlink.ID)
            newlink.D_link.innode.append(newlink.ID)
            
            
            #self.Nodes[int(row[0])-1].innode.append(newlink.ID)
            linkid+=1
            self.Links.append(newlink) 
        f.close()
    def read_od(self,path):#读取OD信息
        f=open(path)
        file=f.readlines() 
        linkid=1
        for i in file:
            row=i.strip('\n').split('\t')
            newnode=self.Nodes[int(row[0])-1]
            if newnode.origin==-1:
                neworigin=OD()
                self.num_Origins+=1
                neworigin.ID=self.num_Origins
                neworigin.origin_node=newnode
                self.Nodes[int(row[0])-1].origin=neworigin.ID
                self.Origins.append(neworigin)
            else:
                neworigin=self.Origins[newnode.origin-1]
                
            self.Origins[newnode.origin-1].Destination.append(self.Nodes[int(row[1])-1].ID)
            self.Origins[newnode.origin-1].odlinks_demand.append(float(row[2]))
                
        f.close()        
    def Dijkstra(self,start,end):#Dijkstra求解最短路(也叫标号法)
        startpos=0
        endpos=1
        path=[]
        checkpath=[None for i in range(len(self.Nodes))]
        boolcheckpath=[]
        self.Djpathcost=[]
        self.Djpath=[None for i in range(len(self.Nodes))]
        bscanStatus=[None for i in range(len(self.Nodes))]
        for i in range(len(self.Nodes)):
            self.Djpath.append(-1)
            self.Djpathcost.append(9999999)#标号法初始路阻最大
            boolcheckpath.append(False)
            
        self.Djpathcost[start-1]=0
        checkpath[0]=start-1
        while startpos!=endpos:
            if startpos>=len(self.Nodes):
                startpos=0
            i=checkpath[startpos]
            startpos+=1
            newnode=self.Nodes[i]
            for j in range(len(newnode.outnode)):
                newlink=self.Links[newnode.outnode[j]-1]
                k=newlink.D_link.ID
                tt=newlink.Traveltime
                if self.Djpathcost[k-1]>self.Djpathcost[i]+tt:
                    self.Djpathcost[k-1]=self.Djpathcost[i]+tt
                    self.Djpath[k-1]=i
                    if endpos>=len(self.Nodes):
                        endpos=0
                    checkpath[endpos]=k-1
                    endpos+=1
                    bscanStatus[k-1]=True
        return self.Djpathcost[end-1]
    def Dijkstra_path(self,start,end):#记录最短路径,与上述函数基本相同,输出结果不同,为了后面调用方便,写了两次
        startpos=0
        endpos=1
        path=[]
        checkpath=[None for i in range(len(self.Nodes))]
        boolcheckpath=[]
        self.Djpathcost=[]
        self.Djpath=[None for i in range(len(self.Nodes))]
        bscanStatus=[None for i in range(len(self.Nodes))]
        for i in range(len(self.Nodes)):
            self.Djpath.append(-1)
            self.Djpathcost.append(9999999)
            boolcheckpath.append(False)         
        self.Djpathcost[start-1]=0
        checkpath[0]=start-1
        while startpos!=endpos:
            if startpos>=len(self.Nodes):
                startpos=0
            i=checkpath[startpos]
            startpos+=1
            newnode=self.Nodes[i]
            for j in range(len(newnode.outnode)):
                newlink=self.Links[newnode.outnode[j]-1]
                k=newlink.D_link.ID
                tt=newlink.Traveltime
                if self.Djpathcost[k-1]>self.Djpathcost[i]+tt:
                    self.Djpathcost[k-1]=self.Djpathcost[i]+tt
                    self.Djpath[k-1]=i
                    if endpos>=len(self.Nodes):
                        endpos=0
                    checkpath[endpos]=k-1
                    endpos+=1
                    bscanStatus[k-1]=True
        Djpathlink=[]
        point_out=end-1
        while True:
            i=0
            point_in=self.Djpath[point_out]
            for j in range(len(self.Links)):
                newlink=self.Links[j]
                if point_in==newlink.O_link.ID-1 and point_out==newlink.D_link.ID-1:
                    Djpathlink.insert(0,newlink.ID)
                    point_out=point_in
            i+=1
            if point_in==start-1:
                break        
        return Djpathlink
  
    def all_none(self):#全有全无分配函数
        all_none_linkflow=[0 for i in range(len(self.Links))]
        for i in range(len(self.Links)):           
            self.Links[i].Traveltime=self.Links[i].FFT*(1+self.Links[i].B*(self.Linkflows[i]/self.Links[i].Capacity)**self.Links[i].Power)#更新路段行驶时间
            all_none_linkflow[i]=0
        for i in range(len(self.Origins)):       
            o_node=self.Origins[i].origin_node.ID
            for j in range(len(self.Origins[i].Destination)):
                d_node=self.Origins[i].Destination[j]
                demand=self.Origins[i].odlinks_demand[j]
                Djpathlink=self.Dijkstra_path(o_node,d_node)#找最短路
                for index in Djpathlink:
                    all_none_linkflow[index-1]+=demand#将流量加载到最短路上                  
        return all_none_linkflow
    
    def getUEerr(self):#计算UE误差
        sum1=0
        for i in range(len(self.Links)):
            
            self.Links[i].Traveltime=self.Links[i].FFT*(1+self.Links[i].B*(self.Linkflows[i]/self.Links[i].Capacity)**self.Links[i].Power)
            sum1+=self.Links[i].Traveltime*self.Linkflows[i]#计算流量与行驶时间的乘积(UE公式中的积分项)
        
        sum2=0
        for i in range(len(self.Origins)):
            for j in range(len(self.Origins[i].Destination)):
                demand=self.Origins[i].odlinks_demand[j]
                cost=self.Dijkstra(self.Origins[i].origin_node.ID,self.Origins[i].Destination[j])
                sum2+=demand*cost#计算需求与行驶时间的乘积
           
        return 1-sum2/sum1
    def Optfunction(self,Descent,Lamuda):#计算最优下行方向函数,返回值越大,下降方向越优
        Sum=0
        for i in range(len(self.Links)):
            x=self.Linkflows[i]+Lamuda*Descent[i]
            Sum+=Descent[i]*(self.Links[i].FFT*(1+self.Links[i].B*(x/self.Links[i].Capacity)**self.Links[i].Power))
        return Sum
 
    def Frank_Wolfe(self):#Frank-Wolfe主函数
        
        iter=0#迭代次数
        self.Linkflows=[0 for i in range(len(self.Links))]
        self.Linkflows=self.all_none()
        while self.err>self.max_err:
            oldlinkflow=self.Linkflows
            newlinkflow=self.all_none()
            Descent=[]
            for i in range(len(self.Links)):
                Descent.append(newlinkflow[i]-self.Linkflows[i])
            Lamuda=0
            left=0#二分法求最优步长
            right=1
            mid=0
            f_left=self.Optfunction(Descent,left)      
            f_right=self.Optfunction(Descent,right)      
            f_mid=0
            if f_left*f_right>0:
                if abs(f_left)>abs(f_right):
                    Lamuda=right
                else:
                    Lamuda=left
                    
            else:              
                while right-left>self.max_err:
                    mid=(left+right)/2
                    f_left=self.Optfunction(Descent,left)
                    f_right=self.Optfunction(Descent,right)
                    f_mid=self.Optfunction(Descent, mid)
                    if f_left*f_mid>0:
                        left=mid
                    else: right=mid
                Lamuda=(left+right)/2
    
            for i in range(len(self.Links)):#更新路段流量
                self.Linkflows[i]+=Lamuda*Descent[i] 
            iter+=1
            self.err=self.getUEerr()        

        #print(self.Linkflows)
        
net=network()
#读取路网文件
net.read_nodes('Node.txt')
net.read_link('Link.txt')
net.read_od('ODPairs.txt')
net.Frank_Wolfe()#F-W主函数
print(net.Linkflows)#输出每条路段上的流量
           
        
        
    
        
        
        
        
    
    
  • 7
    点赞
  • 51
    收藏
    觉得还不错? 一键收藏
  • 3
    评论
求解交通流量分配模型的有效方法 #include "stdafx.h" #include #include #include "os.h" #include "my_types.h" #include "md_alloc.h" #include "my_util.h" #include "message.h" #include "tui.h" #include "meta.h" #include "link_functions.h" #include "od_table.h" #include "od_flow.h" #include "mincostroutes.h" #include "ls_bisect.h" #include "fw_status.h" extern int no_zones, no_nodes, no_links; /* Gloabal variables */ my_float **ODflow, TotalODflow; /* Local Declarations */ /* void FW(void); Should there be a function for fw, or should it be included in main? */ static void Init(char *tuiFileName); static void Close(char *tuiFileName); static void InitODflow(void); static void CloseODflow(void); /* End of local declarations. */ void main(int argc, char **argv ) { my_float *MainVolume, *SubVolume, *SDVolume, Lambda; int **MinPathPredLink; struct fw_status_struct fw_status; char *tuiFileName; StatusMessage("General", "Ready, set, go..."); switch(argc){ case 2: tuiFileName=argv[1]; break; case 1: tuiFileName="control.tui"; break; default: ExitMessage("Wrong number of command line arguments (%d). \n" "Syntax: fw .", argc-1); } Init(tuiFileName); MainVolume = (my_float*)Alloc_1D(no_links, sizeof(my_float) ); SDVolume = SubVolume = (my_float*)Alloc_1D(no_links, sizeof(my_float) ); /* Compute search direction and sub-volume in the same place. */ MinPathPredLink = (int**)Alloc_2D(no_zones,no_nodes, sizeof(int)); InitFWstatus(&fw_status); FindMinCostRoutes (MinPathPredLink, NULL); Assign (ODflow,MinPathPredLink,MainVolume); FirstFWstatus(MainVolume, &fw_status); UpdateLinkCost(MainVolume); for ( fw_status.Iteration = 1; ContinueFW(fw_status); fw_status.Iteration++) { FindMinCostRoutes (MinPathPredLink, NULL); Assign (ODflow,MinPathPredLink,SubVolume); VolumeDifference( SubVolume, MainVolume, SDVolume); /* Which yields the search direction. */ Lambda = LinksSDLineSearch ( MainVolume, SDVolum

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值