matlab 常用算法大全

Matlab 高级算法程序代码汇总

一、灰色预测模型matlab程序

% renkou1=renkou(:,1);%年末常住人口数

% renkou2=renkou(:,2);%户籍人口

% renkou3=renkou(:,3);%非户籍人口

% shjian=1979:2010;

%以上数据自己给

x0=renkou2';

n=length(x0);

lamda=x0(1:n-1)./x0(2:n)

range=minmax(lamda)

x1=cumsum(x0)

for i=2:n

z(i)=0.5*(x1(i)+x1(i-1));

end

B=[-z(2:n)',ones(n-1,1)];

Y=x0(2:n)';

u=B\Y

x=dsolve('Dx+a*x=b','x(0)=x0');

x=subs(x,{'a','b','x0'},{u(1),u(2),x1(1)});

yuce1=subs(x,'t',[0:n-1]);

digits(6),y=vpa(x) %为提高预测精度，先计算预测值，再显示微分方程的解yuce=[x0(1),diff(yuce1)]

epsilon=x0-yuce %计算残差

delta=abs(epsilon./x0) %计算相对误差

rho=1-(1-0.5*u(1))/(1+0.5*u(1))*lamda %计算级比偏差值

%以深圳人口数据得到预测模型及预测误差相关数据

lamda =

Columns 1 through 8

0.9741 0.9611 0.9419 0.8749 0.9311 0.9093

0.9302 0.9254

Columns 9 through 16

0.9245 0.9278 0.9442 0.9376 0.9127 0.9148

0.9332 0.9477

Columns 17 through 24

0.9592 0.9445 0.9551 0.9562 0.9594 0.9461

0.9469 0.9239

Columns 25 through 31

0.9140 0.9077 0.9243 0.9268 0.9312 0.9446

0.9618

range =

0.8749 0.9741

x1 =

1.0e+003 *

Columns 1 through 8

0.0313 0.0634 0.0967 0.1322 0.1727 0.2162 0.2641 0.3155

Columns 9 through 16

0.3711 0.4313 0.4961 0.5647 0.6380 0.7182 0.8059 0.8999

Columns 17 through 24

0.9990 1.1024 1.2119 1.3265 1.4463 1.5712

1.7033 1.8427

Columns 25 through 32

1.9936

2.1588 2.3407 2.5375 2.7499 2.9780

3.2194 3.4705

u =

-0.0665

31.3737

y =

-472.117+503.377*exp(.664533e-1*t)

yuce =

Columns 1 through 8

31.2600 34.5876 36.9641 39.5040 42.2183 45.1192 48.2194 51.5326

Columns 9 through 16

55.0734 58.8576 62.9017 67.2238 71.8428 76.7792 82.0548 87.6928

Columns 17 through 24

93.7183 100.1578 107.0397 114.3945 122.2547 130.6550 139.6324 149.2267

Columns 25 through 32

159.4802 170.4382 182.1492 194.6649 208.0405 222.3352 237.6121 253.9386

epsilon =

Columns 1 through 8

0 -2.4976 -3.5741 -4.0540 -1.6983 -1.5992

-0.3594 -0.0826

Columns 9 through 16

0.5266 1.2824 1.9183 1.4262 1.3772 3.4408

5.6352

6.2772

Columns 17 through 24

5.4417 3.2222 2.4203 0.2055 -2.4047 -5.7350

-7.5924 -9.7767

Columns 25 through 32

-8.5502 -5.3082 -0.2192 2.1651 4.3395 5.7348

3.8379 -2.9086

delta =

Columns 1 through 8

0 0.0778 0.1070 0.1144 0.0419 0.0367

0.0075 0.0016

Columns 9 through 16

0.0095 0.0213 0.0296 0.0208 0.0188 0.0429

0.0643 0.0668

Columns 17 through 24

0.0549 0.0312 0.0221 0.0018 0.0201 0.0459

0.0575 0.0701

Columns 25 through 32

0.0567 0.0321 0.0012 0.0110 0.0204 0.0251

0.0159 0.0116

rho =

Columns 1 through 8

-0.0411 -0.0271 -0.0066 0.0650 0.0049 0.0282

0.0058 0.0110

Columns 9 through 16

0.0119 0.0084 -0.0091 -0.0020 0.0245 0.0223

0.0027 -0.0128

Columns 17 through 24

-0.0251 -0.0094 -0.0208 -0.0219 -0.0254 -0.0111

-0.0119 0.0126

Columns 25 through 31

0.0232 0.0300 0.0122 0.0095 0.0048 -0.0095

-0.0280

二、遗传算法程序代码

% Optimizing a function using Simple Genetic Algorithm with elitist preserved

%Max f(x1,x2)=100*(x1*x1-x2).^2+(1-x1).^2; -2.0480<=x1,x2<=2.0480 % Author: Wang Yonglin (wylin77@https://www.docsj.com/doc/c512823938.html,)

clc;clear all;

format long;%设定数据显示格式

%初始化参数

T=100;%仿真代数

N=80;% 群体规模

pm=0.05;pc=0.8;%交叉变异概率

umax=2.048;umin=-2.048;%参数取值范围

L=10;%单个参数字串长度，总编码长度2L

bval=round(rand(N,2*L));%初始种群

bestv=-inf;%最优适应度初值

%迭代开始

for ii=1:T

%解码，计算适应度

for i=1:N

y1=0;y2=0;

for j=1:1:L

y1=y1+bval(i,L-j+1)*2^(j-1);

end

x1=(umax-umin)*y1/(2^L-1)+umin;

for j=1:1:L

y2=y2+bval(i,2*L-j+1)*2^(j-1);

end

x2=(umax-umin)*y2/(2^L-1)+umin;

obj(i)=100*(x1*x1-x2).^2+(1-x1).^2; %目标函数

xx(i,:)=[x1,x2];

end

func=obj;%目标函数转换为适应度函数

p=func./sum(func);

q=cumsum(p);%累加

[fmax,indmax]=max(func);%求当代最佳个体

if fmax>=bestv

bestv=fmax;%到目前为止最优适应度值

bvalxx=bval(indmax,:);%到目前为止最佳位串

optxx=xx(indmax,:);%到目前为止最优参数

end

Bfit1(ii)=bestv; % 存储每代的最优适应度

%%%%遗传操作开始

%轮盘赌选择

for i=1:(N-1)

r=rand;

tmp=find(r<=q);

newbval(i,:)=bval(tmp(1),:);

end

newbval(N,:)=bvalxx;%最优保留

bval=newbval;

%单点交叉

for i=1:2:(N-1)

cc=rand;

if cc

point=ceil(rand*(2*L-1));%取得一个1到2L-1的整数ch=bval(i,:);

bval(i,point+1:2*L)=bval(i+1,point+1:2*L);

bval(i+1,point+1:2*L)=ch(1,point+1:2*L);

end

end

bval(N,:)=bvalxx;%最优保留

%位点变异

mm=rand(N,2*L)

mm(N,:)=zeros(1,2*L);%最后一行不变异，强制赋0

bval(mm)=1-bval(mm);

end

%输出

plot(Bfit1);% 绘制最优适应度进化曲线

bestv %输出最优适应度值

optxx %输出最优参数

三、种子群算法程序代码

%declare the parameters of the optimization

max_iterations = 1000;

no_of_particles = 50;

dimensions = 1;

delta_min = -0.003;

delta_max = 0.003;

c1 = 1.3;

c2 = 1.3;

%initialise the particles and teir velocity components

for count_x = 1:no_of_particles

for count_y = 1:dimensions

particle_position(count_x,count_y) = rand*10;

particle_velocity(count_x,count_y) = rand;

p_best(count_x,count_y) = particle_position(count_x,count_y); end

end

%initialize the p_best_fitness array

for count = 1:no_of_particles

p_best_fitness(count) = -1000;

end

%particle_position

%particle_velocity

%main particle swrm routine

for count = 1:max_iterations

%find the fitness of each particle

%change fitness function as per equation requiresd and dimensions

for count_x = 1:no_of_particles

%x = particle_position(count_x,1);

%y = particle_position(count_x,2);

%z = particle_position(count_x,3);

%soln = x^2 - 3*y*x + z;

%x = particle_position(count_x);

%soln = x^2-2*x+1;

x = particle_position(count_x);

soln = x-7;

if soln~=0

current_fitness(count_x) = 1/abs(soln); else

current_fitness =1000;

end

end

%decide on p_best etc for each particle

for count_x = 1:no_of_particles

if current_fitness(count_x) >p_best_fitness(count_x)

p_best_fitness(count_x) = current_fitness(count_x);

for count_y = 1:dimensions

p_best(count_x,count_y) = particle_position(count_x,count_y); end

end

end

%decide on the global best among all the particles

[g_best_val,g_best_index] = max(current_fitness);

%g_best contains the position of teh global best

for count_y = 1:dimensions

g_best(count_y) = particle_position(g_best_index,count_y); end

%update the position and velocity compponents

for count_x = 1:no_of_particles

for count_y = 1:dimensions

p_current(count_y) = particle_position(count_x,count_y);

end

for count_y = 1:dimensions

particle_velocity(count_y) = particle_velocity(count_y) + c1*rand*(p_best(count_y)-p_current(count_y)) +

c2*rand*(g_best(count_y)-p_current(count_y));

particle_positon(count_x,count_y) = p_current(count_y)

+particle_velocity(count_y);

end

end

end

g_best

current_fitness(g_best_index)

clear all, clc % pso example

iter = 1000; % number of algorithm iterations

np = 2; % number of model parameters

ns = 10; % number of sets of model parameters

Wmax = 0.9; % maximum inertial weight

Wmin = 0.4; % minimum inertial weight

c1 = 2.0; % parameter in PSO methodology

c2 = 2.0; % parameter in PSO methodology

Pmax = [10 10]; % maximum model parameter value

Pmin = [-10 -10]; % minimum model parameter value

Vmax = [1 1]; % maximum change in model parameter

Vmin = [-1 -1]; % minimum change in model parameter

modelparameters(1:np,1:ns) = 0; % set all model parameter estimates for all model parameter sets to zero

modelparameterchanges(1:np,1:ns) = 0; % set all change in model parameter estimates for all model parameter sets to zero

bestmodelparameters(1:np,1:ns) = 0; % set best model parameter estimates for all model parameter sets to zero

setbestcostfunction(1:ns) = 1e6; % set best cost function of each model parameter set to a large number

globalbestparameters(1:np) = 0; % set best model parameter values for all model parameter sets to zero

bestparameters = globalbestparameters'; % best model parameter values for all model parameter sets (to plot)

globalbestcostfunction = 1e6; % set best cost function for all model parameter sets to a large number

i = 0; % indicates ith algorithm iteration

j = 0; % indicates jth set of model parameters

k = 0; % indicates kth model parameter

for k = 1:np % initialization

for j = 1:ns

modelparameters(k,j) = (Pmax(k)-Pmin(k))*rand(1) + Pmin(k); % randomly distribute model parameters

modelparameterchanges(k,j) = (Vmax(k)-Vmin(k))*rand(1) + Vmin(k); % randomly distribute change in model parameters

end

end

for i = 2:iter

for j = 1:ns

x = modelparameters(:,j);

% calculate cost function

costfunction = 105*(x(2)-x(1)^2)^2 + (1-x(1))^2;

if costfunction

bestmodelparameters(:,j) = modelparameters(:,j);

setbestcostfunction(j) = costfunction;

end

四、模拟退火算法

% for d=1:50 %循环10次发现最小路径为4.115，循环50次有3次出现4.115 T_max=80; %input('please input the start temprature');

T_min=0.001; %input('please input the end temprature');

iter_max=100;%input('please input the most interp steps on the fit temp');

s_max=100; %input('please input the most steady steps ont the fit temp');

T=T_max;

load .\address.txt;

order1=randperm(size(address,1))';%生成初始解。

figure(1);

plot(address(order1,1),address(order1,2),'*r-')

title('随机产生的路径');

totaldis1=distance(address,order1);

for n=1:size(address,1)

text(address(n,1)+0.01,address(n,2),num2str(n))%标号

end

text(0.9,0.9,num2str(totaldis1))

figure(2);

while T>=T_min

iter_num=1;

s_num=1;

plot(T,totaldis1,'r.')

hold on

while iter_num

order2=exhgpath(order1); %随机交换两个城市位置

totaldis2=distance(address,order2);

R=rand;

DeltaDis=totaldis2-totaldis1; %新的距离-原来的距离

if DeltaDis<0;

order1=order2;

totaldis1=totaldis2; %新距离小，无条件接受

elseif (exp((totaldis1-totaldis2)/(T))>R)%%%%本算法最核心的思想：以一定概率接受坏的结果，防止局部最优

order1=order2;

totaldis1=totaldis2;

else s_num=s_num+1;

end

iter_num=iter_num+1;

end

T=T*0.99;

end

set(gca,'xscale','log');%或者使用semilogx，有相同效果

xlabel('退火温度');ylabel('总距离');

order1

totaldis1

figure(3)

plot(address(order1,1),address(order1,2),'*b-')

title('最终路径');

for n=1:size(address,1)

text(address(n,1)+0.01,address(n,2),num2str(n))%标号

end

text(0.9,0.9,num2str(totaldis1))

dstc(d)=totaldis1;

% end

——————————————————————————————

function y=exhgpath(order)

while 1

b=size(order,1);

r=unidrnd(b,1,2);

if r(1)-r(2)~=0

break

end

end

b=order(r(2));

order(r(2))=order(r(1));

order(r(1))=b;

y=order;

------------------------------------------------------------

function y=distance(address,order)

nmb=size(address,1);

y=0;

for i=1:nmb-1

y=y+sqrt((address(order(i+1),1)-address(order(i),1))^2+(address(order (i+1),2)-address(order(i),2))^2);

end

y=y+sqrt((address(order(i+1),1)-address(order(1),1))^2+(address(order (i+1),2)-address(order(1),2))^2);

人工神经网络程序代码

%产生指定类别的样本点，并在图中绘出

X = [0 1; 0 1]; % 限制类中心的范围

clusters = 5; % 指定类别数目

points = 10; % 指定每一类的点的数目

std_dev = 0.05; % 每一类的标准差

P = nngenc(X,clusters,points,std_dev);

plot(P(1,:),P(2,:),'+r');

title('输入样本向量');

xlabel('p(1)');

ylabel('p(2)');

%建立网络

net=newc([0 1;0 1],5,0.1); %设置神经元数目为5

%得到网络权值，并在图上绘出

figure;

plot(P(1,:),P(2,:),'+r');

w=net.iw{1}

hold on;

plot(w(:,1),w(:,2),'ob');

hold off;

title('输入样本向量及初始权值');

xlabel('p(1)');

ylabel('p(2)');

figure;

plot(P(1,:),P(2,:),'+r');

hold on;

%训练网络

net.trainParam.epochs=7;

net=init(net);

net=train(net,P);

%得到训练后的网络权值，并在图上绘出

w=net.iw{1}

plot(w(:,1),w(:,2),'ob');

hold off;

title('输入样本向量及更新后的权值');

xlabel('p(1)');

ylabel('p(2)');

a=0;

p = [0.6 ;0.8];

a=sim(net,p)

－－－－－－－－－－－－－－－－－－－

example8_2

%随机生成1000个二维向量，作为样本，并绘出其分布P = rands(2,1000);

plot(P(1,:),P(2,:),'+r')

title('初始随机样本点分布');

xlabel('P(1)');

ylabel('P(2)');

%建立网络，得到初始权值

net=newsom([0 1; 0 1],[5 6]);

w1_init=net.iw{1,1}

%绘出初始权值分布图

figure;

plotsom(w1_init,https://www.docsj.com/doc/c512823938.html,yers{1}.distances)

%分别对不同的步长，训练网络，绘出相应的权值分布图for i=10:30:100

net.trainParam.epochs=i;

net=train(net,P);

figure;

plotsom(net.iw{1,1},https://www.docsj.com/doc/c512823938.html,yers{1}.distances)

end

%对于训练好的网络，选择特定的输入向量，得到网络的输出结果p=[0.5;0.3];

a=0;

a = sim(net,p)

－－－－－－－－－－－－－－－－－－－－－－－－－－example8_3

%指定输入二维向量及其类别

P = [-3 -2 -2 0 0 0 0 +2 +2 +3;

0 +1 -1 +2 +1 -1 -2 +1 -1 0];

C = [1 1 1 2 2 2 2 1 1 1];

%将这些类别转换成学习向量量化网络使用的目标向量

T = ind2vec(C)

%用不同的颜色，绘出这些输入向量

plotvec(P,C),

title('输入二维向量');

xlabel('P(1)');

ylabel('P(2)');

%建立网络

net = newlvq(minmax(P),4,[.6 .4],0.1);

%在同一幅图上绘出输入向量及初始权重向量

figure;

plotvec(P,C)

hold on

W1=net.iw{1};

plot(W1(1,1),W1(1,2),'ow')

title('输入以及权重向量');

xlabel('P(1), W(1)');

ylabel('P(2), W(2)');

hold off;

%训练网络，并再次绘出权重向量

figure;

plotvec(P,C);

hold on;

net.trainParam.epochs=150;

net.trainParam.show=Inf;

net=train(net,P,T);

plotvec(net.iw{1}',vec2ind(net.lw{2}),'o');

%对于一个特定的点，得到网络的输出

p = [0.8; 0.3];

a = vec2ind(sim(net,p))

五、蚁群算法程序代码

% the procedure of ant colony algorithm for VRP

%

% % % % % % % % % % %

%initialize the parameters of ant colony algorithms

load data.txt;

d=data(:,2:3);

g=data(:,4);

m=31; % 蚂蚁数

alpha=1;

belta=4;% 决定tao和miu重要性的参数

lmda=0;

rou=0.9;%衰减系数

q0=0.95;

% 概率

tao0=1/(31*841.04);%初始信息素

Q=1;%蚂蚁循环一周所释放的信息素

defined_phrm=15.0; % initial pheromone level value

QV=100; % 车辆容量

vehicle_best=round(sum(g)/QV)+1;%所完成任务所需的最少车数

V=40;

% 计算两点的距离

for i=1:32;

for j=1:32;

dist(i,j)=sqrt((d(i,1)-d(j,1))^2+(d(i,2)-d(j,2))^2);

end;

end;

%给tao miu赋初值

for i=1:32;

for j=1:32;

if i~=j;

%s(i,j)=dist(i,1)+dist(1,j)-dist(i,j);

tao(i,j)=defined_phrm;

miu(i,j)=1/dist(i,j);

end;

end;

end;

for k=1:32;

for k=1:32;

deltao(i,j)=0;

end;

end;

best_cost=10000;

for n_gen=1:50;

print_head(n_gen);

for i=1:m;

%best_solution=[];

print_head2(i);

sumload=0;

cur_pos(i)=1;

rn=randperm(32);

n=1;

nn=1;

part_sol(nn)=1;

%cost(n_gen,i)=0.0;

n_sol=0; % 由蚂蚁产生的路径数量

M_vehicle=500;

t=0; %最佳路径数组的元素数为0

while sumload<=QV;

for k=1:length(rn);

if sumload+g(rn(k))<=QV;

gama(cur_pos(i),rn(k))=(sumload+g(rn(k)))/QV;

A(n)=rn(k);

n=n+1;

end;

end;

fid=fopen('out_customer.txt','a+');

fprintf(fid,'%s %i\t','the current position is:',cur_pos(i));

fprintf(fid,'\n%s','the possible customer set is:')

fprintf(fid,'\t%i\n',A);

fprintf(fid,'------------------------------\n');

fclose(fid);

p=compute_prob(A,cur_pos(i),tao,miu,alpha,belta,gama,lmda,i);

maxp=1e-8;

na=length(A);

for j=1:na;

if p(j)>maxp

maxp=p(j);

index_max=j;

end;

end;

old_pos=cur_pos(i);

if rand(1)

cur_pos(i)=A(index_max);

else

krnd=randperm(na);

cur_pos(i)=A(krnd(1));

bbb=[old_pos cur_pos(i)];

ccc=[1 1];

if bbb==ccc;

cur_pos(i)=A(krnd(2));

end;

end;

tao(old_pos,cur_pos(i))=taolocalupdate(tao(old_pos,cur_pos(i)),rou,tao0);%对所经弧进行局部更新

sumload=sumload+g(cur_pos(i));

nn=nn+1;

part_sol(nn)=cur_pos(i);

temp_load=sumload;

if cur_pos(i)~=1;

rn=setdiff(rn,cur_pos(i));

n=1;

A=[];

end;

if cur_pos(i)==1; % 如果当前点为车场,将当前路径中的已访问用户去掉后,开始产生新路径

if setdiff(part_sol,1)~=[];

n_sol=n_sol+1; % 表示产生的路径数,n_sol=1,2,3,..5,6...,超过5条对其费用加上车辆的派遣费用

fid=fopen('out_solution.txt','a+');

fprintf(fid,'%s%i%s','NO.',n_sol,'条路径是:');

fprintf(fid,'%i ',part_sol);

fprintf(fid,'\n');

fprintf(fid,'%s','当前的用户需求量是:');

fprintf(fid,'%i\n',temp_load);

PID算法Matlab仿真程序和C程序

增量式PID控制算法Matlab仿真程序设一被控对象G（s）=50/(0.125s^2+7s),用增量式PID控制算法编写仿真程序（输入分别为单位阶跃、正弦信号，采样时间为1ms，控制器输出限幅：[-5,5],仿真曲线包括系统输出及误差曲线，并加上注释、图例）。程序如下clear all; close all; ts=0.001; sys=tf(50,[0.125,7, 0]); dsys=c2d(sys,ts,'z'); [num,den]=tfdata(dsys,'v'); u_1=0.0;u_2=0.0; y_1=0.0;y_2=0.0; x=[0,0,0]'; error_1=0; error_2=0; for k=1:1:1000 time(k)=k*ts; S=2; if S==1 kp=10;ki=0.1;kd=15; rin(k)=1; % Step Signal elseif S==2 kp=10;ki=0.1;kd=15; %Sin e Signal rin(k)=0.5*sin(2*pi*k*ts); end du(k)=kp*x(1)+kd*x(2)+ki*x(3); % PID Controller u(k)=u_1+du(k); %Restricting the output of controller if u(k)>=5 u(k)=5; end if u(k)<=-5 u(k)=-5; end %Linear model yout(k)=-den(2)*y_1-den(3)*y_2+nu m(2)*u_1+num(3)*u_2; error(k)=rin(k)-yout(k); %Return of parameters u_2=u_1;u_1=u(k); y_2=y_1;y_1=yout(k); x(1)=error(k)-error_1; %C alculating P x(2)=error(k)-2*error_1+error_2; %Calculating D x(3)=error(k); %Calculating I error_2=error_1; error_1=error(k); end figure(1); plot(time,rin,'b',time,yout,'r'); xlabel('time(s)'),ylabel('rin,yout'); figure(2); plot(time,error,'r') xlabel('time(s)');ylabel('error'); 微分先行PID算法Matlab仿真程序%PID Controler with differential in advance clear all; close all; ts=20; sys=tf([1],[60,1],'inputdelay',80); dsys=c2d(sys,ts,'zoh'); [num,den]=tfdata(dsys,'v'); u_1=0;u_2=0;u_3=0;u_4=0;u_5=0;

matlab常用操作命令

matlab常用操作备忘（1）2007-11-30 22:01:06 分类： 北京理工大学 20981 陈罡 帮助朋友做几个数据的卷积的仿真，一用才知道，呵呵，发现对不住偶的导师了。。。好多matlab的关键字和指令都忘记了。特意收集回顾一下： （1）管理命令和函数 addpath :添加目录到MATLAB搜索路径 doc :在Web浏览器上现实HTML文档 help :显示Matlab命令和M文件的在线帮助 helpwin helpdesk :help 兄弟几个 lookfor :在基于Matlab搜索路径的所有M文件中搜索关键字 partialpath:部分路径名 8*) path :所有关于路径名的处理 pathtool :一个不错的窗口路径处理界面 rmpath :删除搜索路径中指定目录 type :显示指定文件的内容 ver :版本信息 version :版本号 web :打开web页 what :列出当前目录吓所有的M文件 Mat文件和 Mex文件 whatsnew :显示readme文件 which :显示文件位置 （2）管理变量和工作区 clear :从内存中删除所有变量，clear x y z是删除某个变量 disp :显示文本或数组内容 length :数组长度(最长维数) load :重新载入变量(从磁盘上) mlock :锁定文件，防止文件被错误删除 munlock :解锁文件 openvar :在数组编辑器中打开变量 pack :整理内存空间 save :保存变量到文件 8*) size :数组维数

who whos :列出内存变量 workspace :显示工作空间窗口 （3）管理命令控制窗口(command窗口) clc :清空命令窗口 echo :禁止或允许显示执行过程 format :设置输出显示格式 home :光标移动到命令窗口左上角 more :设置命令窗口页输出格式 （4）文件和工作环境 cd :改变工作目录 copyfile :复制文件 delete :删除文件和图形对象 diary :把命令窗口的人机交互保存到文件 dir :显示目录 edit :编辑文本文件 fileparts :返回文件的各个部分 fullfile :使用指定部分建立文件全名 inmem :返回内存(伪代码区)的matlab函数名 ls :在unix系统中列出目录(win中亦可) matlabroot :根目录 mkdir :新建目录 open :打开文件 pwd :显示当前目录 tempdir :返回系统临时目录的名字 tempname :随机给出一个临时字符串(可用作文件名) ! :直接调用操作系统command命令 （5）启动和推出matlab matlabrc :Matlab的启动M文件 exit quit :退出Matlab startup :运行matlab启动文件 （6）程序设计 builtin :从可重载方法中调用内置函数 eval :执行包含可执行表达式的字符串

0.618法的matlab实现

实验报告 实验题目: 0.618法的MATLAB实现学生姓名: 学号: 实验时间: 2013-5-13

一．实验名称: 0.618法求解单峰函数极小点 二．实验目的及要求: 1. 了解并熟悉0.618法的方法原理, 以及它的MATLAB 实现. 2. 运用0.618法解单峰函数的极小点. 三．实验内容： 1. 0.618法方法原理: 定理: 设f 是区间],[b a 上的单峰函数, ] ,[ ,)2()1(b a x x ∈, 且)2()1(x x <. 如果)()()2()1(x f x f >, 则对每一个],[)1(x a x ∈, 有)()()2(x f x f >; 如果)()()2()1(x f x f ≤, 则对每一个] ,[) 2(b x x ∈, 有)()()1(x f x f ≥. 根据上述定理, 只需选择两个试探点, 就可将包含极小点的区间缩短. 事实上, 必有 如果)()()2()1(x f x f >, 则],[)1(b x x ∈; 如果)()() 2()1(x f x f ≤, 则][)2(x a x ，∈. 0.618 法的基本思想是, 根据上述定理, 通过取试探点使包含极小点的区间(不确定区间)不断缩短, 当区间长度小到一定程度时, 区间上各点的函数值均接近极小值, 因此任意一点都可作为极小点的近似. 0.618 法计算试探点的公式: ). (618.0),(382.0k k k k k k k k a b a a b a -+=-+=μλ 2. 0.618法的算法步骤: ①置初始区间],[11b a 及精度要求0>L , 计算试探点1λ和1μ, 计算函数值)(1λf 和)(1μf . 计算公式是 ).(618.0 ),(382.011111111a b a a b a -+=-+=μλ 令1=k . ②若L a b k k <-, 则停止计算. 否则, 当)()(k k f f μλ>时, 转步骤③; 当)()(k k f f μλ≤时, 转步骤④. ③置k k a λ=+1, k k b b =+1, k k μλ=+1,)(618.01111++++-+=k k k k a b a μ, 计算函数值)(1+k f μ, 转步骤⑤.

最短路径算法_matlab程序[1]

算法描述： 输入图G，源点v0，输出源点到各点的最短距离D 中间变量v0保存当前已经处理到的顶点集合，v1保存剩余的集合 1.初始化v1,D 2.计算v0到v1各点的最短距离，保存到D for each i in v0;D(j)=min[D(j),G(v0(1),i)+G(i,j)] ,where j in v1 3.将D中最小的那一项加入到v0，并且从v1删除这一项。 4.转到2，直到v0包含所有顶点。 %dijsk最短路径算法 clear,clc G=[ inf inf 10 inf 30 100; inf inf 5 inf inf inf; inf 5 inf 50 inf inf; inf inf inf inf inf 10; inf inf inf 20 inf 60; inf inf inf inf inf inf; ]; %邻接矩阵 N=size(G,1); %顶点数 v0=1; %源点 v1=ones(1,N); %除去原点后的集合 v1(v0)=0; %计算和源点最近的点 D=G(v0,:); while 1 D2=D; for i=1:N if v1(i)==0 D2(i)=inf; end end D2 [Dmin id]=min(D2); if isinf(Dmin),error,end v0=[v0 id] %将最近的点加入v0集合，并从v1集合中删除 v1(id)=0; if size(v0,2)==N,break;end %计算v0（1）到v1各点的最近距离 fprintf('计算v0（1）到v1各点的最近距离\n');v0,v1 id=0; for j=1:N %计算到j的最近距离 if v1(j)

实验1熟悉matlab环境和基本操作

实验1 熟悉Matlb环境及基本操作 实验目的： 1．熟悉Matlab环境，掌握Matlab的主要窗口及功能； 2．学会Matlab的帮助使用； 3．掌握向量、矩阵的定义、生成方法和基本运算； 4．掌握Matlab的基本符号运算； 5．掌握Matlab中的二维图形的绘制和控制。 实验内容： 1．启动Matlab，说明主窗口、命令窗口、当前目录窗口、工作空间窗口、历史窗口、图形窗口、M文件编辑器窗口的功能。 2．实例操作Matlab的帮助使用。 3．实例操作向量、矩阵的定义、生成方法和基本运算。 4．实例操作Matlab的基本符号运算。 5．实例操作Matlab中的二维图形绘制和控制。 实验仪器与软件： 1．CPU主频在2GHz以上，内存在512Mb以上的PC； 2．Matlab 7及以上版本。 实验讲评： 实验成绩： 评阅教师： 年月日

实验1 熟悉Matlab环境及基本操作 一、Matlab环境及主要窗口的功能 运行Matlab安装目录下的matlab.exe文件可启动Matlab环境，其默认布局如下图： 其中， 1．主窗口的功能是：主窗口不能进行任何计算任务操作，只用来进行一些整体的环境参数设置，它主要对6个下拉菜单的各项和10个按钮逐一解脱。 2．命令窗口的功能是：对MATLAB搜索路径中的每一个M文件的注释区的第一行进行扫描，一旦发现此行中含有所查询的字符串，则将该函数名及第一行注释全部显示在屏幕上。 3. 历史窗口的功能是：历史窗口显示命令窗口中的所有执行过的命令，一方面可以查看曾经执行过的命令，另一方面也可以重复利用原来输入的命令行，可以从命令窗口中直接通过双击某个命令行来执行该命令，

实验一 MATLAB基本操作及运算

实验一 MATLAB 基本操作及运算 一、 实验目的 二、 实验的设备及条件 三、 实验内容 1、 建立以下标量： 1） a=3 2） ，（j 为虚数单位） 3） c=3/2πj e 2、 建立以下向量： 1） Vb= 2.71382882????????-???? 2） Vc=[4 3.8 … -3.8 -4 ] （向量中的数值从4到-4，步长为-0.2） 3、 建立以下矩阵： 1） 3 333Ma ????=?????? Ma 为一个7×7的矩阵，其元素全为3. 2） 11191212921020100Mb ??????=??????

Mb 为一个10×10的矩阵. 3） 114525173238Mc ????=?????? 4、 使用题1中的变量计算下列等式的x,y,z 的值： 1） ((15)/6)111a x e --=+ 2） 2x π= 3） 3ln([()()]sin(/3))x b c b c a π=+-R ，其中R 表示复数实部。 5、 求解函数值22/(2.25)ct y e -=，其中c 取值见题1，t 的取值范围为题2中行 向量Vc 。 6、 使用题1和题3中所产生的标量和矩阵计算等式 1()()T Mx a Mc Mc Mc -=?? 其中*为矩阵所对应行列式的值，参考det 。 7、 函数的使用和矩阵的访问。 1） 计算矩阵Mb 每一列的和，结果应为行向量形式。 2） 计算整个矩阵Mb 的平均值。 3） 用向量[1 1…1] 替换Mb 的最上一行的值 4） 将矩阵Mb 的第2~5行，第3到9列的元素所构成的矩阵赋值给矩阵SubMb 。 5） 删除矩阵Mb 的第一行； 6） 使用函数rand 产生一个1×10的向量r ，并将r 中值小于0.5的元素设置为0。 8、 已知CellA （1, 1）＝‘中国’，CellA （1，2）＝‘北京’，CellA （2，1）是一个3乘3的单位阵，CellA （2, 2）＝[1 2 3]，试用MATLAB 创建一个2×2的细胞数组CellA 。 9、 已知结构数组student 中信息包含有姓名，学号，性别，年龄和班级，试用MATLAB 创建相应的结构数组student 。该数组包含有从自己学号开始连续5个同学的信息（如果学号在你后面的同学不足5个则往前排序），创建完成后查看自己的信息。

最优化方法的Matlab实现(公式(完整版))

第九章最优化方法的MatIab实现 在生活和工作中，人们对于同一个问题往往会提出多个解决方案，并通过各方面的论证从中提取最佳方案。最优化方法就是专门研究如何从多个方案中科学合理地提取出最佳方案的科学。由于优化问题无所不在，目前最优化方法的应用和研究已经深入到了生产和科研的各个领域，如土木工程、机械工程、化学工程、运输调度、生产控制、经济规划、经济管理等，并取得了显著的经济效益和社会效益。 用最优化方法解决最优化问题的技术称为最优化技术，它包含两个方面的内容： 1）建立数学模型即用数学语言来描述最优化问题。模型中的数学关系式反映了最优化问题所要达到的目标和各种约束条件。 2）数学求解数学模型建好以后，选择合理的最优化方法进行求解。 最优化方法的发展很快，现在已经包含有多个分支，如线性规划、整数规划、非线性规划、动态规划、多目标规划等。 9.1 概述 利用Matlab的优化工具箱，可以求解线性规划、非线性规划和多目标规划问题。 具体而言，包括线性、非线性最小化，最大最小化，二次规划，半无限问题，线性、非线性方程（组）的求解，线性、非线性的最小二乘问题。另外，该工具箱还提供了线性、非线性最小化，方程求解，曲线拟合，二次规划等问题中大型课题的求解方法，为优化方法在工程中的实际应用提供了更方便快捷的途径。 9.1.1优化工具箱中的函数 优化工具箱中的函数包括下面几类： 1 ?最小化函数

2.方程求解函数 3.最小—乘（曲线拟合）函数

4?实用函数 5 ?大型方法的演示函数 6.中型方法的演示函数 9.1.3参数设置 利用OPtimSet函数，可以创建和编辑参数结构；利用OPtimget函数，可以获得o PtiOns优化参数。 ? OPtimget 函数 功能：获得OPtiOns优化参数。 语法：

实验一 Matlab基本操作(2016)

实验一 MATLAB 基本操作 一、实验目的 1. 学习和掌握MA TLAB 的基本操作方法 2. 掌握命令窗口的使用 3. 熟悉MATLAB 的数据表示、基本运算 二、实验内容和要求 1. 实验内容 1) 练习MATLAB7.0或以上版本 2) 练习矩阵运算与数组运算 2. 实验要求 1) 每位学生独立完成，交实验报告 2) 禁止玩游戏！ 三、实验主要软件平台 装有MATLAB7.0或以上的PC 机一台 四、实验方法、步骤及结果测试 1. 实验方法：上机练习。 2. 实验步骤： 1) 开启PC ，进入MA TLAB 。 2) 使用帮助命令，查找sqrt 函数的使用方法 答: help sqrt 3) 矩阵、数组运算 a) 已知 ??????????=987654321A ，???? ??????=963852741B ，求)2()(A B B A -?+ 答： A=[1, 2, 3; 4, 5, 6; 7, 8, 9]; B=[1, 4, 7; 2, 5, 8; 3, 6, 9]; (A+B)*(2*B-A) b) 已知?? ????-=33.1x ，??????=π24y ，求T xy ，y x T c) 已知??????????=987654321A ，???? ??????=300020001B ，求A/B, A\B. d) 已知???? ??????=987654321A ，求：(1) A 中第三列前两个元素；(2) A 中所有第二行元素；(3) A 中四个角上的元素；(4) 交换A 的第1、3列。(5) 交换A 的第1、2行。(6) 删除A 的第3列。

