天天看點

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

點選藍字

關注我們

提高知識水準

Transportation

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

導語

交通配置設定(Traffic Assignment)是交通規劃四階段法的最後一步,其中随機配置設定(Stochastic Assignment)方法是其中一種重要方法。随機配置設定主要基于“随機使用者平衡”(Stochastic User Equilibrium,即SUE)理論實作,其中Dial算法是求解SUE問題的經典方法,具體理論可參考《交通規劃》教材中随機交通配置設定方法。本文在簡單回顧Dial理論方法基礎上,通過C++進行程式設計實作,同時設計了求解交通問題簡單的并行計算結構,為求解大規模問題提供思路與方法。

1.首先回顧Dial方法的基本理論與基本步驟。

2.其次較為詳盡解釋Dial方法的C++程式設計實作。

3.最後提供程式完整解決方案,供大家交流學習。

一、Dial算法具體步驟

步驟1:初始化。确定有效路段和有效徑路.

•計算從起點r到所有節點的最小阻抗,記為r(i);

•計算從所有節點到終點s的最小阻抗,記為s(i);

•定義Q(i)為路段起點為i的路段終點的集合;

•定義D(i)為路段終點為i的路段起點的集合;

•對每個路段(i,j),根據下式計算路段似然值(Link Likelyhood)L(i,j)(此時通常假定參數b=1):

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

由該公式可知,凡是似然值等于0的路段都是不合理路段,不應該考慮包含它們的徑路;凡是似然值大于0的路段都可以考慮包含在有效徑路中;當某徑路包含的所有路段的似然值都等于1時,該徑路必然是最小阻抗徑路。

Dial算法中的初始化過程,實質就是確定出行量配置設定在使其有效地遠離其起始節點的 徑路上,那些“走回頭路”的徑路将被剔除掉。

步驟2:從起點r開始按照r(i)升的順序,向前計算路段權重。

從起點r開始,按照r(i)的上升順序依次考慮每個節點,對每個節點,計算離開它的所有路段的權重值,對于節點i,其權重W(i,j), j∈O(i)的計算公式為:

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

當達到終點s,即i=s時就停止權重的計算。

