图论模型和算法

最新推荐文章于 2023-04-14 21:03:32 发布

爱笑的Gary哥

最新推荐文章于 2023-04-14 21:03:32 发布

阅读量584

点赞数 1

分类专栏：数学建模

本文链接：https://blog.csdn.net/Garyboyboy/article/details/113485130

版权

数学建模专栏收录该内容

6 篇文章 8 订阅

订阅专栏

图论模型和算法

Gary哥哥的哥哥的哥哥 2021.1.30

本博客聚焦于matlab解决与图论相关的数学建模问题
读者在阅读之前，最好对数据结构或者算法设计与分析的图论知识有所了解

图论基本概念

假设读者有相关基础，因此图论最基本的概念这里省略相关知识可见此博客

在这里插入图片描述

说明：

接近中心度： v顶点到所有顶点的最短路径的总和的倒数，越大证明越“中心”
中间中心度（表示v的频率）： s到t顶点的最短路径数量，分子为路径中经过了v顶点的最短路径数量
特征向量中心度：（求矩阵特征向量可得）

graphallshortestpaths	求图中所有顶点对之间的最短距离
graphconnredcomp	找无 (有) 向图的 (强/弱) 连通分支
graphisreddag	测试有向图是否含有圈
graphisomorphism	确定一个图是否有生成树
graphmaxflow	计算有向图的最大流
graphminspantree	在图中找最小生成树
graphpred2path	把前驱顶点序列变成路径的顶点序列
graphshortestpath	求指定一对顶点间的最短距离和路径
graphtopoorder	执行有向无圈图的拓扑排序
graphtraverse	求从一顶点出发, 所能遍历图中的顶点

网络分析工具箱

为了方便计算一些更加复杂的图论中的变量，我给出一个功能强大的工具箱
非Matlab官方的，是MIT开发的一个工具箱，地址（那个zip文件就是了）

函数名	功能
degrees	求图中所有顶点的度，入度和出度
ave_neighbor_deg	求图中所有顶点的相邻顶点平均度
closeness	求图中所有顶点的接近中心度
node_betweenness_faster	求图中所有顶点的中间中心度
edge_betweenness	求图中所有边的中间中心度
…	…

案例

在这里插入图片描述

n = 10; %顶点总数

% 给每个人都编号
Andre = 1; Betty = 2; Carol = 3; Dave  = 4; Ed    = 5;
Fanny = 6; Garth = 7; Hale  = 8; Ike   = 9; Jane  =10;

% 根据图构造邻接矩阵
A = zeros(10);
A(Andre, [Betty, Carol,  Dave, Fanny]) = 1;% Andre跟哪几个人有边
A(Betty, [Andre,  Dave,    Ed, Garth]) = 1;
A(Carol, [Andre,  Dave, Fanny]) = 1;
A( Dave, [Andre, Betty, Carol,    Ed, Fanny, Garth]) = 1;
A(   Ed, [Betty,  Dave, Garth]) = 1;
A(Fanny, [Andre, Carol,  Dave, Garth,  Hale]) = 1;
A(Garth, [Betty,  Dave,    Ed, Fanny,  Hale]) = 1;
A( Hale, [Fanny, Garth,   Ike]) = 1;
A(  Ike, [ Hale,  Jane]) = 1;
A( Jane, [  Ike]) = 1;

Cd = degrees(A)' /(n-1)   % 计算点度中心度并且标准化. 
Cc = closeness(A)*(n-1)   % 计算接近中心度并标准化

最小哈密尔顿路径问题

function [order, totdist] = minhamiltonpath(D)
%Step 1: Setup linear progrm
N = size(D, 1);
idxs = nchoosek(1:N, 2);
dist = D(sub2ind([N, N], idxs(:, 1), idxs(:, 2)));
lendist = length(dist);

%There need to be two "trips" atached to each "stop"
Aeq = spalloc(N, length(idxs), N*(N-1));
for i = 1:N
    whichIdxs = idxs(:,1) == i | idxs(:,2) == i;
    Aeq(i,whichIdxs) = 1;
end
beq = 2*ones(N, 1);

intcon = 1:lendist;
lb = zeros(lendist, 1);
ub = ones(lendist, 1);

opts = optimoptions('intlinprog','Display','off');

%Step 2: Solve linear program, adding inequality constraints until
%subtours are eliminated from found solution
x = intlinprog(dist,intcon,[],[],Aeq,beq,lb,ub,opts);

tours = detectSubtours(x, idxs);
ntours = length(tours);

A = spalloc(0, lendist, 0);
b = [];
while ntours > 1
    fprintf(1, '%i subtours\n', ntours);
    for i = 1:ntours
        rowIdx = size(A, 1)+1;
        touri = tours{i}; % Extract the current subtour
        variations = nchoosek(1:length(touri), 2);
        for j = 1:length(variations)
            whichVar = (sum(idxs==touri(variations(j,1)),2)) & ...
                       (sum(idxs==touri(variations(j,2)),2));
            A(rowIdx, whichVar) = 1;
        end
        b(rowIdx) = length(touri)-1;
    end
    x = intlinprog(dist,intcon,A,b,Aeq,beq,lb,ub,opts);
    tours = detectSubtours(x, idxs);
    ntours = length(tours);