e) 已知[]321=x ，[]654=y ，求：y x *.，y x /.，y x \.，y x .^， 2.^x ，x .^2。 f) 给出x=1,2,…,7时，x x sin 的值。 3）常用的数学函数 a ）随机产生一个3x3的矩阵A ，求：(1) A 每一行的最大、最小值，以及最大、最小值所在的列；(2) A 每一列的最大、最小值，以及最大、最小值所在的行；(3) 整个矩阵的最大、最小值；(4) 每行元素之和；(5) 每列元素之和；(6) 每行元素之积；(7) 每列元素之积。 b) 随机产生两个10个元素的向量x ，y 。（1） 求x 的平均值、标准方差。（2） 求x ，y 的相关系数。（3）对x 排序，并记录排序后元素在原向量中的位置。 4) 字符串操作函数 建立一个字符串向量（如‘ABc123d4e56Fg9’），然后对该向量做如下处理： (1) 取第1～5个字符组成的子字符串。 (2) 将字符串倒过来重新排列。 (3) 将字符串中的小写字母变成相应的大写字母，其余字符不变。 (4) 统计字符串中小写字母的个数。

最优化算法-Matlab程序

CG程序代码 function [x,y] = cg(A,b,x0) %%%%%%%%%%%%%%%%%CG算法%%%%%%%%%%%% r0 = A*x0-b; p0 = -r0; k = 0; r = r0; p = p0; x = x0; while r~=0 alpha = -r'*p/(p'*A*p); x = x+alpha*p; rold = r; r = rold+alpha*A*p; beta = r'*r/(rold'*rold); p = -r+beta*p; plot(k,norm(p),'.--'); hold on k = k+1; end y.funcount = k; y.fval = x'*A*x/2-b'*x;

