fleury MATLAB,数学建模部分源码分享
写在前面
数学建模竞赛已经结束了。由于种种原因最终无缘国奖,还是很难受的。结合我的经历,有一些想法想告诉所有正在准备数学建模竞赛的同学。
数学建模竞赛它不会绝对相信实力,运气真的是必须考虑的一部分。
无论数学建模是否拿奖都受益终生,它会极大地提高工科生思维和建模能力。
如果没有必要千万不要熬夜,状态会很差很差很差,得不偿失。
心态爆炸是导致数学建模竞赛失败的大部分原因。
找个好队友真的extremely重要
正题
最近准备换系统,把我们在历次比赛中的代码拿出。反正删了也是删了,还不如传上来以后还能回味一下。如果中间有bug或error欢迎反馈。原谅我不写注释。。。
在此真诚感谢数学建模32种方法的整理者,这里我只提供一个下载链接
层次分析法AHP
clear, clc;
prin = xlsread('AHP.xlsx', 'A1:D4');
price = xlsread('AHP.xlsx', 'A6:C8');
oc = xlsread('AHP.xlsx', 'A10:C12');
conf = xlsread('AHP.xlsx', 'A14:C16');
surf = xlsread('AHP.xlsx', 'A18:C20');
[Vpric, Dpric] = eig(prin);
DpricMax = max(max(Dpric));
CIpric = (DpricMax - length(prin)) ./ (length(prin) - 1);
RI = 0.90;
CRpric = CIpric ./ RI;
if CRpric < 0.10 %#ok
disp(['准则层CR = ', num2str(CRpric), ',接受误差']);
end
wCR = sum(prin, 2)';
indexwCRmax = find(wCR == max(wCR));
switch indexwCRmax
case 1
disp('最大影响因素:price');
case 2
disp('最大影响因素:oil consumption');
case 3
disp('最大影响因素:conformance');
case 4
disp('最大影响因素:surface');
end
wPrice = sum(wCR(1) .* price, 2);
wOc = sum(wCR(2) .* oc, 2);
wConf = sum(wCR(3) .* conf, 2);
wSurf = sum(wCR(4) .* surf, 2);
wScheme = sum([wPrice, wOc, wConf, wSurf], 2);
disp(['最倾向第',num2str(find(wScheme == max(wScheme))) , '辆车']);
单源最短路Dijkstra算法
function [pathMat, pathWeightMat] = Dijkstra(weightmat, startpoint)
%
%function Dijkstra receives two arguments, the first one is the matrix of path weight,
%the second one is the startpoint.
%function Dijkstra returns two results, pathMat is the matrix of the
%shorest path, pathWeightMat is the matrix of weight
%an example:
% >> m. at = [
% 0 1 5 Inf Inf Inf Inf Inf Inf
% 1 0 3 7 5 Inf Inf Inf Inf
% 5 3 0 Inf 1 7 Inf Inf Inf
% Inf 7 Inf 0 2 Inf 3 Inf Inf
% Inf 5 1 2 0 3 6 9 Inf
% Inf Inf 7 Inf 3 0 Inf 5 Inf
% Inf Inf Inf 3 6 Inf 0 2 7
% Inf Inf Inf Inf 9 5 2 0 4
% Inf Inf Inf Inf Inf Inf 7 4 0]
% >> [path, weight] = Dijkstra(mat, 1)
% path =
%
% 1 0 0 0 0 0 0 0 0
% 1 2 0 0 0 0 0 0 0
% 1 2 3 0 0 0 0 0 0
% 1 2 3 5 4 0 0 0 0
% 1 2 3 5 0 0 0 0 0
% 1 2 3 5 6 0 0 0 0
% 1 2 3 5 4 7 0 0 0
% 1 2 3 5 4 7 8 0 0
% 1 2 3 5 4 7 8 9 0
%
%
weightMat = weightmat;
if ismember(1, weightMat ~= weightMat')
[row, col] = find((weightMat ~= weightMat') == 1);
error(['Matrix input illegal, ', '(', num2str(row(1)), ',', num2str(col(1)),...
')and(', num2str(row(2)) , ',', num2str(col(2)), ')']);
end
unfinishMat = ones(1, length(weightMat(1, :)));
pathWeightMat = inf .* ones(length(weightMat(1, :)), 1);
pathMat = zeros(size(weightMat));
startPoint = startpoint;
presPoint = startPoint;
presPathWeight = 0;
path = presPoint;
pathMat(presPoint, 1) = presPoint;
while ismember(1, unfinishMat)
i = 1;
unfinishMat(presPoint) = 0;
tempPathMat = weightMat(presPoint, :) .* unfinishMat;
tempPathMat(tempPathMat == 0) = inf;
nextPoint = find(tempPathMat == min(tempPathMat), 1, 'last');
if presPoint == nextPoint
break;
end
while i <= length(unfinishMat)
if tempPathMat(i) ~= inf && ~isnan(tempPathMat(i)) && i ~= presPoint
tempPath = [path, i];
tempWeight = presPathWeight + tempPathMat(i);
if tempWeight < pathWeightMat(i)
pathWeightMat(i) = tempWeight;
pathMat(i, :) = [tempPath, zeros(1, length(pathMat(i, :)) - length(tempPath))];
end
end
i = i + 1;
end
path = pathMat(nextPoint, pathMat(nextPoint, :) ~= 0);
presPathWeight = pathWeightMat(nextPoint);
presPoint = nextPoint;
end
pathWeightMat(startPoint) = 0;
欧拉回路Fleury算法
参考后封装的,原文太久忘了,抱歉
function EulerCircuit = Fleury(G, st)
G(G == inf) = 0;
if nargin == 2
presPoint = st;
circuit = st;
else
presPoint = 1;
circuit = 1;
end
while sum(any(G))
altDeg = (sum(G) - G(presPoint, :)) .* (G(presPoint, :) > 0)
if ~any(altDeg)
[~, nextPoint] = find(G ~= 0)
circuit = [circuit, nextPoint(1)]
break
end
nextPoint = find(altDeg == max(altDeg), 1)
circuit = [circuit, nextPoint]
G(presPoint, nextPoint) = G(presPoint, nextPoint) - 1
G(nextPoint, presPoint) = G(presPoint, nextPoint)
presPoint = nextPoint
G
altDeg = (sum(G, 2) - G(:, presPoint)) .* (G(:, presPoint) > 0)
if ~any(altDeg)
[nextPoint, ~] = find(G ~= 0)
circuit = [circuit, nextPoint(1)]
break
end
nextPoint = find(altDeg == max(altDeg), 1)
circuit = [circuit, nextPoint]
G(presPoint, nextPoint) = G(presPoint, nextPoint) - 1
G(nextPoint, presPoint) = G(presPoint, nextPoint)
presPoint = nextPoint
G
end
EulerCircuit = circuit;
完全最短路Flory算法
function [pathMat, pathWeightMat] = Floyd(weightMat, startpoint, endpoint)
% Floyd find out the shorest path between two given vertex
%
% [a, b] = Floyd(G, s, e)
%
% G: the graph
% s: start point
% e: end point
% a: the shorest path between s and e
% b: the path length of a
%
% Example:
% >> G = [ 0 12 19 28 40 59
% 12 0 13 20 29 41
% 19 13 0 14 21 30
% 28 20 14 0 15 22
% 40 29 21 15 0 16
% 59 41 30 22 16 0];
% s = 1; e = 6;
% >> [a, b] = Floyd(G, s, e)
% a =
% 1 3 6
% b =
% 49
%
if ismember(1, weightMat ~= weightMat')
[row, col] = find((weightMat ~= weightMat') == 1);
error(['Matrix input illegal, ', '(', num2str(row(1)), ',', num2str(col(1)),...
')and(', num2str(row(2)) , ',', num2str(col(2)), ')']);
end
n=length(weightMat);
U=weightMat;
m=1;
%Floyd algorithm
while m<=n
for i=1:n
for j=1:n
if U(i,j)>U(i,m)+U(m,j)
U(i,j)=U(i,m)+U(m,j);
end
end
end
m=m+1;
end
pathWeightMat = U(startpoint,endpoint);
P1=zeros(1,n);
k=1;
P1(k)=endpoint;
V=ones(1,n)*inf;
kk=endpoint;
while kk~=startpoint
for i=1:n
V(1,i)=U(startpoint,kk)-weightMat(i,kk);
if V(1,i)==U(startpoint,i)
P1(k+1)=i;
kk=i;
k=k+1;
end
end
end
k=1;
wrow=find(P1~=0);
for j=length(wrow):(-1):1
pathMat(k)=P1(wrow(j));
k=k+1;
end
遗传算法(求解哈密顿图的一种优化算法)
注意,此算法可能主要部分不是遗传算法,而是其中的改良圈算法。学校教材代码,无法验证但结果不错,注释部分可删
clc,clear;
close all;
load data;
data = [data; [data(1, 1), data(1, 2)]]; % 数据部分,其中第一行为起点和终点(环)
x = data(:, 1);
y = data(:, 2);
fullData = data;
codeLen = length(x);
libSize = 50;
dist = zeros(codeLen);
for i = 1: codeLen - 1
for j = i + 1: codeLen
pathLen = sqrt((x(i) - x(j)) ^ 2 + (y(i) - y(j)) ^ 2);
dist(i,j) = pathLen;
end
end
dist = dist + dist';
%通过改良圈算法选取优良父代A
% rec = [];
for k = 1: libSize
c = randperm(codeLen - 2);
c1 = [1, c + 1, codeLen];
flag = 1;
while flag > 0
flag = 0;
for m = 1: codeLen - 3
for n = m + 2: codeLen - 1
if dist(c1(m), c1(n)) + dist(c1(m + 1), c1(n + 1)) < dist(c1(m), c1(m + 1)) + dist(c1(n), c1(n + 1))
flag = 1;
c1(m + 1:n) = c1(n: - 1:m + 1);
end
end
end
% temp = 0;
% for i = 1: length(c1) - 1
% temp = temp + dist(c1(i), c1(i + 1));
% end
% rec = [rec, temp];
end
lib(k,c1) = 1: codeLen;
end
[~,index] = sort(lib,2);
pathLen = zeros(length(lib(:, 1)), 1);
for j = 1:length(lib(:, 1))
for i = 1: codeLen - 1
pathLen(j) = pathLen(j) + dist(index(j, i),index(j, i + 1));
end
end
path = index(pathLen == min(pathLen),:);
len = min(pathLen)
X = fullData(path,1);
Y = fullData(path,2);
figure;
plot(X, Y, ' - o')
% % figure;
% % plot(rec)
数据表data.mat
1304 2312
3639 1315
4177 2244
3712 1399
3488 1535
3326 1556
3238 1229
4196 1004
4312 790
4386 570
3007 1970
2562 1756
2788 1491
2381 1676
1332 695
3715 1678
3918 2179
4061 2370
3780 2212
3676 2578
4029 2838
4263 2931
3429 1908
3507 2367
3394 2643
3439 3201
2935 3240
3140 3550
2545 2357
2778 2826
2370 2975
最小生成树Kruskal算法
function MST = Kruskal(G, method)
%Kruskal Get the MST of graph.
%
% MST = Kruskal(G, method)
%
% Kruskal has two arguments. G is necessary while method can be omitted.
% G is the matrix of the graph, in which the elements stands for the path
% weight between the vertexes, and G(i, i) = 0.
% if vertex(i) and vertex(j) are non-adjacent, G(i, j) and G(j, i) could
% be both 0 or inf.
% method points out whether the function will find out the Minimum Spanning Tree
% or the Maximum Spanning Tree. if method = 'min' or omitted, function
% will returns Minimum Spanning Tree, if method = 'max', Maximum Spanning
% Tree will be returned.
% Kruskal returns a Ax3 matrix ,A is the number of edges, column 1 and 2
% are the vertexes of edges, column 3 is the weight of the relevant edge.
%
% Example:
% if G = [0 5 Inf Inf
% 5 0 2 3
% Inf 2 0 4
% Inf 3 4 0]
% then Kruskal(G) is
% 3 2 2
% 4 2 3
% 2 1 5
% ## [4 2 3] means from vertex4 to vertex2 and the distance is 3.
% Kruskal(G, 'max') is
% 2 1 5
% 4 3 4
% 4 2 3
if nargin == 1 || strcmp(method, 'min')
weightMat = tril(G);
weightMat(weightMat == 0) = inf;
checkMat = zeros(size(weightMat));
path = [];
edgeNum = 0;
while ~all(all(weightMat == inf))
[row, col] = find(weightMat == min(min(weightMat)), 1);
checkMat(row, col) = 1;
checkMat(col, row) = 1;
tempMat = checkMat;
while 1
index = find(sum(tempMat, 2) == 1 | sum(tempMat, 2) == 0, 1);
if isempty(index)
checkMat(row, col) = 0;
break;
end
tempMat(index, :) = [];
tempMat(:, index) = [];
if isempty(tempMat)
path = [path; row, col, weightMat(row, col)];
edgeNum = edgeNum + 1;
break;
end
end
weightMat(row, col) = inf;
end
elseif nargin == 2 || strcmp(method, 'max')
weightMat = tril(G);
weightMat(weightMat == inf) = -inf;
weightMat(weightMat == 0) = -inf;
checkMat = zeros(size(weightMat));
path = [];
edgeNum = 0;
while ~all(all(weightMat == -inf))
[row, col] = find(weightMat == max(max(weightMat)), 1);
checkMat(row, col) = 1;
checkMat(col, row) = 1;
tempMat = checkMat;
while 1
index = find(sum(tempMat, 2) == 1 | sum(tempMat, 2) == 0, 1);
if isempty(index)
checkMat(row, col) = 0;
break;
end
tempMat(index, :) = [];
tempMat(:, index) = [];
if isempty(tempMat)
path = [path; row, col, weightMat(row, col)];
edgeNum = edgeNum + 1;
break;
end
end
weightMat(row, col) = -inf;
end
end
MST = path;
最大流问题
这个代码没有验证过,因为这是培训的时候讲的,回来就顺便仿了
function f=maxFlow(u)
f = zeros(size(u));
n = length(u);
maxf(n) = 1;
while maxf(n)>0
maxf = zeros(1,n);
pred = zeros(1,n);
list = 1;
record = list;
maxf(1) = inf;
%list是未检查邻接点的标号点,record是已标号点
while (~isempty(list))&&(maxf(n) == 0)
flag = list(1);
list(1) = [];
label1 = find(u(flag, :) - f(flag, :));
label1 = setdiff(label1, record);
list = union(list, label1);
pred(label1) = flag;
maxf(label1) = min(maxf(flag), u(flag,label1) - f(flag, label1));
record = union(record, label1);
label2 = find(f(:, flag));
label2 = label2';
label2 = setdiff(label2, record);
list = union(list, label2);
pred(label2) = - flag;
maxf(label2) = min(maxf(flag), f(label2, flag));
record = union(record, label2);
end
if maxf(n) > 0
v2 = n;
v1 = pred(v2);
while v2 ~= 1
if v1 > 0
f(v1,v2) = f(v1, v2) + maxf(n);
else
v1 = abs(v1);
f(v2, v1) = f(v2, v1) - maxf(n);
end
v2 = v1;
v1 = pred(v2);
end
end
end
主成分分析PCA
这个算法好像matlab有自带的函数,自行help pca。这个代码只作为原理步骤展示。原文参考地址见代码
%example refer to https://blog.csdn.net/zhongkelee/article/details/44064401
clear, clc;
score = xlsread('PCA.xlsx', 'B2:G10');
%第一步,分别求每列平均值,然后对于所有的样例,都减去对应的均值。
for i = 1:length(score(1,:))
score(:, i) = score(:, i) - mean(score(:, i));
end
%第二步,求特征协方差矩阵
covMatrix = zeros(length(score(1,:)), length(score(1,:)));
for i = 1:length(score(1,:))
for j = 1: length(score(1,:))
temp = cov(score(:, i), score(:, j));
covMatrix(i ,j) = temp(2, 1);
end
end
%第三步,求协方差的特征值和特征向量
[VcovMatrix, DcovMatrix] = eig(covMatrix);
%第四步,将特征值按照从大到小的顺序排序,选择其中最大的k个,然后将其对应的k个特征向量分别作为列向量组成特征向量矩阵
projectDimension = 3; %提取的主成分个数
projectMatrix = zeros(6, projectDimension);
for k = 1:projectDimension
index = find(sum(DcovMatrix == max(max(DcovMatrix))) == 1);
projectMatrix(:, k) = VcovMatrix(:, index);
DcovMatrix(DcovMatrix == max(max(DcovMatrix))) = 0;
end
%第五步,做乘积
dataProjected = score * projectMatrix
排队论/生灭过程
这个代码好像只适合生灭过程,其他类型自行查询
t_arr = zeros(1,21); %每位顾客到达时间
t_leave = zeros(1,21); %每位顾客离去时间
time_wait_c = zeros(1,21); %顾客等待时间累加
time_free_s = zeros(1,21); %收款台空闲时间累加
t_exp = exprnd(2,1,21); %服从指数分布的随机数
t_norm = normrnd(1,1/3,1,21); %服从正态分布的随机数
for i = 1:1:20
t_arr(i+1) = t_arr(i)+t_exp(i);
if t_arr(i+1) >= t_leave(i)
time_wait_c(i+1) = time_wait_c(i);
t_leave(i+1) = t_arr(i+1)+t_norm(i+1);
time_free_s(i+1) = t_arr(i+1)-t_leave(i)+time_free_s(i);
else
time_wait_c(i+1) = t_leave(i)-t_arr(i+1)+time_wait_c(i);
t_leave(i+1) = t_leave(i)+t_norm(i+1);
time_free_s(i+1) = time_free_s(i);
end
end
b = [t_arr',t_leave',time_wait_c',time_free_s']; %队伍情况
[brow,bcol] = size(b);
b = [b,zeros(brow,1)];
for j = 2:brow
l = 0;
if j-l-1 > 0
while b(j,1) <= b(j-l-1,2) l = l+1;
end;
b(j,bcol+1) = l;
end;
end;
ave_wait_time = time_wait_c(end)/length(time_wait_c) %平均等待时间
ave_stay_time = sum(t_leave-t_arr)/length(t_arr) %平均逗留时间
served_per_min = t_leave(end)/length(t_leave) %平均每分钟服务的顾客人数
画六面体的函数
这个函数完全是闲得蛋疼写的,当初做了个堆箱子的题的可视化,没别的办法了就只好上matlab,比较吃内存。所以能找到别的方法尽量还是别用这个吧
function simplot(arg)
x0 = arg(1);
y0 = arg(2);
z0 = arg(3);
x1 = arg(4);
y1 = arg(5);
z1 = arg(6);
d = 1;
color = {'blue', 'red', 'green', 'black', 'yellow'};
color = color{randi(5)};
x = x0: d: x1; % x起始坐标: 间隔: 终止坐标
y = y0: d: y1; % y起始坐标: 间隔: 终止坐标
z = z0: d: z1; % z起始坐标: 间隔: 终止坐标
hold on;
% xoy
[X, Y] = meshgrid(x, y);
Z = z0 * ones(size(X));
surf(X, Y, Z, 'EdgeColor', color)
Z = (z1) * ones(size(X));
surf(X, Y, Z, 'EdgeColor', color)
% yoz
[Y, Z] = meshgrid(y, z);
X = x0 * ones(size(Y));
surf(X, Y, Z, 'EdgeColor', color)
X = x1 * ones(size(Y));
surf(X, Y, Z, 'EdgeColor', color)
% zox
[X, Z] = meshgrid(x, z);
Y = y0 * ones(size(X));
surf(X, Y, Z, 'EdgeColor', color)
Y = y1 * ones(size(X));
surf(X, Y, Z, 'EdgeColor', color)
xlabel('x');
ylabel('y');
zlabel('z');
最后祝各位在数学建模的路上一路顺风?!
第一次写这个玩意儿,markdown不熟,排版见谅。
fleury MATLAB,数学建模部分源码分享相关推荐
- matlab均线程序化交易,【策略分享】Matlab量化交易策略源码分享
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 策略名称:横盘突破策略 策略思路: 日内交易策略,收盘平仓: 横盘突破在过去30根K线的高低点围绕中轴上下0.5%的范围内波动时: 上轨=过去30根K线的 ...
- matlab 量化 策略,【策略分享】Matlab量化交易策略源码分享
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 策略名称:横盘突破策略 策略思路: 日内交易策略,收盘平仓: 横盘突破在过去30根K线的高低点围绕中轴上下0.5%的范围内波动时: 上轨=过去30根K线的 ...
- 量化交易策略matlab交易方案,【策略分享】Matlab量化交易策略源码分享
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 策略名称:横盘突破策略 策略思路: 日内交易策略,收盘平仓: 横盘突破在过去30根K线的高低点围绕中轴上下0.5%的范围内波动时: 上轨=过去30根K线的 ...
- Matlab期货量化交易特征选取,【策略分享】Matlab量化交易策略源码分享
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 策略名称:横盘突破策略 策略思路: 日内交易策略,收盘平仓: 横盘突破在过去30根K线的高低点围绕中轴上下0.5%的范围内波动时: 上轨=过去30根K线的 ...
- matlab量化策略代码,【策略分享】Matlab量化交易策略源码分享
该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 策略名称:横盘突破策略 策略思路: 日内交易策略,收盘平仓: 横盘突破在过去30根K线的高低点围绕中轴上下0.5%的范围内波动时: 上轨=过去30根K线的 ...
- matlab布林算法代码,MATLAB量化交易策略源码分享之 布林通道+高低点
策略原理: 通过布林带以及突破后的高低点的形成产生交易信号 采取跟踪止损出场 回测曲线: 2017-2-27 10:09:43 上传 下载附件 (65.41 KB) 策略代码: function S ...
- 多目标人工秃鹫优化算法(MATLAB源码分享,智能优化算法) 提出了一种多目标版本的人工秃鹫优化算法(AVOA)
多目标人工秃鹫优化算法(MATLAB源码分享,智能优化算法) 提出了一种多目标版本的人工秃鹫优化算法(AVOA),用于多目标优化问题. AVOA的灵感来源于非洲秃鹫的生活方式. 档案.网格和领导者选择 ...
- 【MATLAB第11期】#源码分享 |时间序列数据绘图,横坐标更改为时间轴 横坐标轴参数更改 日期间隔设置 日期标签或格式更改
[MATLAB第11期]#源码分享 |时间序列数据绘图,横坐标更改为时间轴 横坐标轴参数更改 绘图问题解决方案 1.导入数据方式(识别时间数据) 2.案例展示 (1)打开数据 (2)定义时间和数据 ( ...
- 数学建模matlab imread,《matlab数学建模方法与实践》第三章 数据的处理
3.1数据的获取 3.1.1从Excel中读取数据 a= xlsread('路径',n(sheetn),'A(n-1):Bn(读入的数据范围)') 3.1.2从TXT中读取数据 3.1.2.1.1 同 ...
最新文章
- 【预告】腾讯移动分析 MTA 即刻登陆 2017 GMTC 全球移动技术大会
- Windows 7防火墙设置详解(三)
- 二进制安全与MySQL的关系_《高可用MySQL》节选 -- 安全和二进制日志
- 【渝粤教育】国家开放大学2018年秋季 0033-21T工程数学 参考试题
- iOS统计项目的代码总行数
- C++之运算符重载(上)
- 【BZOJ3489】A simple rmq problem kd-tree
- 第四章 自上而下分析
- 计算机操作系统详细学习笔记(一):计算机操作系统概述
- 敏感词过滤算法Aho-Corasick
- Mac上搭建IOS自动化测试环境
- 易支付源码 28k支付第四方支付源码-Oreo支付系统
- Linux第7章Gdk及Cairo基础,Linux第7章Gdk及Cairo基础概要1.ppt
- 软件测试技术 实验一:黑盒测试1
- LinuxOracle R6安装oracle 11gR2详细步骤
- 基于stm32的mpu6050传感器实验
- 联想YOGA14s(2021) i5-11300H MX450双系统安装Ubuntu20.04全过程
- 【Arcgis】绘制甘肃省行政区划地图
- 由被WebInspect***引发的php header()使用问题
- 用C语言实现高铁客运订票系统