end
totdist = dist'*x;
order = tours{:};

function subTours = detectSubtours(x,idxs)
% 返回子旅途元胞数组，即图中的子圈
 
x = round(x); % 纠正不确切的整数
r = find(x);
substuff = idxs(r,:); % 旅行线的节点对
unvisited = ones(length(r),1); % 跟踪未访问的旅行
curr = 1; % 正在评价的子旅途
startour = find(unvisited,1); % 第一个未访问的旅行
while ~isempty(startour)
    home = substuff(startour,1);
    nextpt = substuff(startour,2);
    visited = nextpt;
    unvisited(startour) = 0;
    
    while nextpt ~= home
        % 找以nextpt为起点的旅行
        [srow,scol] = find(substuff == nextpt);
        
        % 确定相应旅行的节点
        trow = srow(srow ~= startour);
        scol = 3-scol(trow == srow); % 1变2，2变1
        startour = trow;
        nextpt = substuff(startour,scol); % 子旅途的下一节点位置
        
        visited = [visited,nextpt]; % 将节点加入子旅途
        unvisited(startour) = 0; % 更新访问过的位置
    end
    subTours{curr} = visited; % 保存找到的子旅途
    
    curr = curr + 1;
    startour = find(unvisited,1);
end
end

tic
lat = [ 39.54  31.14  39.09  29.32  45.45 ...	   
        43.52  41.50  40.49  38.02  37.52 ...	
        36.38  34.48  34.16  36.03  38.20 ...	
        36.38  43.48  31.51  32.02  30.14 ...	
        28.11  28.41  30.37  30.39  26.35 ...	
        26.05  23.08  20.02  22.48  25.00 ...
        29.39  22.18  22.14  25.03          ];

lon = [116.28 121.29 117.11 106.32 126.41 ...	  
       125.19 123.24 111.48 114.28 112.34 ...
       117.00 113.42 108.54 103.49 106.16 ...	
       101.45  87.36 117.18 118.50 120.09 ...
       113.00 115.52 114.21 104.05 106.42 ...	
       119.18 113.15 110.20 108.20 102.41 ...
        90.08 114.10 113.35 121.31          ];

n = length(lat);

R = 6378.137;

dist = zeros(n);
for i = 1:n
    for j = i+1:n
        dist(i,j) = distance(lat(i),lon(i), lat(j),lon(j), R);
    end
end

dist = dist + dist';
[order,totdist] = minhamiltonpath(dist);
plot(lon(order([1 : end ,1])), lat(order([1:end,1])),'o-')

灾情巡视路径问题 (CUMCM1998B)

教程（pdf最后）

代码

三个县长把所有乡村跑一遍

兵分三路，总路程短一点，大家差不多时间一起停

分三个子图

分的方案很重要
先得出以o为顶点的最小生成树
再去分配

function [ncomp, icomp] = roadtree(ifplot)
if nargin<1; ifplot = 'plot'; end

[w, x,y, n, O] = road2graph('original', 'noplot');
[dist, path, pred] = graphshortestpath(w, O, 'directed', false);

I = setdiff(1:n, O);
J = pred(I);

tree = zeros(n);
for k = 1:length(I)
    i = I(k); j = J(k);
    tree(i,j) = w(i,j); tree(j,i) =w(j,i);
end

treeO = tree; treeO(:,O) = 0; treeO(O,:) = 0;

%Find strongly or weakly connected components in graph
[ncomp, icomp] = graphconncomp(sparse(treeO),'directed',false);