步驟3:從終點$開始,按照s(Q上升的順序,向後計算路段交通量。

從終點s開始,按照s(j)的上升順序依次考慮每個節點,對每個節點,計算進入它的所 有路段的交通量,對于節點j,其交通量x(i,j), i∈D(i)的計算公式為:

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

當達到起點r,即j=r時停止計算。同時整個算法結束。式中qn表示從起點r到終點s的OD交通量,方括号内的交通量之和是指節點j的所有下遊路段上的交通量之和,它們應該先于路段(i,j)上的交通量事先已被計算出來,是以節點的考慮順序是十分重要的。

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

二、資料檔案準備

下面以邵春福老師《交通規劃》中的路網為案例基礎進行資料檔案準備。檔案包括路網節點node.csv,路網連接配接邊link.csv,待配置設定流量demand.csv。

(1)node.csv

Node資料存儲于g_node_vector中,字段資訊如下圖所示:

node_id zone_id x_coord y_coord
結點編号 O-D小區編号 X坐标 Y坐标

(2)link.csv

Link資料存儲于g_node_vector中,字段資訊如下圖所示:

road_link_id from_node_id to_node_id length
路段編号 路段起點編号 路段終點編号 路段長度

(3)demand.csv

Demand資料存儲與g_demand_vector中,字段資訊如下圖所示:

from_zone_id to_zone_id number_of_agents
交通量産生地點編号 交通量吸引地點編号 交通量數值

路網拓撲結構如圖所示:

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...
bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

三、簡單C++程式設計實作

大規模交通網絡具有較高的計算挑戰性,C++在計算效率上具有相當優勢,故本次使用C++語言實作,不熟悉C++語言也可通過Python實作算法過程。鑒于篇幅原因,這裡僅介紹程式整體架構及Dial核心算法部分,I/O等可參見源程式,若有一定C++基礎更易了解。

(1)主程式架構

int main(){// 1.計時開始clock_t startTime, endTime;startTime = clock();// 2.讀取資料,包括node,link,demandg_ReadInputData(assignment);// 3.準備所有node間最短徑路,使用Floyd方法實作(Dijstra亦可)float** adj_matrix;adj_matrix = floyd_SPP(assignment);// 4.重置所有link流量為0for (int i = 0; i < assignment.g_number_of_links; i++) {          g_link_vector[i].volume = 0.0;}// 5.準備并行計算工作台(将所有OD均分至工作台)g_assign_computing_tasks_to_memory_blocks_SUE(assignment);// 6.每個工作台中分别執行Dial算法#pragma omp parallel forfor (int ProcessID = 0; ProcessID < g_NetworkForSUE_vector.size(); ProcessID++) {          g_NetworkForSUE_vector[ProcessID]->Dial_subnetwork(adj_matrix);}
           
bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...
// 7.輸出結果g_output_RLSUE_result(assignment);// 8. 釋放Floyd矩陣記憶體空間for (int i = 0; i < assignment.g_number_of_nodes; i++){          delete[] adj_matrix[i];}delete[] adj_matrix;// 9.計時停止endTime = clock();cout << "The run time is: " << (double)(endTime - startTime) / CLOCKS_PER_SEC << "s" << endl;system("pause");return 1;}
           

(2)基本函數及類的聲明與定義

(按程式中順序介紹)

1.class Assignment {……}: 該項計算任務的全局資訊,node數量、link數量、demand數量、并行計算線程數量等。

2.class CNode, class Clink, class CDemand: 分别存儲node(編号、坐标等),link(兩端節點,長度等),demand(O-D點,需求量等)的相關資訊。

3. class CCSVParser {……}:讀取CSV檔案工具類。(這個可以不看,主要作用是打開檔案,找到字段等)。

4.void g_ReadInputData() {……}:讀取node.csv, link.csv, demand.csv檔案,将資料存至相應vector中。

5.class NetworkForSUE {……}: Dial算法核心類(劃重點!),算法實作的核心步驟在這裡面,後面會詳細介紹。

6. void g_assign_computing_tasks_to_memory_blocks_SUE(){……}: 實作并行計算核心步驟,将不同O-D對平均分在不同工作台運算,提高計算機利用效率,如下圖所示(以4核為例):

7. void fopen_ss()/void g_output_RLSUE_result():結果輸出函數,輸出result.csv,主要包括每條link的流量資訊。

8. float** floyd_SPP():計算最短路徑函數,通過Floyd方法。

(3)Dial算法實作程式

Dial算法程式可在NetworkForSUE類下找到,即:void Dial_subnetwork()。代碼有些長,希望你可以看完。

第1步:初始化。

第1.1步:計算從起點r到所有節點的最小阻抗,記為r(i);計算從所有節點到終點s的最小阻抗,記為s(i)。程式中分别使用數組r_distance[i], s_distance[i]表示。

for (int i = 0; i < assignment.g_number_of_nodes; i++) {        r_distance[i] = adj_matrix[origin_node][i]; //阻抗值從flyod矩陣獲得。}for (int i = 0; i < assignment.g_number_of_nodes; i++) {        s_distance[i] = adj_matrix[i][destination_node];}
           

第1.2步:計算路段似然值,計算公式見上文。

for (int i = 0; i < assignment.g_number_of_links; i++) {        int from_node = g_link_vector[i].from_node_seq_no;        int to_node = g_link_vector[i].to_node_seq_no;        if ((r_distance[from_node] < r_distance[to_node]) & (s_distance[from_node] > s_distance[to_node])) {                 link_likelihood[i] = exp((float)b_Dial * (r_distance[to_node] - r_distance[from_node] - (float)g_link_vector[i].lenth));        if (link_likelihood[i] > 65535) //由于部分最短路不存在(距離無窮大),造成似然值過大,此時該路不可行,似然值為零                 link_likelihood[i] = 0;        }}
           

第2步:從起點r開始按照r(i)上升的順序,向前計算路段權重。

第2.1步:按照r(i)上升的順序形成排序數組,記錄于r_index_list中。按照s(i)上升的順序形成排序數組,記錄于s_index_list中。

vector<float> r_distance_array(assignment.g_number_of_nodes);for (int i = 0; i < r_distance_array.size(); i++)    r_distance_array[i] = r_distance[i];vector<size_t> r_index_list;r_index_list = sort_indexes_e(r_distance_array); //排序函數,傳回值為關鍵字排序後的下标值vector<float> s_distance_array(assignment.g_number_of_nodes);for (int i = 0; i < s_distance_array.size(); i++)    s_distance_array[i] = s_distance[i];vector<size_t> s_index_list;s_index_list = sort_indexes_e(s_distance_array);
           

第2.2步:按照r_index_list的順序,向前逐漸計算路段權重。

for each (int node_list_no in r_index_list) { //按照r_index_list的順序周遊    CNode node = g_node_vector[node_list_no];    if (node.node_seq_no == destination_node) //周遊至終點停止        break;    if (node.node_seq_no == origin_node) { //周遊至起點開始,起點的label為0.        for each (int next_link in g_node_vector[node_list_no].m_outgoing_link_seq_no_vector) {            link_weights_list[next_link] = link_likelihood[next_link];            link_label_list[next_link] = 1.0;        }    }    else { //周遊至中間點時,計算相應label及label_sum資訊。        float current_label = 1.0;        for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) { //周遊該點的入弧及前繼結點。            if (link_likelihood[pre_link] > 0) {                current_label = current_label * link_label_list[pre_link];            }        }        if (current_label == 1.0) {            for each (int next_link in g_node_vector[node_list_no].m_outgoing_link_seq_no_vector) {//計算離開該店所有路段的權重值                float sum_w = 0.0;                for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) {                    sum_w += link_weights_list[pre_link];                }                link_weights_list[next_link] = link_likelihood[next_link] * sum_w;                link_label_list[next_link] = 1.0;            }        }    }}
           

第3步:從終點s開始按照s(i)上升的順序,向後計算路段交通量。

for each (int node_list_no in s_index_list) { //按照s_index_list的順序周遊    CNode node = g_node_vector[node_list_no];    if (node.node_seq_no == origin_node) //遇到起點時停止        break;    if (node.node_seq_no == destination_node) { //周遊至終點時開始,初始化路段flow,label等資訊,并在終點進行第一次配置設定。        float sum_w_pre = 0.0;        for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) {            sum_w_pre += link_weights_list[pre_link];        }        for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) {//計算進入終點所有路段交通量            if (link_weights_list[pre_link] == 0) {                link_flow_list[pre_link] = 0.0;                link_label_list[pre_link] = 1.0;            }            else {                link_flow_list[pre_link] = (link_weights_list[pre_link] / sum_w_pre) * demand_volume;                link_label_list[pre_link] = 1.0;            }        }    }    else { //周遊至中間點時,累計流量和并配置設定流量。        float current_label = 1.0;        float sum_flow = 0.0;        for each (int next_link in g_node_vector[node_list_no].m_outgoing_link_seq_no_vector) {            if (link_weights_list[next_link] > 0) {                current_label = current_label * link_label_list[next_link];                sum_flow += link_flow_list[next_link];            }        }        if (current_label == 1.0) {            float sum_w_pre = 0.0;            for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) {                sum_w_pre += link_weights_list[pre_link];            } //計算路段權重和            for each (int pre_link in g_node_vector[node_list_no].m_incoming_link_seq_no_vector) {//計算進入該點所有路段交通量                if (link_weights_list[pre_link] == 0) {                    link_flow_list[pre_link] = 0.0;                    link_label_list[pre_link] = 1.0;                }                else {                    link_flow_list[pre_link] = (link_weights_list[pre_link] / sum_w_pre) * sum_flow;                    link_label_list[pre_link] = 1.0;                }            }        }    }}
           

(4)程式運作結果

程式運作結果在result.csv中,包括每條link上的流量資訊,字段如下圖所示,進一步可以根據BPR函數計算相應的Travel Time,可自行嘗試。

link_id from_node_id to_node_id volume
1 1 2 647.948
……

(5)大規模路網案例測試

在較大網絡Chicago測試中,node數量933,link數量2950,demand數量23210。網絡拓撲結構如圖所示:

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

在i5低壓雙核筆記本,x64 release模式下不同線程計算結果如下:

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

受限于本人程式設計水準及裝置,代碼效率仍有一定提升優化空間,歡迎大佬或有高性能計算機小夥伴修改與測試。

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

結束語

文中提到的源程式已上傳至Github平台,網址:https://github.com/EntaiWang99/Learning-Network/tree/master/Dial_SUE。大家可以直接打開解決方案運作,注意需打開OpenMP功能,歡迎感興趣的小夥伴讨論與交流。

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

參考

[1]  邵春福. 交通規劃原理[M]. 中國鐵道出版社, 2004.

[2]  STALite: https://github.com/xzhou99/STALite

編輯:莊桢

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

“交通科研Lab”:分享學習點滴,期待科研交流!

bucket sort sample sort 并行_随機交通配置設定經典Dial算法基本原理及C++并行計算實作...

如果覺得還不錯

點這裡!???

繼續閱讀