常用算法的matlab的源程序

参考目录：

微分方程的数值解见F盘微分方程的数值解.ppt

优化问题111优化问题matlab27.doc

优化工具箱简洁.ppt

辅助优化计算与设计.pdf

检索文献：题目---英语关键词和中文关键词----中文数据库---修正关键词----外文数据库中文关键词----翻译英语关键词

外国网站搜索

###################################################### 曲线拟合及拟合后相应点的值

clear

x=0:0.1:1;

y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];

n=2;

p=polyfit(x , y , n);

xi=linspace(0 , 1 , 100);

z=polyval(p , xi);

plot(x , y , '-.ro ' , xi , z , ' : ')

xlabel('x') ;

ylabel('y=f(x)')

title('Second Order Curve Fitting')

legend('points','line')%图例

###################################################### >> %plot画图的例子

>> x=0:pi/20:2*pi;

>> y1=x.*sin(x);

>> plot(x,y1,'-.r*')%显示每一点的用*，并用虚线连起来

>> hold on %还要在这个加其他的图

>> y2=exp(x/100).*sin(x-pi/2);

>> plot(x,y2,'g') %用实线条画图

>> plot(x,y2,'-mo') %用实线条画图

>> plot(x,y2,'-mo')

>> y3=sin(x-pi);

>> plot(x,y3,'gs')

>> title('y1=..;y2=.;y3=..;')%加标题

>> legend('y1','y2','y3')%加上图利

>> xlabel('x axis biaoshi....');ylabel('y axis biaoshi..')

>> axis([0 2*pi -5 2])%调整坐标轴

>> grid on %加上网格

%要想加上箭头可以手动

###################################################### ###############

1.三维立体网格图mesh

>> clear

>> x=-5:0.01:5;

>> y=-5:0.01:5;

>> z=x.^2+y.^2;

>> [X,Y]=meshgrid(x,y);

>> Z=X.^2+Y.^2;

>> figure,mesh(X,Y,Z)

>> title('san wei surface')

>> xlabel('x axis')

>> ylabel('y axis')

>> zlabel('z axis')

>> axis([-5 5 -5 5 0 100])

grid on %画坐标方格

colorbar%加上彩色比照图

%%%%%%%%%%%%

其中将mesh改为meshz （后加colormap([1 0 0]) axis tight)和waterfall(后家colormap axis tight) 分别可以得到设立一个账莫平台

其中将mesh改为meshc画立体图提取登高线contour提取登高线contour3(Z)立体的登高线

###################################################### #################################

2.三维线和点

clear

t=1:0.01:15;

x=2*(exp(-0.2.*t).*cos(8.*t));

y=2*(exp(-0.4.*t).*cos(8.*t));

plot3(t,x,y)

%标注

>> xlabel('x axis')

>> ylabel('y axis')

>> zlabel('z axis')

>> axis([-5 5 -5 5 0 100])

###################################################### %灰色预测的matlab程序

%x0为原始值（输入六个原始数据）

%f为预测值得到8个数据（后两个为预测数据）精确度很高

%e为残差

function [f,e]=GM(x0)

x1 = cumsum(x0);

B = zeros(5,1);

for i =1:5

B(i) = -(x1(i) + x1(i+1)) /2;

end

for i =1:5

y(i) =x0(i+1);

end

y=y';

B = [B ones(5,1)];