if strcmp(ifplot,'plot')
    color = {'r','g','b','m','k','c'};
    for k = 1:length(I)
        i = I(k); j = J(k);
        plot([x(i),x(j)], [y(i),y(j)]',['o-', color{icomp(i)}])
        hold on
    end
end

function [w, x,y, n, O] = road2graph(gtype, ifplot)
if nargin<1; 
    gtype = 'original'; ifplot = 'plot'; 
elseif nargin<2; 
    ifplot = 'plot';
end


nvil  = 35;   % number of villages
ntown = 18;   % number of towns
n = nvil + ntown;

ID = num2cell(nvil+1:nvil+ntown); % index of towns
[A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R] = deal(ID{:});

%% coordinates of villages & towns for plot
x = [61.7 51.3 52.0 51.3 44.4 35.9 33.3 35.5 22.9 20.9 ... %  1-10
     18.9  9.4  9.1  1.6  0.0  8.2 13.9 15.1 21.1 26.1 ... % 11-20
     24.7 20.3 28.6 29.9 30.5 43.5 44.1 50.6 59.0 60.8 ... % 21-30
     64.4 68.9 68.1 75.1 78.5 67.2 65.2 58.4 45.9 28.6 ... % 31-35, A-E
     19.9 13.1  1.9  8.1 14.9 20.9 26.5 39.3 34.9 56.3 ... %  F- O
     50.5 56.2 61.1                                    ]'; %  P- R
 
y = [28.9 23.6 19.0  6.9 24.6 25.4 19.2 13.3 12.7  0.0 ... %  1-10
     19.2 11.8 24.3 26.7 39.7 47.2 45.9 35.6 28.3 33.6 ... % 11-20
     39.8 48.6 46.7 54.2 37.4 41.7 47.7 45.5 45.7 53.7 ... % 21-30
     45.8 50.5 40.5 36.5 42.7 34.2 25.9 20.9 15.9 15.3 ... % 31-35, A-E
      8.8 17.6 17.6 37.7 28.6 41.2 28.7 32.6 44.3 29.6 ... %  F- O
     35.9 50.8 38.8                                    ]'; %  P- R

%% distance matrix 
w = zeros(n);

w(1, A) = 10.3; w(1, B) =  5.9; w(1, C) = 11.2; w(1, O) =  6.0; 
w(2, 3) =  4.8; w(2, 5) =  8.3; w(2, O) =  9.2;
w(3, C) =  7.9; w(3, D) =  8.2;
w(4, 8) = 20.4; w(4, D) = 12.7;
w(5, 6) =  9.7; w(5, D) = 11.3; w(5, M) = 11.4;

w(6, 7) =  7.3; w(6, L) = 11.8; w(6, M) =  9.5;
w(7, D) = 15.1; w(7, E) =  7.2; w(7, L) = 14.5;
w(8, E) =  8.0;
w(9, E) =  7.8; w(9, F) =  5.6;
w(10,F) = 10.8;

w(11,E) = 14.2; w(11,G) =  6.8; w(11,J) = 13.2;
w(12,F) = 12.2; w(12,G) =  7.8; w(12,H) = 10.2;
w(13,14)=  8.6; w(13,G) =  8.6; w(13,I) = 16.4; w(13,J) = 9.8;
w(14,15)= 15.0; w(14,H) =  9.9;
w(15,I) =  8.8;

w(16,17)=  6.8; w(16,I) = 11.8;
w(17,22)=  6.7; w(17,K) =  9.8;
w(18,I) =  8.2; w(18,J) =  8.2; w(18,K) =  9.2;
w(19,20)=  9.3; w(19,J) =  8.1; w(19,L) =  7.2;
w(20,21)=  7.9; w(20,25)=  6.5; w(20,L) =  5.5;

w(21,23)=  9.1; w(21,25)=  7.8; w(21,K) =  4.1;
w(22,23)= 10.0; w(22,K) = 10.1;
w(23,24)=  8.9; w(23,N) =  7.9;
w(24,27)= 18.8; w(24,N) = 13.2;
w(25,M) = 12.0; w(25,N) =  8.8;

w(26,27)=  7.8; w(26,N) = 10.5; w(26,P) = 10.5;
w(27,28)=  7.9; 
w(28,P) = 12.1; w(28,Q) =  8.3;
w(29,P) = 15.2; w(29,Q) =  7.2; w(29,R) =  7.9;
w(30,32)= 10.3; w(30,Q) =  7.7;

w(31,32)=  8.1; w(31,33)=  7.3; w(31,R) =  9.2;
w(32,33)= 19.0; w(32,35)= 14.9;
w(33,35)= 20.3; w(33,A) =  7.4;
w(34,35)=  8.2; w(34,A) = 11.5; w(34,B) = 17.6;

w(A, B) = 12.2; w(A, R) =  8.8;
w(B, C) = 11.0; 
w(I, J) = 15.8;
w(C, O) = 11.5;
w(M, N) = 14.2; w(M, O) = 19.8;

w(O, P) = 10.1; w(O, R) = 12.9;

w = sparse(w'); % upper triangular matrix to lower trianglular matrix

%% plot graph
if strcmp(ifplot,'plot')
    [i,j] = find(w>0);
    plot([x(i),x(j)]', [y(i),y(j)]','o:b','linewidth',4);
end
%%
if strcmp(gtype,'complete')
   w = graphallshortestpaths(w, 'directed', false);
end

总路程：三个子图最小哈密尔顿回路最小
均衡化最小哈密尔顿回路三个中max-min小
具体过程见上面的pdf和代码

爱笑的Gary哥

关注

1
点赞
踩
5

收藏

觉得还不错? 一键收藏
1
评论
图论模型和算法

图论模型和算法Gary哥哥的哥哥的哥哥 2021.1.30本博客聚焦于matlab解决与图论相关的数学建模问题读者在阅读之前，最好对数据结构或者算法设计与分析的图论知识有所了解图论基本概念假设读者有相关基础，因此图论最基本的概念这里省略相关知识可见此博客说明：接近中心度： v顶点到所有顶点的最短路径的总和的倒数，越大证明越“中心”中间中心度（表示v的频率）： s到t顶点的最短路径数量，分子为路径中经过了v顶点的最短路径数量特征向量中心度：（求矩阵特征向量可得）相关
复制链接

扫一扫