function [x,y] = cg_FR(fun,dfun,x0) %%%%%%%%%%%%%%%CG_FR算法%%%%%%%%%%%%%%% error = 10^-5; f0 = feval(fun,x0); df0 = feval(dfun,x0); p0 = -df0; f = f0; df = df0; p = p0; x = x0; k = 0; while ((norm(df)>error)&&(k<1000)) f = feval(fun,x); [alpha,funcNk,exitflag] = lines(fun,0.01,0.15,0.85,6,f,df'*p,x,p);%%用线搜索找下降距离%% if exitflag == -1 disp('Break!!!'); break; end x = x+alpha*p; dfold = df; df = feval(dfun,x); beta = df'*df/(dfold'*dfold); p = -df+beta*p; plot(k,norm(df),'.--'); hold on k = k+1; end y.funcount = k; y.fval = feval(fun,x); y.error = norm(df);

王能超 计算方法——算法设计及MATLAB实现课后代码

第一章插值方法 1.1Lagrange插值 1.2逐步插值 1.3分段三次Hermite插值 1.4分段三次样条插值 第二章数值积分 2.1 Simpson公式 2.2 变步长梯形法 2.3 Romberg加速算法 2.4 三点Gauss公式 第三章常微分方程德差分方法 3.1 改进的Euler方法 3.2 四阶Runge-Kutta方法 3.3 二阶Adams预报校正系统 3.4 改进的四阶Adams预报校正系统 第四章方程求根 4.1 二分法 4.2 开方法 4.3 Newton下山法 4.4 快速弦截法 第五章线性方程组的迭代法 5.1 Jacobi迭代 5.2 Gauss-Seidel迭代 5.3 超松弛迭代 5.4 对称超松弛迭代 第六章线性方程组的直接法 6.1 追赶法 6.2 Cholesky方法 6.3 矩阵分解方法 6.4 Gauss列主元消去法

第一章插值方法 1.1Lagrange插值 计算Lagrange插值多项式在x=x0处的值. MATLAB文件：（文件名：Lagrange_eval.m）function [y0,N]= Lagrange_eval(X,Y,x0) %X,Y是已知插值点坐标 %x0是插值点 %y0是Lagrange插值多项式在x0处的值 %N是Lagrange插值函数的权系数 m=length(X); N=zeros(m,1); y0=0; for i=1:m N(i)=1; for j=1:m if j~=i; N(i)=N(i)*(x0-X(j))/(X(i)-X(j)); end end y0=y0+Y(i)*N(i); end 用法》X=[…];Y=[…]; 》x0= ; 》[y0,N]= Lagrange_eval(X,Y,x0) 1.2逐步插值 计算逐步插值多项式在x=x0处的值. MATLAB文件：（文件名：Neville_eval.m）function y0=Neville_eval(X,Y,x0) %X,Y是已知插值点坐标 %x0是插值点 %y0是Neville逐步插值多项式在x0处的值 m=length(X); P=zeros(m,1); P1=zeros(m,1); P=Y; for i=1:m P1=P; k=1; for j=i+1:m k=k+1;

MATLAB操作手册常用快捷键大全

常用对象操作：除了一般windows窗口的常用功能键外。 1、!dir 可以查看当前工作目录的文件。!dir& 可以在dos状态下查看。 2、who 可以查看当前工作空间变量名，whos 可以查看变量名细节。 3、功能键： 功能键快捷键说明 方向左键Ctrl+B 光标向后移一个字符 方向右键Ctrl+F 光标向前移一个字符 Ctrl+方向右键Ctrl+R 光标向右移一个字符 Ctrl+方向左键Ctrl+L 光标向左移一个字符 home Ctrl+A 光标移到行首 End Ctrl+E 光标移到行尾 Esc Ctrl+U 清除一行 Del Ctrl+D 清除光标所在的字符 Backspace Ctrl+H 删除光标前一个字符 Ctrl+K 删除到行尾 Ctrl+C 中断正在执行的命令 4、clc可以命令窗口显示的内容，但并不清除工作空间。 二、函数及运算 1、运算符： ＋：加，－：减，*：乘，/：除，\：左除^：幂，‘：复数的共轭转置，（）：制定运算顺序。 2、常用函数表： sin( ) 正弦（变量为弧度） Cot( ) 余切（变量为弧度） sind( ) 正弦（变量为度数） Cotd( ) 余切（变量为度数） asin( ) 反正弦（返回弧度） acot( ) 反余切（返回弧度） Asind( ) 反正弦（返回度数） acotd( ) 反余切（返回度数） cos( ) 余弦（变量为弧度） exp( ) 指数 cosd( ) 余弦（变量为度数） log( ) 对数 acos( ) 余正弦（返回弧度） log10( ) 以10为底对数 acosd( ) 余正弦（返回度数） sqrt( ) 开方 tan( ) 正切（变量为弧度） realsqrt( ) 返回非负根 tand( ) 正切（变量为度数） abs( ) 取绝对值 atan( ) 反正切（返回弧度） angle( ) 返回复数的相位角

MVDR算法matlab程序

clc clear all close all %% 常量定义 Freqs=1.6e9; %工作频率 c=3e8; %光速 lamda=c/Freqs; %波长 d=0.5*lamda; %单元间距 M=16; %天线阵元数 fs=2e6; %采样频率 pd=10; %快拍数 %% 模型建立 %--------------第一个干扰模型-------------------- thetaJ1=20*pi/180; %干扰方向 FreqJ1=5e5; %第一个干扰的频率 J1NR=sqrt(10^(60/10)); J1one=J1NR*exp(j*(2*pi*FreqJ1*(1:1:pd)/fs)); %1*pd %--------------第二个干扰模型-------------------- ThetaJ2=60*pi/180; %干扰方向 FreqJ2=6e5; %第二个干扰的频率 J2NR=sqrt(10^(60/10)); J2one=J2NR*exp(j*(2*pi*FreqJ2*(1:1:pd)/fs)); %1*pd %--------------信号模型-------------------- ThetaS=30*pi/180; FreqS=7e5; SNR=sqrt(10^(40/10)); Sone=SNR*exp(j*(2*pi*FreqS*(1:1:pd)/fs)); %1*pd %--------------空域处理------------------------- I1=zeros(M,1); I2=zeros(M,1); IS=zeros(M,1); for n=1:M I1(n)=exp(j*2*pi*(n-1)*d*sin(thetaJ1)/lamda); I2(n)=exp(j*2*pi*(n-1)*d*sin(ThetaJ2)/lamda); IS(n)=exp(j*2*pi*(n-1)*d*sin(ThetaS)/lamda); end J1=I1*J1one; J1=J1.'; J2=I2*J2one; J2=J2.'; %------------产生噪声------------------------- noise=sqrt(1/2)*randn(pd,M)+j*sqrt(1/2)*randn(pd,M);

实验一Matlab基本操作-新版.pdf

实验一Matlab基本操作 题目： 1．利用基本矩阵产生 3x3 和15x8 的单位阵，全 1 阵，全0 阵，均匀分布的随 机阵（[-1，1]之间），正态分布随机阵（方差4，均值1） 2．利用diag()函数和rot90()产生下列矩阵： 然后求解a 阵的逆矩阵aa 及b 阵的特征值和对应特征向量，并利用reshape 将 aa 阵变换成行向量。 3．产生一均匀分布在(-5,5)随机阵(50x2)，精确到小数点后一位。 4．编程实现当α∈[-π，π]，间隔为1o 时，求解正弦和余弦的值，并利用plot() 函数绘制正弦，余弦曲线。 5．利用rand 函数产生(0,1)间均匀分布的10x10 随机矩阵a,然后统计a 中大于等于0.6 的元素个数。 6．利用randn 函数产生均值为0，方差为 1 的10x10 正态分布随机阵，然后统计其 中大于-0.5，小于0.5 的元素个数。 7．编程实现下表功能： 8．有一矩阵a,找出矩阵中其值大于 1 的元素，并将他们重新排列成列向量b。9．在一保定市区9 月份平均气温变化测量矩阵temp_Baoding_sep 中(48x30)，存在有奇异值（大于42o C，小于0o C）,编程实现删除奇异值所在的行。 10．在给定的100x100 矩阵中，删除整行内容全为0 的行，删除整列内容全为0 的列。 程序： 1. %3X3矩阵 a1=eye(3) a2=ones(3) a3=zeros(3) a4=1-2*rand(3) a5=2*randn(3)+1 %15X8矩阵 b1=eye(15,8) b2=ones(15,8) b3=zeros(15,8) b4=1-2*rand(15,8) b5=2*randn(15,8)+1 运行结果：

[设计]罚函数法MATLAB程序

[设计]罚函数法MATLAB程序 一、进退法、0.618法、Powell法、罚函数法的Matlab程序设计罚函数法(通用) function y=ff(x,k) y=-17.86*0.42*x(1)/(0.8+0.42*x(1))*(1-exp(- 2*(0.8+0.42*x(1))/3))*exp(-1.6)*x(2)-22. 99*x(1)/(0.8+x(1))*(1-exp(-2*(0.8+x(1))/3))*x(3)+k*(x(2)- (1.22*10^2*(9517.8*exp(-1 .6-2*0.42*x(1)/3)*x(2)+19035.6*exp(- 2*x(1)/3)*x(3)))/(1.22*10^2+9517.8*exp(-1.6-2 *0.42*x(1)/3)*x(2)+19035.6*exp(-2*x(1)/3)*x(3)))^2+k*(x(3)-exp(-0.8-2*x(1)/3)*x(3) -exp(-2.4-2*0.42*x(1)/3)*x(2))^2; % 主函数，参数包括未知数的个数n，惩罚因子q，惩罚因子增长系数k，初值x0，以及允许的误差r function G=FHS(x0,q,k,n,r,h,a) l=1; while (l) x=powell(x0,n,q,r(1),h,a); %调用powell函数 g(1)=ff1(x),g(2)=ff2(x) . . . g(p)=ffp(x); %调用不等式约束函数，将其值 %存入数组g h(1)=hh1(x),h(2)=hh2(x) . . . h(t)=hht(x); %调用等式约束函数，将其值%存入数组h for i=1:p

图论算法及matlab程序的三个案例

图论实验三个案例 单源最短路径问题 Dijkstra 算法 Dijkstra 算法是解单源最短路径问题的一个贪心算法。其基本思想是，设置一个顶点集合S 并不断地作贪心选择来扩充这个集合。一个顶点属于集合S 当且仅当从源到该顶点的最短路径长度已知。设v 是图中的一个顶点，记()l v 为顶点 v 到源点v 1的最短距离， ,i j v v V ?∈，若 (,)i j v v E ?，记i v 到j v 的权ij w =∞。 Dijkstra 算法： ① 1{}S v =,1()0l v =；1{}v V v ??-,()l v =∞,1i =,1{}S V v =-； ② S φ=，停止，否则转③； ③ ()min{(),(,)} j l v l v d v v =， j v S ∈，v S ?∈； ④ 存在 1 i v +，使 1()min{()} i l v l v +=，v S ∈； ⑤ 1{} i S S v +=， 1{} i S S v +=-，1i i =+，转②； 实际上，Dijkstra 算法也是最优化原理的应用：如果12 1n n v v v v -是从1v 到 n v 的最短路径，则 12 1 n v v v -也必然是从1v 到 1 n v -的最优路径。 在下面的MATLAB 实现代码中，我们用到了距离矩阵，矩阵第i 行第j 行元 素表示顶点i v 到j v 的权ij w ，若i v 到j v 无边，则realmax ij w =，其中realmax 是 MATLAB 常量，表示最大的实数+308)。 function re=Dijkstra(ma)

龙格库塔方法matlab实现

龙格库塔方法matlab实现~ function ff=rk(yy,x0,y0,h,a,b)%yy为y的导函数,x0,y0,为初值,h为步长,a,b为区间 c=(b-a)/h+1;i1=1; %c为迭代步数；i1为迭代步数累加值 y=y0;z=zeros(c,6); %z生成c行，5列的零矩阵存放结果； %每行存放c次迭代结果，每列分别存放k1~k4及y的结果 for x=a:h:b if i1<=c k1=feval(yy,x,y); k2=feval(yy,x+h/2,y+(h*k1)/2); k3=feval(yy,x+h/2,y+(h*k2)/2); k4=feval(yy,x+h,y+h*k3); y=y+(h/6)*(k1+2*k2+2*k3+k4); z(i1,1)=x;z(i1,2)=k1;z(i1,3)=k2;z(i1,4)=k3;z(i1,5)=k4;z(i1,6)=y; i1=i1+1; end end fprintf(‘结果矩阵，第一列为x(n),第二列~第五列为k1~k4，第六列为y（n+1）的结果') z %在命令框输入下列语句 %yy=inline('x+y'); %>> rk(yy,0,1,0.2,0,1) %将得到结果 %结果矩阵，第一列为x(n),第二列~第五列为k1~k4第六列为y（n+1）的结果 %z = % 0 1.0000 1.2000 1.2200 1.4440 1.2428 % 0.2000 1.4428 1.6871 1.7115 1.9851 1.5836 % 0.4000 1.9836 2.2820 2.3118 2.6460 2.0442 % 0.6000 2.6442 3.0086 3.0451 3.4532 2.6510 % 0.8000 3.4510 3.8961 3.9407 4.4392 3.4365 % 1.0000 4.4365 4.9802 5.0345 5.6434 4.4401

三个遗传算法matlab程序实例

遗传算法程序（一）: 说明: fga.m 为遗传算法的主程序; 采用二进制Gray编码,采用基于轮盘赌法的非线性排名选择, 均匀交叉,变异操作,而且还引入了倒位操作! function [BestPop,Trace]=fga(FUN,LB,UB,eranum,popsize,pCross,pMutation,pInversion,options) % [BestPop,Trace]=fmaxga(FUN,LB,UB,eranum,popsize,pcross,pmutation) % Finds a maximum of a function of several variables. % fmaxga solves problems of the form: % max F(X) subject to: LB <= X <= UB % BestPop - 最优的群体即为最优的染色体群 % Trace - 最佳染色体所对应的目标函数值 % FUN - 目标函数 % LB - 自变量下限 % UB - 自变量上限 % eranum - 种群的代数,取100--1000(默认200) % popsize - 每一代种群的规模；此可取50--200(默认100) % pcross - 交叉概率,一般取0.5--0.85之间较好(默认0.8) % pmutation - 初始变异概率,一般取0.05-0.2之间较好(默认0.1) % pInversion - 倒位概率,一般取0.05－0.3之间较好(默认0.2) % options - 1*2矩阵,options(1)=0二进制编码(默认0),option(1)~=0十进制编 %码,option(2)设定求解精度(默认1e-4) % % ------------------------------------------------------------------------ T1=clock; if nargin<3, error('FMAXGA requires at least three input arguments'); end if nargin==3, eranum=200;popsize=100;pCross=0.8;pMutation=0.1;pInversion=0.15;options=[0 1e-4];end if nargin==4, popsize=100;pCross=0.8;pMutation=0.1;pInversion=0.15;options=[0 1e-4];end if nargin==5, pCross=0.8;pMutation=0.1;pInversion=0.15;options=[0 1e-4];end if nargin==6, pMutation=0.1;pInversion=0.15;options=[0 1e-4];end if nargin==7, pInversion=0.15;options=[0 1e-4];end if find((LB-UB)>0) error('数据输入错误,请重新输入(LB

0计算方法及MATLAB实现简明讲义课件PPS8-1欧拉龙格法

第8章 常微分方程初值问题数值解法 8.1 引言 8.2 欧拉方法 8.3 龙格-库塔方法 8.4 单步法的收敛性与稳定性 8.5 线性多步法

8.1 引 言 考虑一阶常微分方程的初值问题 00(,),[,],(). y f x y x a b y x y '=∈=（1.1） （1.2） 如果存在实数 ，使得 121212(,)(,).,R f x y f x y L y y y y -≤-?∈（1.3） 则称 关于 满足李普希茨(Lipschitz ）条件， 称为 的李普希茨常数（简称Lips.常数）. 0>L f y L f (参阅教材386页)

计算方法及MATLAB 实现 所谓数值解法，就是寻求解 在一系列离散节点 )(x y <<<<<+121n n x x x x 上的近似值 . ,,,,,121+n n y y y y 相邻两个节点的间距 称为步长. n n n x x h -=+1 如不特别说明，总是假定 为定数， ),2,1( ==i h h i 这时节点为 . ) ,2,1,0(0 =+=i nh x x n 初值问题(1.1),(1.2)的数值解法的基本特点是采取 “步进式”. 即求解过程顺着节点排列的次序一步一步地向前推进. 00(,),[,], (). y f x y x a b y x y '=∈=

描述这类算法，只要给出用已知信息 ,,,21--n n n y y y 计算 的递推公式. 1+n y 一类是计算 时只用到前一点的值 ，称为单步法. 1+n y n y 另一类是用到 前面 点的值 ， 1+n y k 11,,,+--k n n n y y y 称为 步法. k 其次，要研究公式的局部截断误差和阶，数值解 与 精确解 的误差估计及收敛性，还有递推公式的计算 稳定性等问题. n y )(n x y 首先对方程 离散化，建立求数值解的递推 公式. ),(y x f y ='

(完整版)纵横向拉开档次法的MATLAB实现

简介：本文档为《纵横向拉开档次法的MATLAB实现》，可适用于工程科技领域，主题内容包含globalxystdszxystdxy定义全局变量loadshuju原始数据xystd=zscore(shuju)数据无量纲处理xystdrow,符等。 global xystdsz xystd x y %定义全局变量 load shuju %原始数据 xystd= zscore (shuju); %数据无量纲处理 [xystdrow,xystdcol]=size(xystd); %----------区域知识创造能力评价---------- for tt=1:xystdcol xystdsz{tt}(:,:)=xystd{tt}(:,1:10); %提取区域知识创造能力指标无量纲值 end [xystdszrow,xystdszcol]=size(xystdsz); [xyrow,xycol]=size(xystdsz{1}); w0=zeros(1,xycol); for i=1:xycol w0(1,i)=1/xycol; % 优化初始值 end Aeq=[]; beq=[]; lb=zeros(1,xycol);ub=ones(1,xycol); %zeros生成零矩阵；ones生成全1阵。 options =optimset('largescale','off'); %优化函数，largescale大规模算法 [w,faval]=fmincon(@YHQU,w0,[],[],Aeq,beq,lb,ub,@fun,options ); %优化求权重；fmincon用来求解非线性多元函数最小值。 wqz1=w./sum(w); %权重归一化 for tt=1:xystdszcol z{tt}(:,1)=xystd{tt}(:,1:10)*wqz1'; % 求评价值 pxacz(:,tt)=px(z{tt}(:,1)) ; % 对评价值排序 end clear w0 w lb ub faval ; clear global xystdsz; %--------区域知识流动能力评价------------ for tt=1:xystdszcol xystdsz{tt}(:,:)=xystd{tt}(:,11:16); %提取区域知识流动能力指标无量纲值 end global xystdsz; [xystdszrow,xystdszcol]=size(xystdsz); [xyrow,xycol]=size(xystdsz{1}); w0=zeros(1,xycol); for i=1:xycol w0(1,i)=1/xycol; % 优化w初始值 end Aeq=[]; beq=[]; lb=zeros(1,xycol);ub=ones(1,xycol); options =optimset('largescale','off'); [w,faval]=fmincon(@YHQU,w0,[],[],Aeq,beq,lb,ub,@fun,options );