au = (inv(B'*B))*B'*y;

au = au';

a=au(1);

u=au(2);

k = x0(1)-u/a;

t =1:8;

t1=-a*(t-1); t2=-a*(t-2);

f = k*(exp(t1) - exp(t2));

for i = 1:6

e(i)=100*(x0(i)-f(i))/x0(i);

end

###################################################### #################

插值问题预测温度

一维

>> hours=1:1:12;

>> tem=[5 8 9 15 25 29 31 30 22 25 27 24];

>> inh=1:0.1:12;

>> intem=interp1(hours,tem,inh,'spline');

>> plot(hours,tem,'-.ro',inh,intem,'-.g*')

二维

x=1:5;

>> y=1:3;

>> z=[82 81 80 82 84;79 63 61 65 81;84 84 82 85 86];

>> mesh(x,y,z)

>> inx=1:0.2:5;

>> iny=1:0.2:3;

inz=interp2(x,y,z,inx',iny,'cubic');

mesh(inx,iny,inz)

散点插值

具体格式见插值.ppt

gtext('string') 当光标位于一个图形窗口内时，等待用户单击鼠标或键盘。若按下鼠标或键盘，则在光标的位置放置给定的文字“string”

legend 在图形上添加图例。该命令对有多种图形对象类型（线条图，条形图，饼形图等）的窗口中显示一个图例

############################

stem 画二维离散数据的柄形图。该图用线条显示数据点与x轴的距离

x = linspace(0,2,10);

stem(exp(-x.^2),'fill','-.')

stem3 画三维离散数据的柄形图。该图用一线段显示数据离开xy平面的高度

[X,Y,Z] = peaks(20);

stem3(X,Y,Z,’r*')

#################

直方图

>>x = -5:0.1:5;

>>y = randn(1000,1);

>>hist(y,x)

#######################333

feather羽毛图

t=0:0.1:2*pi;

x=sin(t).*t;

feather(t,x)

############################################

pie 饼形图pie3为三维饼型图

x=[1 2 3 4];

explode=[1 1 1 0];

pie(x,explode)

legend('s','k','d','g')

title('')

##################################333

bar条形图

1、行向量

x=[1 2 3 4]; bar(x)

2、矩阵

y=[1 2 3;4 5 6];bar(y,1)

3、函数条形图

x=-2.9:0.2:2.9;

bar(x,exp(x.*sin(x)))

barh二维水平条形图。

可以加上一切修饰

###################################################### ########

fill

fill3用指定的颜色填充三维多边形。

x=0:0.1:5;

>> y=sin(x);

>> z=sin(x).^2;

>> c=1;

>> fill3(x,y,z,c)

神经网络预测

close all

clear

echo on

clc

%BP建模

%原始数据归一化

%就业人口预测4年预测一年从94年开始

m_data=[ 1729.55 1709.26 1674.90 1679.91 1686.71;

1709.26 1674.90 1679.91 1686.71 1688.73;

1674.90 1679.91 1686.71 1688.73 1690.00;

1679.91 1686.71 1688.73 1690.00 1697.00;

1686.71 1688.73 1690.00 1697.00 1710.50;

1688.73 1690.00 1697.00 1710.50 1726.36;

1690.00 1697.00 1710.50 1726.36 1715.62;

1697.00 1710.50 1726.36 1715.62 1700.50];

%定义网络输入p和期望输出t

pause

clc

p1=m_data(:,1:4);

t1=m_data(:,5);

p=p1';t=t1';

[pn,minp,maxp,tn,mint,maxt]=premnmx(p,t)

%设置网络隐单元的神经元数(5~30验证后5个最好）

n=5;

%建立相应的BP网络

pause

clc

net=newff(minmax(pn),[n,1],{'tansig','purelin'},'traingdm');

inputWeights=net.IW{1,1};

inputbias=net.b{1};

layerWeights=net.IW{1,1};

layerbias=net.b{2};

pause

clc

% 训练网络

net.trainParam.show=50;

net.trainParam.lr=0.05;

net.trainParam.mc=0.9;

net.trainParam.epochs=20000;

net.trainParam.goal=1e-3;

pause

clc

%调用TRAINGDM算法训练BP网络

net=train(net,pn,tn);

%对BP网络进行仿真

A=sim(net,pn);

E=A-tn;

M=sse(E)

N=mse(E)

pause

clc

p2=[1690.00 1697.00 1710.50 1726.36;

1697.00 1710.50 1726.36 1715.62;

1710.50 1726.36 1715.62 1700.50];

p2=p2';

p2n=tramnmx(p2,minp,maxp);

a2n=sim(net,p2n);

a2=postmnmx(a2n,mint,maxt)

随机数的产生；

R=rand(m,n) 生成(0,1)上均匀分布的m行n列随机矩阵

P=randperm(n) 生成1,2…,n的一个随机排列

normrnd(mu,sigma,m,n).产生m*‰n阶的正态分布随机数矩阵.

指数分布exprnd

泊松分布poissrnd

Matlab提供了两种除法运算：左除（\）和右除（/）。一般情况下，x=a\b是方程a*x =b的解，而x=b/a是方程x*a=b的解。

1]画图问题：

1对于正余玄等，其实离散的，若将t取得非常小，是连续的，如果t的间隔大，plot(t1,y1,'r.')则明显离散，plot(t1,y1,t1,y1,'r.') 则就会有一条折线将离散的点连接起来了，

2、subplot（a，b，c），将图分为ab快；

grid on 画坐标方格

title （‘加标题’）

xlabel('x'), ylabel('y')

axis([0,pi,-1,1])%控制轴的范围

axis square %使坐标轴长度相同

legend('1','2','3')%图例

text （x ，y ，‘’）在图片上（）

3、 t=(0:pi/50:2*pi);k=0.4:0.1:1;y=cos(t)*k;plot(t,y)这是错误语句； t=(0:pi/50:2*pi)'才对，另外关于幂的问题> Y=sin(th)*sqrt(25-a^2); ??? Error using ==> mpower Matrix must be square.

>> Y=sin(th)*sqrt(25-a.^2); 4、对函数加包络线，如y=sin （x ）[-1，1]

5]数据拟合经典

线形拟合的经典：（详见拟合PPT ）

x=0:0.1:1

y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];

P=polyfit(x,y,2);

xi=linspace(0,1,100);

z=polyval(P,xi);

hold on

plot(x,y,'r.','MarkerSize',20);

plot(xi,z,'b','LineWidth',2),hold off

legend('points','line')%图例

非线性拟合两例

1．0.0.2()kt c t a be =+

function f=nonlinefit(x,t)

f=x(1)+x(2).*exp(0.02*x(3)*t);

t=100:100:1000;

y=1e-03*[4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59];

x0=[0.2 0.05 0.05];

x=lsqcurvefit('nonlinefit',x0,t,y)

f=nonlinefit(x,t)

2．()kt d

c t e v -=

function f=nonlinefit2(x,t)

f=300*exp(-x(1)*t)./x(2);

t=[0.25 0.5 1 1.5 2 3 4 6 8];

c=[19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01];

x0=[0.01 0.02];

x=lsqcurvefit('nonlinefit2',x0,t,c)

f=nonlinefit2(x,t)

插值经典：

hours=1:12;

temps=[5 8 9 15 25 29 31 30 22 25 27 24];

h=1:0.1:12;%（将时间细化，）

t=interp1(hours,temps,h,'spline'); (每个细化后的时间通过插值得到的值，直接输出数据将是很多的)

plot(hours,temps,'+',h,t,hours,temps,'r:') %作图依次将散点图，插值后的图，及散点的折线图画出

xlabel('Hour'),ylabel('Degrees Celsius’)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 6、k=sym('k')% k是一个符号变量;

symsum(S) 对通项S求和,其中k为变量且从0变到k-1.

symsum(S,v) 对通项S求和,指定其中v为变量且v从0变到v-1.

symsum(S,a,b) 对通项S求和,其中k为变量且从a变到b.

symsum(S,v,a,b)对通项S求和,指定其中v为变量且v从a变到b.

7、解微分方程：

解析解

dsolve('fangcheng1 ','fangcheng 2'…..,'chu shi tiao jian ','bian liang' )

数值解

[自变量，函数]=ode45/ode23/ode15s（‘保存的函数名’，[自变量的范围]，[函数的初始值]）

高化底，

8、回归分析：

经典例制：x=0.1:0.01:0.18;x=[x,0.2,0.21,0.23]';

y=[42,41.5,45,45.5,45,47.5,49,55,50,55,55.5,60.5]';

X=[ones(12,1),x];

>> polyfit(x,y,1)

[b,bint,r,rint,stats]=regress(y,X,0.05);

>> b,bint,stats,rcoplot(r,rint)

9、数组的输入：

x=linspace(a,b,n);

x=rand(1,n);//0-1之间的n的随即数

x=randperm(n);//1到n的随即整数

10、将硬盘的数据读入变量

读入txt：a=load(' shu ju lu jing ')

读入excel:a=xlsread(' shu ju lu jing ')

11、矩阵函数

eye(m,n);m*n的矩阵；

magic(n);n为方阵；

ones(m,n);

rand(m,n);是m*n的（0，1）分布随即数矩镇

randn(m,n);m* n矩阵，元素为正态分布的随即数珠；

zeros(m,n);

meshgrid（y,x）产生x-y平面的网格点数据x,y；

det 计算矩镇的行列市，

eig 特征值

inv 逆

rank 秩

trace 迹

norm 摸或范数

orth 求出可将矩镇化为对角阵的正交化

poly 矩镇的特征多项式

12、点乘与矩镇运算：

A.' A'

A.^n 矩镇的每个元素自乘n方A^n矩镇自乘n方

A.^p 矩镇的每个元素自求p次方A^p矩镇的p次幂

A.*B A*B

A./B A/B

B.\A B\A

13、（具体见统计工具）

统计工具箱有关命令：

1）数据的输入，保存，调用；data=[];save data;load data; 2）基本统计量:

均值mean(x),

中位数median(x);

标准差std(x);

方差var(x);

偏度skewness(x);

峰度kurtosis(x);

偏度：反映分布的对称性，〉0表示右偏态，此时数据位于均值得右边比位于左边的多，峰度：如果

3，常见的分布函数：

常见的几种分布的命令字符为：

正态分布：norm 指数分布：exp

帕松分布：poiss β分布：beta

威布尔分布：weib 2

χ分布：chi2

t分布：t F分布：F

Matlab工具箱对每一种分布都提供五类函数，其命令字符为：

概率密度：pdf 概率分布：cdf

逆概率分布：inv 均值与方差：stat

随机数生成：rnd

（当需要一种分布的某一类函数时，将以上所列的分布命令字符与函数命令字符接起来，并输入自变量（可以是标量、数组或矩阵）和参数即可.）

4、

统计图形的作图

normplot（）正态分布检验图；

hist（y,x）做直方图；

chi2rnd（5，1，100）产生自由度为5的1000个卡方随即数；

pie（X，explode）//X一为数组，作图时将元素求和，然后计算出每个元素的的比绿，explode：数组X 的同维数组，元素为1表示响应的突出

scatter（nd，pop，‘rs’）；

scatter3（X，Y，Z）立体散点图

MA TLAB中灰度直方图的显示

MATLAB图象处理工具箱提供了imhist函数来计算和显示图象的直方图，imhist函数的语法格式为：

imhist(I,n)

imhist(X,map)

其中imhist(I,n)计算和显示灰度图象I的直方图，n为指定的灰度级数目，默认值为256。imhist(X,map)计算和显示索引色图象X的直方图，map为调色板。

例：

I = imread('rice.tif');

imshow(I)

figure, imhist(I)

5、参数估计：

1、正态总体的参数估计

设总体服从正态分布，则其点估计和区间估计可同时由以下命令获得：

[muhat,sigmahat,muci,sigmaci] = normfit(X,alpha)

此命令在显著性水平alpha下估计数据X的参数（alpha缺省时设定为0.05），返回值muhat是X的均值的点估计值，sigmahat是标准差的点估计值, muci是均值的区间估计,sigmaci是标准差的区间估计.

（1）[muhat, muci] = expfit(X,alpha)----- 在显著性水平alpha下，求指数分布的数据X的均值的点估计及其区间估计.

（2）[lambdahat, lambdaci] = poissfit(X,alpha)----- 在显著性水平alpha

下，求泊松分布的数据X 的参数的点估计及其区间估计.

（3）[phat, pci] = weibfit(X,alpha)----- 在显著性水平alpha下，求Weibull分布的数据X 的参数的点估计及其区间估计.

6、假设检验：

在总体服从正态分布的情况下，可用以下命令进行假设检验、

总体方差sigma2已知时，总体均值的检验使用z-检验

[h,sig,ci] = ztest(x,m,sigma,alpha,tail)

检验数据x 的关于均值的某一假设是否成立，其中sigma 为已知方差，alpha 为显著性水平，究竟检验什么假设取决于tail 的取值：tail = 0，检验假设“x 的均值等于m ”

tail = 1，检验假设“x 的均值大于m ”

tail =-1，检验假设“x 的均值小于m ”

tail的缺省值为0，alpha的缺省值为0.05.

返回值h 为一个布尔值，h=1 表示可以拒绝假设，h=0 表示不可以拒绝假设，sig 为假设成立的概率，ci 为均值的1-alpha 置信区间. 总体方差sigma2未知时，总体均值的检验使用t-检验

两总体均值的假设检验使用t-检验

非参数检验：总体分布的检验

7、Matlab工具箱提供了两个对总体分布进行检验的命令:

1）h = normplot(x)

此命令显示数据矩阵x的正态概率图.如果数据来自于正态

分布，则图形显示出直线性形态.而其它概率分布函数显示出曲线形态.

h = weibplot(x)

14、

两自变量作图：

ezsurf（f）

ezsurf（f，domain）

ezsurf（x，y，z）

ezsurf（x，y，z，[smin，smax，tmin，tmax]）

ezsurf（。。。，n）；

ezsurf（。。。，‘circ’）

domain两自变量的定义域。x，y，z：参数方程

[smin，smax，tmin，tmax]自定义绘图区域，

n 网格点的数，circ’作图定义于为圆

15、优化问题：(具体解释可以参考1、优化问题)

线形规划：linprog

非线性规划问题：有约束的一元函数的最小值：fminbnd

无约束多元函数最小值:fminsearch fminunc（阶数达时）

有约束的多元函数最小值：fmincon

function y=objfun(x);

M文件：

y=exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);

%编写返回约束之的函数

function [c ceq]=confun(x)

%非线性不等式约束

c=[1.5+x(1)*x(2)-x(1)-x(2):-x(1)*x(2)-10];

%线形等式约束、

ceq=[];

x0=[-1,1];

options=optimset('largescale','off');

>> [x.fval]=fmincon('objfun',x0,[],[],[],[],[],[],'confun',options)

二次规划问题quadprog

“半无限”有约束的多元函数最优解fseminf

极小化极大（Minmax）问题

多目标规划问题

最小二乘最优问题: 约束线性最小二乘lsqlin()

非线性数据（曲线）拟合: 在MA TLAB5.x中，使用函数curvefit解决这类问题。现在用函数lsqcurvefit

非线性最小二乘lsqnonlin

非负线性最小二乘lsqnonneg

非线性方程(组)求解fzero

非线性方程组的解fsolve

16、图论工具箱的使用说明：下面是从网上搜索到的图论工具箱，但有些英文确实不能理解，请教各位牛人

GrTheory - Graph Theory Toolbox.

Functions:

一些例子 grTheoryTest - test program for all functions.

grBase - find all bases of digraph;

grCoBase - find all contrabases of digraph;

寻找独立的割集grCoCycleBasis - find all independent cut-sets for a connected graph;

边着色问题 grColEdge - solve the color problem for graph edges;

点着色问题grColVer - solve the color problem for graph vertexes;

寻找所有的连通分枝grComp - find all components of graph;

连通图的独立圈grCycleBasis - find all independent cycles for a connected graph;

有向图分解为强连通分枝grDecOrd - solve the problem about decomposition（分解） of the digraph

to the sections with mutually accessed vertexes (strongly connected components强连通的分枝); 最短路问题grDistances - find the distances between any vertexes of graph;

grEccentricity - find the (weighted) eccentricity（离心率） of all vertexes, radius, diameter, center vertexes and the periphery（外围） vertexes;

寻找欧拉图grIsEulerian - find the Eulerian cycle of graph;

寻找最大完全子图grMaxComSu - solve the maximal complete sugraph problem for the graph;

网络中的最大流grMaxFlows - solve the maximal flow problem for the digraph;

寻找最大匹配grMaxMatch - solve the maximal matching problem for the graph;

grMaxStabSet - solve the maximal stable set（稳定集） problem for the graph; grMinAbsEdgeSet - solve the minimal absorbant set problem for the graph edges; grMinAbsVerSet - solve the minimal absorbant set problem for the graph vertexes;

寻找最小割集grMinCutSet - solve the minimal cut-set problem for the digraph;

寻找最小边覆盖grMinEdgeCover - solve the minimal edge cover problem for the graph;

寻找最小生成树grMinSpanTree - solve the minimal spanning tree problem for the graph;

寻找最小点覆盖grMinVerCover - solve the minimal vertex cover problem for the graph;

grPERT - solve the project evaluation research task;

画图（有向或无向）grPlot - draw the plot of the graph (digraph);

寻找有向图的最短路问题grShortPath - solve the shortest path problem for the digraph;

寻找不对称的旅行者问题（TSP问题）grTravSale - solve the nonsymmetrical traveling salesman problem; grValidation - auxiliary function (the data validation);

还可以参考2、图论算法

17、回归分析：（详细可见回归工具）

经典例之：x=[143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164]';

>> X=[ones(16,1) x];

>> Y=[88 85 88 91 92 93 93 95 96 98 97 96 98 99 100 102]'; >> [b,bint,r,rint,stats]=regress(X,Y)

??? Error using ==> regress

Y must be a vector.

>> [b,bint,r,rint,stats]=regress(Y,X)

>> rcoplot(r,rint)%%%可以看出第2个可以删取

>> X(2,:)=[];Y(2)=[];

>> [b,bint,r,rint,stats]=regress(Y,X)

多元二项回归：有两种方法：

1）>> x1=[1000 600 1200 500 300 400 1300 1100 1300 300]; x2=[5 7 6 6 8 7 5 4 3 9];

y=[100 75 80 70 50 65 90 100 110 60]';

x=[x1' x2'];

>> rstool(x,y,'purequadratic')

2）略

非线性拟合包括函数拟合

逐步回归stepwise（x，y，inmodel，alpha）例：

>> x1=[4.46 4.11 3.58 4.42 4.62]';

>> x2=[4.42 3.82 3.31 4.37 4.47]';

>> x3=[4.23 3.29 3.24 4.34 4.53]';

>> X=[x1 x2 x3];

>> y=[4.41 3.38 3.17 4.39 4.69]'; >> stepwise(X,y)

区分回归与拟合！！

18

二、网格法（半穷举法）

算法略。（类似下面的代码格式）

网格法实现代码参考：

s1=0.01;%变量x1搜索步长

s2=0.01;%变量x2搜索步长

opt_value =inf;%预置一个足够大的数

%二重嵌套循环进行“穷举”

for x1=-5:s1:5,

for x2=0:s2:10,

v = 100*(x2-x1^2)^2+(1-x1)^2;

if v< opt_value,

opt_value = v;

opt_x = [x1,x2];

end?????

end

end

%显示结果

opt_value

opt_x

计算机模拟：

追逐模拟：

v=1;

dt=0.05;

x=[0 0 10 10];

y=[0 10 10 0];

d=20;

while (d>0.1)

x(5)=x(1);y(5)=y(1);

for i=1:4

d=sqrt((x(i+1)-x(i))^2+(y(i+1)-y(i))^2);

x(i)=x(i)+v*dt*(x(i+1)-x(i))/d;

y(i)=y(i)+v*dt*(y(i+1)-y(i))/d;

plot(x(i),y(i),'.'),hold on

end

end

数值函数

N[expr]表达式的机器精度近似值

N[expr, n] 表达式的n位近似值，n为任意正整数 NSolve[lhs==rhs, var] 求方程数值解

NSolve[eqn, var, n] 求方程数值解，结果精度到n位

NDSolve[eqns, y, {x, xmin, xmax}]微分方程数值解

NDSolve[eqns, {y1,y2,...}, {x, xmin, xmax}]

微分方程组数值解

FindRoot[lhs==rhs, {x,x0}] 以x0为初值，寻找方程数值解 FindRoot[lhs==rhs, {x, xstart, xmin, xmax}]

NSum[f, {i,imin,imax,di}] 数值求和，di为步长

NSum[f, {i,imin,imax,di}, {j,..},..] 多维函数求和

NProduct[f, {i, imin, imax, di}]函数求积

NIntegrate[f, {x, xmin, xmax}] 函数数值积分

优化函数：

FindMinimum[f, {x,x0}] 以x0为初值，寻找函数最小值

FindMinimum[f, {x, xstart, xmin, xmax}]

ConstrainedMin[f,{inequ},{x,y,..}]

inequ为线性不等式组，f为x,y..之线性函数，得到最小值及此时的x,y..取值

ConstrainedMax[f, {inequ}, {x, y,..}]同上

LinearProgramming[c,m,b] 解线性组合c.x在m.x>=b&&x>=0约束下的

最小值，x,b,c为向量,m为矩阵

LatticeReduce[{v1,v2...}] 向量组vi的极小无关组

数据处理：

Fit[data,funs,vars]用指定函数组对数据进行最小二乘拟和

data可以为{{x1,y1,..f1},{x2,y2,..f2}..}多维的情况

emp: Fit[{10.22,12,3.2,9.9}, {1, x, x^2,Sin[x]}, x]

Interpolation[data]对数据进行差值,

data同上，另外还可以为{{x1,{f1,df11,df12}},{x2,{f2,.}..}指定各阶导数

InterpolationOrder默认为3次，可修改

ListInterpolation[array]对离散数据插值，array可为n维

ListInterpolation[array,{{xmin,xmax},{ymin,ymax},..}]

FunctionInterpolation[expr,{x,xmin,xmax}, {y,ymin,ymax},..] 以对应expr[xi,yi]的为数据进行插值

Fourier[list] 对复数数据进行付氏变换

InverseFourier[list] 对复数数据进行付氏逆变换

Min[{x1,x2...},{y1,y2,...}]得到每个表中的最小值

Max[{x1,x2...},{y1,y2,...}]得到每个表中的最大值

Select[list, crit] 将表中使得crit为True的元素选择出来 Count[list, pattern] 将表中匹配模式pattern的元素的个数

Sort[list] 将表中元素按升序排列

Sort[list,p] 将表中元素按p[e1,e2]为True的顺序比较list 的任两个元素e1,e2,实际上Sort[list]中默认p=Greater

集合论：

Union[list1,list2..] 表listi的并集并排序

Intersection[list1,list2..] 表listi的交集并排序

Complement[listall,list1,list2...]从全集listall中对listi的差集

九、虚数函数

Re[expr] 复数表达式的实部

Im[expr] 复数表达式的虚部

Abs[expr] 复数表达式的模

Arg[expr] 复数表达式的辐角

Conjugate[expr] 复数表达式的共轭

王能超 计算方法——算法设计及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;

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实现结构可靠度计算 口徐华…朝泽刚‘u刘勇‘21 。 (【l】中国地质大学(武汉工程学院湖北?武汉430074; 12】河海大学土木工程学院江苏?南京210098 摘要:Matlab提供了各种矩阵的运算和操作,其中包含结构可靠度计算中常用的各种数值计算方法工具箱,本文从基本原理和相关算例分析两方面,阐述利用Matlab,编制了计算结构可靠度Matlab程.序,使得Matlab-语言在可靠度计算中得到应用。 关键词:结构可靠度Matlab软件最优化法 中图分类号:TP39文献标识码:A文章编号:1007-3973(200902-095-Ol 1结构可靠度的计算方法 当川概率描述结构的可靠性时,计算结构可靠度就是计算结构在规定时问内、规定条件F结构能够完成预定功能的概率。 从简单到复杂或精确稃度的不同,先后提出的可靠度计算方法有一次二阶矩方法、二次二阶矩方法、蒙特卡洛方法以及其他方法。一次■阶矩方法又分为。I-心点法和验算点法,其中验算点法足H前可靠度分析最常川的方法。 2最优化方法计算可靠度指标数学模型 由结构111n个任意分布的独立随机变量一,x:…以表示的结构极限状态方程为:Z=g(■.托…t=0,采用R-F将非正念变量当罱正态化,得到等效正态分布的均值o:和标准差虹及可靠度指标B,由可靠度指标B的几何意义知。o;辟

开始时验算点未知,把6看成极限状态曲面上点P(■,爿:---37,的函数,通过优化求解,找到B最小值。求解可靠皮指标aJ以归结为以下约束优化模型: rain睁喜t华,2 s.,.Z=g(工i,x2’,…,工:=0 如极限状态方栉巾某个变最(X。可用其他变量表示,则上述模型jfIJ‘转化为无约束优化模型: 。。B!:手f生丛r+阻:坚:坠:盐尘}二剐 t∞oY?’【叫,J 3用MATLAB实现结构可靠度计算 3.1Matlab简介 Matlab是++种功能强、效率高、便.丁.进行科学和工程计算的交互式软件包,汇集了人量数学、统计、科学和工程所需的函数,MATI.AB具有编程简甲直观、用户界mf友善、开放性强等特点。将MATLAB用于蒙特卡罗法的一个显著优点是它拥有功能强大的随机数发生器指令。 3.2算例 3.2.I例:已知非线形极限状态方程z=g(t r'H=567f r-0.5H2=0’f、r服从正态分布。IIf=0.6,o r=0.0786;la|_ 2.18,o r_0.0654;H服从对数正态分布。u H= 3218,O。 =0.984。f、r、H相互独立,求可靠度指标B及验算点(,,r’,H‘。 解:先将H当量正念化:h=ln H服从正态分布,且 ,‘-““了:等专虿’=,。49?口二-、『五ir面_。。3

计算方法及其MATLAB实现第二章作业

作者：夏云木子 1、 >> syms re(x) re(y) re(z) >> input('计算相对误差：'),re(x)=10/1991,re(y)=0.0001/1.991,re(y)=0.0000001/0.0001991 所以可知re(y)最小，即y精度最高 2、 >> format short,A=sqrt(2) >> format short e,B=sqrt(2) >> format short g,C=sqrt(2)

>> format long,D=sqrt(2) >> format long e,E=sqrt(2) >> format long g,F=sqrt(2) >> format bank,H=sqrt(2) >> format hex,I=sqrt(2) >> format +,J=sqrt(2) >> format,K=sqrt(2)

3、 >> syms A >> A=[sqrt(3) exp(7);sin(5) log(4)];vpa(pi*A,6) 4、1/6251-1/6252=1/6251*6252 5、（1）1/（1+3x）-(1-x)/(1+x)=x*(3*x-1)/[(1+3*x)*(1+x)] (2) sqrt(x+1/x)-sqrt(x-1/x)=2/x/[sqrt(x-1/x)+sqrt(x+1/x)] (3) log10(x1)-log(x2)=log10(x1/x2) (4) [1-cos(2*x)]/x =x^2/factorial(2)-x^4/factorial(4)+x^6/factorial(6)-…

matlab用于计算方法的源程序

1、Newdon迭代法求解非线性方程 function [x k t]=NewdonToEquation(f,df,x0,eps) %牛顿迭代法解线性方程 %[x k t]=NewdonToEquation(f,df,x0,eps) %x:近似解 %k:迭代次数 %t:运算时间 %f:原函数，定义为内联函数 ?:函数的倒数，定义为内联函数 %x0:初始值 %eps:误差限 % %应用举例： %f=inline('x^3+4*x^2-10'); ?=inline('3*x^2+8*x'); %x=NewdonToEquation(f,df,1,0.5e-6) %[x k]=NewdonToEquation(f,df,1,0.5e-6) %[x k t]=NewdonToEquation(f,df,1,0.5e-6) %函数的最后一个参数也可以不写。默认情况下，eps=0.5e-6 %[x k t]=NewdonToEquation(f,df,1) if nargin==3 eps="0".5e-6; end tic; k=0; while 1 x="x0-f"(x0)./df(x0); k="k"+1; if abs(x-x0) < eps || k >30 break; end x0=x; end t=toc; if k >= 30 disp('迭代次数太多。'); x="0"; t="0"; end

2、Newdon迭代法求解非线性方程组 function y="NewdonF"(x) %牛顿迭代法解非线性方程组的测试函数 %定义是必须定义为列向量 y(1,1)=x(1).^2-10*x(1)+x(2).^2+8; y(2,1)=x(1).*x(2).^2+x(1)-10*x(2)+8; return; function y="NewdonDF"(x) %牛顿迭代法解非线性方程组的测试函数的导数 y(1,1)=2*x(1)-10; y(1,2)=2*x(2); y(2,1)=x(2).^+1; y(2,2)=2*x(1).*x(2)-10; return; 以上两个函数仅供下面程序的测试 function [x k t]=NewdonToEquations(f,df,x0,eps) %牛顿迭代法解非线性方程组 %[x k t]=NewdonToEquations(f,df,x0,eps) %x:近似解 %k:迭代次数 %t:运算时间 %f:方程组（事先定义） ?:方程组的导数（事先定义） %x0:初始值 %eps:误差限 % %说明：由于虚参f和df的类型都是函数，使用前需要事先在当前目录下采用函数M文件定义% 另外在使用此函数求解非线性方程组时，需要在函数名前加符号“@”，如下所示 % %应用举例： %x0=[0,0];eps=0.5e-6; %x=NewdonToEquations(@NewdonF,@NewdonDF,x0,eps) %[x k]=NewdonToEquations(@NewdonF,@NewdonDF,x0,eps) %[x k t]=NewdonToEquations(@NewdonF,@NewdonDF,x0,eps) %函数的最后一个参数也可以不写。默认情况下，eps=0.5e-6 %[x k t]=NewdonToEquations(@NewdonF,@NewdonDF,x0,eps)

层次分析法计算权重在matlab中的实现

信息系统分析与设计作业 层次分析法确定绩效评价权重在matlab中的实现 小组成员：孙高茹、王靖、李春梅、郭荣1 程序简要概述 编写程序一步实现评价指标特征值lam、特征向量w以及一致性比率CR的求解。 具体的操作步骤是：首先构造评价指标，用专家评定法对指标两两打分，构建比较矩阵，继而运用编写程序实现层次分析法在MATLAB中的应用。 通过编写MATLAB程序一步实现问题求解，可以简化权重计算方法与步骤，减少工作量，从而提高人力资源管理中绩效考核的科学化电算化。 2 程序在matlab中实现的具体步骤 function [w,lam,CR] = ccfx(A) %A为成对比较矩阵，返回值w为近似特征向量 % lam为近似最大特征值λmax，CR为一致性比率 n=length(A(:,1)); a=sum(A); B=A %用B代替A做计算 for j=1:n %将A的列向量归一化 B(:,j)=B(:,j)./a(j); end s=B(:,1); for j=2:n s=s+B(:,j); end c=sum(s);%计算近似最大特征值λmax w=s./c; d=A*w lam=1/n*sum((d./w)); CI=(lam-n)/(n-1);%一致性指标 RI=[0,0,0.58,0.90,1.12,1.24,1.32,1.41,1.45,1.49,1.51];%RI为随机一致

性指标 CR=CI/RI(n);%求一致性比率 if CR>0.1 disp('没有通过一致性检验'); else disp('通过一致性检验'); end end 3 案例应用 我们拟构建公司员工绩效评价分析权重，完整操作步骤如下： 3.1构建的评价指标体系 我们将影响员工绩效评定的指标因素分为：打卡、业绩、创新、态度与品德。 3.2专家打分，构建两两比较矩阵 A = 1.0000 0.5000 3.0000 4.0000 2.0000 1.0000 5.0000 3.0000 0.3333 0.2000 1.0000 2.0000 0.2500 0.3333 0.5000 1.0000 3.3在MATLAB中运用编写好的程序实现 直接在MATLAB命令窗口中输入 [w,lam,CR]=ccfx(A) 继而直接得出 d = 1.3035 2.0000 0.5145 0.3926 w = 0.3102 0.4691 0.1242 0.0966 lam =4.1687

计算方法及其MATLAB实现第一章作业

计算方法作业（作者：夏云木子） 1、help linspace type linspace 2、a1=[5 12 47;13 41 2;9 6 71];a2=[12 9;6 15;7 21];B=a1*a2, C=a1(:,1:2).*a2, D=a1.^2,

E=a1(:).^2 3、a1=[5 12 47;13 41 2;9 6 71];a2=[12 9;6 15;7 21];a1(4:5,1:3)=a2.';a1([4 5],:)=a1([5 4],:);b1=a1

c1=b1(4,1),c2=b1(5,3),D=b1(3:4,:)*a2 4、a1=[5 12 47;13 41 2;9 6 71]; E=eye(3,3); S = a1 + 5*a1' - E, S1=a1^3-rot90(a1)^2+6*E 5、a1=[5 12 47;13 41 2;9 6 71];s=5;A=s-a1,B=s*a1,C=s.*a1,D=s./a1,E=a1./s

6、c=[1 2 3 4;5 6 7 8;9 10 11 12;13 14 15 16];A=c^-4,B=(c^3)^-1,C=(3*c+5*c^-1)/5

7、a=[1 i 3;9i 2-i 8;7 4 8+i];A=a.' 8、abc=[-2.57 8.87;-0.57 3.2-5.5i];m1=sign(abc),m2=round(abc),m3=floor(abc) Sign为符号函数，round表示四舍五入取整，floor表示舍去小数部分取整

9、x=[1 4 3 2 0 8 10 5]';y=[8 0 0 4 2 1 9 11]';A=dot(x,y) 10、a=[3.82 5.71 9.62];b=[7.31 6.42 2.48];A=dot(a,b),B=cross(a,b) 11、P=[5 7 8 0 1];Pf=poly(P);Px=poly2str(Pf,'x') 12、P=[3 0 9 60 0 -90];K1=polyval(P,45),K2=polyval(P,-123),K3=polyval(P,579) 13、P1=[13 55 0 -17 9];P2=[63 0 26 -85 0 105];PP=conv(P1,P2);P1P2=poly2str(PP,'x'),[Q,r]=deconv(P2,P1)

数值计算方法matlab程序

function [x0,k]=bisect1(fun1,a,b,ep) if nargin<4 ep=1e-5; end fa=feval(fun1,a); fb=feval(fun1,b); if fa*fb>0 x0=[fa,fb]; k=0; return; end k=1; while abs(b-a)/2>ep x=(a+b)/2; fx=feval(fun1,x); if fx*fa<0 b=x; fb=fx; else a=x; fa=fx; k=k+1; end end x0=(a+b)/2; >> fun1=inline('x^3-x-1'); >> [x0,k]=bisect1(fun1,1.3,1.4,1e-4) x0 = 1.3247 k = 7 >> 简单迭代法 function [x0,k]=iterate1(fun1,x0,ep,N) if nargin<4 N=500; end if nargin<3 ep=1e-5; end x=x0; x0=x+2*ep;

while abs(x-x0)>ep & k> fun1=inline('(x+1)^(1/3)'); >> [x0,k]=iterate1(fun1,1.5) x0 = 1.3247 k = 7 >> fun1=inline('x^3-1'); >> [x0,k]=iterate1(fun1,1.5) x0 = Inf k = 9 >> Steffesen加速迭代（简单迭代法的加速）function [x0,k]=steffesen1(fun1,x0,ep,N) if nargin<4 N=500; end if nargin<3 ep=1e-5; end x=x0; x0=x+2*ep; k=0; while abs(x-x0)>ep & k

(完整word版)自己编写算法的功率谱密度的三种matlab实现方法

功率谱密度的三种matlab 实现方法 一：实验目的： (1)掌握三种算法的概念、应用及特点； (2)了解谱估计在信号分析中的作用； (3) 能够利用burg 法对信号作谱估计，对信号的特点加以分析。 二;实验内容： （1）简单说明三种方法的原理。 （2）用三种方法编写程序，在matlab 中实现。 （3）将计算结果表示成图形的形式，给出三种情况的功率谱图。 （4）比较三种方法的特性。 （5）写出自己的心得体会。 三：实验原理： 1.周期图法： 周期图法又称直接法。它是从随机信号x(n)中截取N 长的一段，把它视为能量有限x(n)真实功率谱)(jw x e S 的估计)(jw x e S 的抽样. 认为随机序列是广义平稳且各态遍历的，可以用其一个样本x(n)中的一段)(n x N 来估计该随机序列的功率谱。这当然必然带来误差。由于对)(n x N 采用DFT ，就默认)(n x N 在时域是周期的，以及)(k x N 在频域是周期的。这种方法把随机序列样本x(n)看成是截得一段)(n x N 的周期延拓，这也就是周期图法这个名字的来历。

2.相关法（间接法): 这种方法以相关函数为媒介来计算功率谱，所以又叫间接法。这种方法的具体步骤是： 第一步：从无限长随机序列x(n)中截取长度N 的有限长序列列 )(n x N 第二步：由N 长序列)(n x N 求（2M-1）点的自相关函数)(m R x ∧ 序列。 )()(1)(1 m n x n x N m R N n N N x += ∑-=∧ (2-1) 这里，m=-(M-1)…,-1,0,1…,M-1，M N ，)(m R x 是双边序列，但是由自相关函数的偶对称性式，只要求出m=0,。。。,M-1的傅里叶变换，另一半也就知道了。 第三步：由相关函数的傅式变换求功率谱。即 jwm M M m X jw x e m R e S ----=∧∧ ∑= )()(1) 1( 以上过程中经历了两次截断，一次是将x(n)截成N 长，称为加数据窗，一次是将x(n)截成（2M-1）长，称为加延迟窗。因此所得的功率谱仅是近似值，也叫谱估计，式中的)(jw x e S 代表估值。一般取M<

实验2追赶法算法设计及MATLAB实现

数值计算方法 实 验 报 告 实验序号：实验二 实验名称：追赶法算法设计及MATLAB实现 实验人： 专业年级： 教学班： 学号： 实验时间：

实验二追赶法算法设计及MATLAB实现 一、实验目的 １．初步掌握算法设计规则； ２．初步掌握MATLAB程序设计规则． 二、实验内容 1．构造利用追赶法求解三对角线性方程组的算法； 2．在MATLAB环境下编写追赶法的程序（函数）； 3．自由选择若干个三对角线性方程组求解。 三、实验步骤 1．追赶法算法： 算法名称：thomas 输入参数：向量a,b,c,f 输出参数：输出解信息x 算法的自然语言： Step1：u 1=b 1 ,y 1 =b 1 ; Step2：对于i=2,3,….n; Step2.1：当u 1-i ≠，否则转step5 l i =a i /u 1-i ; u i =b i -l i *c 1-i ； y i =f i -l i *y 1-i ； Step3：当u n ≠时，x n =y n /u n ，否则转step5 Step4：对于：i=n-1,n-2,…..,2,1，转step6 x i =(y i -c i *x 1+i )/u i Step5：无解信息，转step7

Step6：输出x Step7：关机 2．MATLAB程序 function [x,L,U]=thomas(a,b,c,f) n=length(b); % 对A进行分解 u(1)=b(1); for i=2:n if(u(i-1)~=0) l(i-1)=a(i-1)/u(i-1); u(i)=b(i)-l(i-1)*c(i-1); else break; end end L=eye(n)+diag(l,-1); U=diag(u)+diag(c,1); x=zeros(n,1); y=x; % 求解Ly=b y(1)=f(1); for i=2:n y(i)=f(i)-l(i-1)*y(i-1); end % 求解Ux=y if(u(n)~=0) x(n)=y(n)/u(n); end for i=n-1:-1:1 x(i)=(y(i)-c(i)*x(i+1))/u(i);

利用MATLAB实现信号DFT的计算

07级电信(2)班 刘坤洋 24 实验一 利用MATLAB 实现信号DFT 的计算 一、实验目的: 1、熟悉利用MATLAB 计算信号DFT 的方法 2、掌握利用MATLAB 实现由DFT 计算线性卷积的方法 二、实验设备:电脑、matlab 软件 三、实验内容: 1、练习用matlab 中提供的内部函数用于计算DFT （1） fft(x),fft(x,N),ifft(x),ifft(x,N)的含义及用法 （2） 在进行DFT 时选取合适的时域样本点数N 请举例,并编程实现 题目: 源程序: >> N=30; %数据的长度 >>L=512; %DFT 的点数 >>f1=100; f2=120; >>fs=600; %抽样频率 >>T=1/fs; %抽样间隔 >>ws=2*pi*fs; >>t=(0:N-1)*T; >>f=cos(4*pi*f1*t)+cos(4*pi*f2*t); >>F=fftshift(fft(f,L)); >>w=(-ws/2+(0:L-1)*ws/L)/(2*pi); >>hd=plot(w,abs(F)); >>ylabel('幅度谱') >> xlabel('频率/Hz') 的频谱 分析利用)π4cos()π4cos()(DFT 21t f t f t x +=Hz 600,Hz 120,Hz 10021===s f f f

>> title('my picture') 结果图 : （3） 在对信号进行DFT 时选择hamming 窗增加频率分辨率 请举例,并编程实现 题目: 源程序:>> N=50; %数据的长度 >>L=512; %DFT 的点数 >>f1=100;f2=150; >>fs=600; %抽样频率 >>T=1/fs; %抽样间隔 >>ws=2*pi*fs; >>t=(0:N-1)*T; >>f=cos(4*pi*f1*t)+0、15*cos(4*pi*f2*t); 的频谱分析利用)π4cos(15.0)π4cos()(DFT 21t f t f t x +=Hz 600,Hz 150,Hz 10021===s f f f

matlab实现有限差分法计算电场强度(最新)

实验一：有限差分法研究静电场边值问题 实验报告人：年级和班级：学号： 1. 实验用软件工具： Matlab 2. 实验原理：电磁场课本P36-38 1)差分方程 2)差分方程组的解 简单迭代法 高斯-赛德尔迭代法 逐次超松弛法 3. 实验步骤： 1）简单迭代法 程序： hx=41;hy=21; v1=zeros(hy,hx); v1(hy,:)=zeros(1,hx); v1(1,:)=ones(1,hx)*100; v1(:,1)=zeros(hy,1); v1(:,hx)=zeros(hy,1); v1 v2=v1;maxt=1;t=0; k=0; while(maxt>1e-5) k=k+1; maxt=0; for i=2:hy-1 for j=2:hx-1 v2(i,j)=(v1(i,j+1)+v1(i+1,j)+v1(i-1,j)+v1(i,j-1))/4; t=abs(v2(i,j)-v1(i,j)); if(t>maxt) maxt=t;end end end v1=v2; end v2 k clf subplot(1,2,1),mesh(v2) axis([0,41,0,21,0,100]) subplot(1,2,2),contour(v2,15) hold on

axis([-1,42,-1,25]) plot([1,1,hx,hx,1],[1,hy+1,hy+1,1,1],'r') text(hx/2,0.3,'0V','fontsize',11); text(hx/2-0.5,hy+0.5,'100V','fontsize',11); text(-0.5,hy/2,'0V','fontsize',11); text(hx+0.3,hy/2,'0V','fontsize',11); hold off 当W=1e-5, 迭代次数：1401次 2）高斯-赛德尔迭代法 程序： hx=41;hy=21; v1=ones(hy,hx); v1(hy,:)=zeros(1,hx); v1(1,:)=ones(1,hx)*100; v1(:,1)=zeros(hy,1); v1(:,hx)=zeros(hy,1); v2=v1;maxt=1;t=0; k=0; while(maxt>1e-5) k=k+1; maxt=0; for i=2:hy-1 for j=2:hx-1 v2(i,j)=(v1(i,j+1)+v1(i+1,j)+v2(i-1,j)+v2(i,j-1))/4; t=abs(v2(i,j)-v1(i,j)); if(t>maxt) maxt=t;end end end v1=v2; end v2 k clf subplot(1,2,1),mesh(v2) axis([0,41,0,21,0,100]) subplot(1,2,2),contour(v2,15) hold on axis([-1,42,-1,25]) plot([1,1,hx,hx,1],[1,hy+1,hy+1,1,1],'r') text(hx/2,0.3,'0V','fontsize',11); text(hx/2-0.5,hy+0.5,'100V','fontsize',11); text(-0.5,hy/2,'0V','fontsize',11); text(hx+0.3,hy/2,'0V','fontsize',11); hold off

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

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

表9-1最小化函数表 函数描述 fgoalattain多目标达到问题 fminbnd有边界的标量非线性最小化 fmincon有约束的非线性最小化 fminimax最大最小化 fminsearch,fminunc无约束非线性最小化 fseminf半无限问题 linprog线性课题 quadprog二次课题 2．方程求解函数 表9-2方程求解函数表 函数描述 \线性方程求解 fsolve非线性方程求解 fzero标量非线性方程求解 3．最小二乘（曲线拟合）函数 表9-3最小二乘函数表 函数描述 \线性最小二乘 lsqlin有约束线性最小二乘 lsqcurvefit非线性曲线拟合 lsqnonlin非线性最小二乘 lsqnonneg非负线性最小二乘

《应用计算方法教程》matlab作业二解析

6-1 试验目的 计算特征值，实现算法 试验内容：随机产生一个10阶整数矩阵，各数均在-5和5之间。 （1） 用MATLAB 函数“eig”求矩阵全部特征值。 （2） 用幂法求A 的主特征值及对应的特征向量。 （3） 用基本QR 算法求全部特征值（可用MATLAB 函数“qr ”实现矩阵的QR 分解）。 原理 幂法：设矩阵A 的特征值为12n ||>||||λλλ≥???≥并设A 有完全的特征向量系12,,,n χχχ???(它们线性无关)，则对任意一个非零向量0n V R ∈所构造的向量序列1k k V AV -=有11()lim ()k j k k j V V λ→∞ -=， 其中()k j V 表示向量的第j 个分量。 为避免逐次迭代向量k V 不为零的分量变得很大（1||1λ> 时）或很小（1||1λ< 时），将每一步的k V 按其模最大的元素进行归一化。具体过程如下： 选择初始向量0V ，令1max(),,,1k k k k k k k V m V U V AU k m +===≥，当k 充分大时1111,max()max() k k U V χλχ+≈ ≈。 QR 法求全部特征值： 111 11222 111 ,1,2,3,k k k k k A A Q R R Q A Q R k R Q A Q R +++==????==??=???? ??????==?? 由于此题的矩阵是10阶的，上述算法计算时间过长，考虑采用改进算法——移位加速。迭 代格式如下： 1 k k k k k k k k A q I Q R A R Q q I +-=?? =+? 计算k A 右下角的二阶矩阵() () 1,1 1,() (),1 ,k k n n n n k k n n n n a a a a ----?? ? ??? 的特征值()()1,k k n n λλ-，当()()1,k k n n λλ-为实数时，选k q 为()()1,k k n n λλ-中最接近(),k n n a 的。 程序