数模培训1(MATLAB基础)
数模培训
一、曲线插值与拟合
二、数值微分与积分
三、微分方程数值解
四、优化问题
五、回归分析
1.一维插值
对表格给出的函数，求出没有给出的函数值。
在实际工作中，经常会遇到插值问题。
例1：表1是待加工零件下轮廓线的一组数据，现需要得到x坐标每改变0.1时所对应的y的坐标.
下面是关于插值的两条命令(专门用来解决这类问题)：
y=interp1(x0,y0,x) 分段线性插值
y=spline(x0,y0,x) 三次样条插值
其中x0,y0是已知的节点坐标，是同维向量。y对应于x处的插值。y与x是同维向量。解决上述问题,我们可分两步:
一用原始数据绘图作为选用插值方法的参考.
二确定插值方法进行插值计算
对于上述问题,可键入以下的命令:
x0=[0,3,5,7,9,11,12,13,14,15]';
y0=[0,1.2,1.7,2.0,2.1,2.0,1.8,1.2,1.0,1.6]'
plot(x0,y0) %完成第一步工作
x=0:0.1:15;
y=interp1(x0,y0,x'); %用分段线性插值完成第二步工作
plot(x,y)
y=spline(x0,y0,x');
plot(x,y) %用三次样条插值完成第二步工作
练习：对y=1/(1+x2),-5≤x≤5，用n（=11）个节点（等分）作上述两种插值，用m（=21）个插值点（等分）作图，比较结果。
解：键入并运行如下命令
n=11;m=21;x=-5:10/(m-1):5;y=1./(1+x.^2);
xo=-5:10/(n-1):5;yo=1./(1+xo.^2);
y1=interp1(xo,yo,x);
y2=spline(xo,yo,x);
plot(x,y,'r',x,y1,'b',x,y2,'k')
练习：在某处测得海洋不同深度处水温如下：
解：输入程序：
D=[446,714,950,1422,1634];
T=[7.04,4.28,3.40,2.54,2.13];
Di=[500,1000,1500];
Ti=interp1(D,T,Di)
MATLAB的命令interp1(X,Y,Xi,’method’)用于一元插值.其中Method可选’nearest’(最近邻插值),’linear’(线性插值),’spline’(三次样条插值),’cubic’(三次多项式插值)
2.二维插值
MATLAB中二维插值的命令是：
z=interp2（x0,y0,z0,x,y,'meth'）
例2：在一个长为5个单位，宽为3个单位的金属薄片上测得15个点的温度值，试求出此薄片的温度分布，并绘出等温线图。（数据如下表）
程序：temps=[82,81,80,82,84;79,63,61,65,87;84,84,82,85,86];
mesh(temps) %根据原始数据绘出温度分布图，可看到此图的粗造度。
%下面开始进行二维函数的三阶插值。
width=1:5; depth=1:3; di=1:0.2:3; wi=1:0.2:5;
[WI,DI]=meshgrid(wi,di);%增加了节点数目
ZI=interp2(width,depth,temps,WI,DI,'cubic');% 对数据（width,depth,temps）进 % 行三阶插值拟合。
surfc(WI,DI,ZI)
contour(WI,DI,ZI)
3．曲线拟合
假设一函数g(x)是以表格形式给出的，现要求一函数f(x)，使f(x)在某一准则下与表格函数（数据）最为接近。
由于与插值的提法不同，所以在数学上理论根据不同，解决问题的方法也不同。
此处，我们总假设f(x)是多项式。
例3：弹簧在力F的作用下伸长x厘米。F和x在一定的范围内服从虎克定律。试根据下
在MATLAB中，用以下命令拟合多项式。
polyfit（x0,y0,n）
一般，也需先观察原始数据的图像，然后再确定拟和成什么曲线。
对于上述问题,可键入以下的命令:
x=[1,2,4,7,9,12,13,15,17]';
F=[1.5,3.9,6.6,11.7,15.6,18.8,19.6,20.6,21.1]';
plot(x,F,'.')
从图像上我们发现：前5个数据应与直线拟合，后5个数据应与二次曲线拟合。于是键入a=polyfit(x(1:5),F(1:5),1); a=polyfit(x(5:9),F(5:9),2)
得
注意：有时，面对一个实际问题，究竟是用插值还是用拟合不好确定，还需大家在实际中仔细区分。同时，大家（包括学过计算方法的同学）注意去掌握相应的理论知识。
4．数值积分
先看一个例子：
例4．现要根据瑞士地图计算其国土面积。于是对地图作如下的测量：以西东方向为横轴，以南北方向为纵轴。（选适当的点为原点）将国土最西到最东边界在x轴上的区间划取足够多的分点x i，在每个分点处可测出南北边界点的对应坐标y1，y2。用这样的方
根据地图比例知18mm相当于40km，试由上表计算瑞士国土的近似面积。（精确值为41288km2）。
解题思路：数据实际上表示了两条曲线，实际上我们要求由两曲线所围成的图形的面积。解此问题的方法是数值积分的方法。具体解时我们遇到两个问题：1。数据如何输入；2。没有现成的命令可用。
解：对于第一个问题，我们可把数据考备成M文件（或纯文本文件）。
然后，利用数据绘制平面图形。键入
load mianji.txt
A=mianji';
plot(A(:,1),A(:,2),'r',A(:,1),A(:,3),'g')
接下来可以计算面积。键入：
a1=trapz(A(:,1)*40/18,A(:,2)*40/18);
a2=trapz(A(:,1)*40/18,A(:,3)*40/18);
d=a2-a1
d = 4.2414e+004
至此，问题可以说得到了解决。之所以说还有问题，是我们觉得误差较大。但计算方法的理论给了我们更精确计算方法。只是MATLAB没有相应的命令。想得到更理想的结果，我们可以自己设计解决问题的方法。（可以编写辛普森数值计算公式的程序，或用拟合的方法求出被积函数，再利用MA TLAB的命令quad,quad8）
5．数值微分
实际的例：已知20世纪美国人口统计数据如下，根据数据计算人口增长率。（其
解题思路：设人口是时间的函数x(t).于是人口的增长率就是x(t)对t 的导数.如果计
算出人口的相关变化率x
x
t r =)(。那么人口增长满足)()(t x t r x
= ，它在初始条件x(0)=x0下的解为?=t
du
u r e x t x 0
)(0)(.（用以检查计算结果的正确性）
解：此问题的特点是以离散变量给出函数x(t),所以就要用差分来表示函数x(t)的导数.
2)
()()(' , )()()(' , )()()('h
h a f h a f a f h h a f a f a f h a f h a f a f --+≈--≈-+≈
常用后一个公式。（因为，它实际上是用二次插值函数来代替曲线x(t)）即常用三点公式来代替函数在各分点的导数值：
h y y y x f h
y y y x f h
y y x f n
n n n k
k k 234)('243)('1-n 1,2,k 2)('12
2
1001+-≈-+-≈
=-≈
--+ MATLAB 用命令diff 按两点公式计算差分;此题自编程序用三点公式计算相关变化率.编程如下: for i=1:length(x) if i==1
r(1)=(-3*x(1)+4*x(1+1)-x(1+2))/(20*x(1)); elseif i~=length(x)
r(i)=(x(i+1)-x(i-1))/(20*x(i)); else
r(length(x))=(x(length(x)-2)-4*x(length(x)-1)+3*x(length(x)))/(20*x(length(x))); end end r=r;
保存为diff3.m 文件听候调用.再在命令窗内键入
X=[1900,1910,1920,1930,1940,1950,1960,1970,1980,1990];
x=[76.0, 92.0, 106.5, 123.2, 131.7, 150.7, 179.3, 204.0, 226.5, 251.4]; diff3;
由于r 以离散数据给出,所以要用数值积分计算.键入 x(1,1)*exp(trapz(X(1,1:9),r(1:9)))
数值积分命令：trapz(x),trapz(x,y),quad(‘fun ’,a,b)等.
6：微分方程数值解（单摆问题）
单摆问题的数学模型是
θθ
sin g l -= 在初始角度不大时,问题可以得到很好地解决,但如果初始角较大,此方程无法求出解
析解.现问题是当初始角为100和300时,求出其解,画出解的图形进行比较。 解:若θ0
较小,则原方程可用0=+θθ
g l 来近似.其解析解为θ(t)= θ0cos ωt,g
l
=
ω. 若不用线性方程来近似,那么有两个模型:
30)0(sin 10)0(sin 00??
?
?
?=-=?????=-=θθθθθθl g l g 取g=9.8,l=25, 100=0.1745, 300=0.5236.用MA TLAB 求这两个模型的数值解,先要作如下
的处理:令x1=θ,x2=θ’,则模型变为
0)0(,5236
.0)0(sin l g - x x x 0)0(,1745.0)0(sin l g - x x x 211221211
221??
?????====???????====x x x x x x 再编函数文件
function xdot=danbai(t,x) xdot=zeros(2,1);
xdot(1)=x(2);xdot(2)=-9.8/25*sin(x(1)); 在命令窗口键入
[t,x]=ode45(‘danbai ’,[0:0.1:20],[0.1745,0]); [t,y]=ode45(‘danbai ’,[0:0.1:20],[0.5236,0]); plot(t,x(:,1),’r ’,t,y(:,1),’k ’);
参考书:数学实验,高等教育出版社
7．解优化问题
线性规划有约束极小问题 模型
UB
x LB , , . .,min 11≤≤=≤=b x A b Ax t s cx z
用命令
[x, fval]= linprog(f,A,b,A1,b1,lb,ub) 例1：Find x that minimizes f(x )=-5x 1-4x 2-6x 3 subject to x 1-x 2+x 3≦20 3x 1+2x 2+4x 3≦42 3x 1+2x 2≦30
0≦x 1, 0≦x 2,0≦x 3
First, enter the coefficients: f = [-5; -4; -6] A = [1 -1 1 3 2 4 3 2 0]; b = [20; 42; 30]; lb = zeros(3,1);
Next, call a linear programming routine:
[x,fval,exitflag,output,lambda] = linprog(f,A,b,[],[],lb)； Entering x, fval,lambda.ineqlin, and lambda.lower gets x =
0.0000 15.0000 3.0000 fval =
-78.0000 和其它信息。 例2：解问题
,, x ,1052x ,7 x ..532x z max 321321321321≥≥+-=++-+=x x x x x x t s x x
把问题极小化并将约束标准化
3
21321321321,,x 0 ,7 x ,1052x - ..532x z min x x x x x x t s x x ≤=++-≤-++--=
键入c=[-2,-3,5];a=[-2,5,-1];b=-10;a1=[1,1,1];b1=7;LB=[0,0,0];
[x ，y]=linprog （c,a,b,a1,b1,LB ）得当X=(6.4286,0.5714,0.0000)时,z=-14.5714最大. 例3: 解问题
.
5,0,0 x ,1222x ,44 x ,62 x ..2x z min 321321321321321≤≥≥≤+-≤-+=+++--=x x x x x x x x t s x x
解:键入
c=[-2,-1,1];a=[1,4,-1;2,-2,1];b=[4;12];a1=[1,1,2];b1=6;lb=[0;0;-inf];ub=[inf;inf;5]; [x,z]=linprog(c,a,b,a1,b1,lb,ub) 得当X=(4.6667,0.0000,0.6667)时,z=-8.6667最小.
非线性规划有约束极小问题 模型1
x g ≤∈)(..),(min t s R x x f n
用命令x=constr('f ',x0)。 Examples
Find values of x that minimize f(x)=-x1x2x3, starting at the point x = [10; 10; 10] and subject to the constraints 0≤x1+2x2+2x3≤72. -x1-2x2-2x3≤0,x1+2x2+2x3≤72, 第一步：编写M 文件 function [f,g]=myfun(x) f=-x(1)*x(2)*x(3);
g(1)=-x(1)-2*x(2)-2*x(3); g(2)=x(1)+2*x(2)+2*x(3)-72;
第二步：求解
在MA TLAB 工作窗中键入 x0=[10,10,10];
x=constr('myfun',x0)即可
模型2：
UB.
x LB ,0)(c 0,c(x) ,A b,Ax ..),(min 111≤≤=≤=≤x b x t s x f
MATLAB 求解此问题的命令是：
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(‘fun ’,x0,A,b,A1,b1,LB,UB,’nonlcon ’, options,p1,p2,…)
fun 是目标函数的m_文件名.nonlcon 是约束函数C(x)和C1(x)的m_文件名.文件输出为[C,C1].
例:求解最优化问题
.
010x - ,0x 1.5 ,
0 x ..),
12424()(min 212121212212
2211≤-≤--+=+++++=x x x x x t s x x x x x e x f x
第1步: 建立目标函数和非线性约束的m_文件.
Function y=e1511(x)% 目标函数的m_文件
Y=exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);
Function [c1,c2]=e1511b(x)% 非线性约束的m_文件 C1=[1.5+x(1)*x(2)-x(1)-x(2);-x(1)*x(2)-10]; C2=0;
第2步: 运行程序.键入 x0=[-1,1];a1=[1,1];b1=0;
[x,f,exitflab,output]=fmincon(‘e1511’,x0,[],[],a1,b1,[],[],’e1511b ’) 得结果.
输出结果的意义:经过4次迭代(iterations:4)收敛到了(exitfag=1)最优解x(1)=-1.2247,x(2)=1.2247,目标函数最优值为1.8951.
非线性无约束极小问题
R x x f ∈),(min
用命令x=fmin('f',x0)。
或用命令x=fminu('f ',x0)，或用命令x=fmins('f ',x0)。
非线性最小二乘问题
)()(min x f x f T
用命令x=leastsq('f ',x0)，或用命令x=curvefit('f ',x0)。
二次规划
b
Ax t s x c Hx x T T ≤+..2/min
用命令x=qp(H,c,A,b)。
关于这些命令的详细使用规则和例子，用借助help 进行查阅。 参考书:科学计算技术与MATLAB,科学出版社
8．回归分析
前面我们曾学过拟合。但从统计的观点看，对拟合问题还需作回归分析。例如：有描述问题甲和问题乙的两组数据(x,y)和(x,z)。设x=[1，2，3，4]；y=[1.0，1.3，1.5,2.02.3];z=[0.6,1.95,0.9,2.85,1.8]。如果在平面上画出散点图，那么问题甲的四个点基本在一条直线上而问题乙的四个点却很散乱。如果都用命令polyfit(x,y,1)，polyfit(x,z,1)来拟合，将得到同一条直线。自然对问题甲的信任程度会高于对问题乙的信任程度。所以有必要对所得结果作科学的评价分析。回归分析就是解决这种问题的科学方法。下面结合三个具体的例子介绍MATLAB 实现回归分析的命令。
中有无异常点、由x 的取值对y 作出预测。
2． 将17至19岁的运动员每两岁一组分为7组，每组两人测量其旋转定向能力，以考察年龄(x)对这种运动能力(y)的影响。现得到一组数据如下表
3． 某厂生产的某产品的销售量与竞争对手的价格x1和本厂的价格x2有关。下表是该
170（元），预测此产品在该城市的销售量。
解1：
在x-y平面上画散点图，直观地知道y与x大致为线性关系。
用命令polyfit(x,y,1)可得y=140.6194x+27.0269。
作回归分析用命令
[b,bint,r,rint,ststs]=regress(y,x,alpha) 可用help查阅此命令的具体用法
残差及置信区间可以用rcoplot(r,rint)画图
设回归模型为y=β0+β1x,
在MATLAB命令窗口中键入下列命令进行回归分析
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];
[b,bint,r,rint,stats]=regress(y,X,0.05);
b,bint,stats,rcoplot(r,rint)
得结果和图
b =
27.0269
140.6194
bint =
22.3226 31.7313
111.7842169.4546
stats =
0.9219 118.0670 0.0000
结果含义为
β0=27.0269 β1=140.6194
β0的置信区间是 [22.3226,31.7313]
β1的置信区间是 [111.7842,169.4546]
R2=0.9219 F=118.0670, p<10-4.
R是衡量y与x的相关程度的指标，称为相关系数。R越大，x与y关系越密切。
通常R大于0.9才认为相关关系成立。
F是一统计指标
p是与F对应的概率，当 p<0.05时，回归模型成立。
此例中 p=0 <10-4<0.05，所以，所得回归模型成立。
观察所得残差分布图，看到第8个数据的残差置信区间不含零点，此点视为异常点，剔除后重新计算。
此时键入：
X(8,:)=[];
y(8)=[];
[b,bint,r,rint,stats]=regress(y,X);
b,bint,stats,rcoplot(r,rint)
得：
b =
27.0992
137.8085
bint =
23.8563 30.3421
117.8534 157.7636
stats =
0.9644 244.0571 0.0000
可以看到：置信区间缩小；R2、F变大，所以应采用修改后的结果。
解2：
在x-y平面上画散点图，直观地知道y与x大致为二次函数关系。
设模型为y=a1x2+a2x+a3
此问题可以利用命令polyfit(x,y,2)来解，也可以像上题一样求解。下面介绍用命令polytool来解。
首先在命令窗口键入
x=17:2:29;x=[x,x];
y=[20.48,25.13,26.15 30,26.1,20.3,19.35,24.35,28.11,26.3,31.4,26.92,25.7,21.3];
polytool（x,y,2）
得到一个交互式窗口
182022242628
10
15
20
25
30
35
窗口中绿线为拟合曲线、红线为y 的置信区间、可通过移动鼠标的十字线或通过在窗口下方输入来设定x 值，窗口左边则输出与x 对应的y 值及y 的置信区间。通过左下方的Export 下拉菜单可输出回归系数等。更详细的解释可通过help 查阅。
解3.这是一个多元回归问题。若设回归模型是线性的，即设y=β0+β1x 1+β2x 2 那么依然用regress(y,x,alpha)求回归系数。键入 x1=[120,140,190,130,155,175,125,145,180,150]; x2=[100,110,90,150,210,150,250,270,300,250]; y=[102,100,120,77,46,93,26,69,65,85]'; x=[ones(10,1),x1',x2'];
[b,bint,r,rint,stats]=regress(y,x);b,bint,stats, 得 b =
66.5176 0.4139 -0.2698 bint =
-32.5060 165.5411 -0.2018 1.0296 -0.4611 -0.0785 stats =
0.6527 6.5786 0.0247
p=0.0247,若显著水平取0,01，则模型不能用；R2=0.6527较小;β0，β1的置信区间包含零点。因此结果不理想。于是设模型为二次函数。此题设模型为纯二次函数：
y=β0+β1x 1+β2x 2+β11x 12+β22x 22
MATLAB 提供的多元二项式回归命令为rstool(x,y,model,alpha).其中alpha 为显著水平、model 在下列模型中选一个： linear （线性）
purequadratic （纯二次）
interaction （交叉）
quadratic （完全二次）
对此例，在命令窗中键入 x(:,1)=[];
rstool(x,y,'purequadratic') 得到一个对话窗：
140160-100
150200
其意义与前面的对话窗意义类似。若要回答“本厂售价160，对手售价170，预测该市销售量”的问题，只需在下方窗口中分别肩入160和170，就可在左方窗口中读到答案及其置信区间。下拉菜单Export 向工作窗输出数据具体操作为：
∑=+
++=m
j j
ij m m x
x x y 1
2110ββββ m
m x x y βββ ++=110k
m
k j j jk m m x x x x y ∑≤≠≤+
++=1110β
βββ ∑≤≤+
++=m
k j k
j jk m m x x x x y ,1110β
βββ
弹出菜单，选all，点击确定。此时可到工作窗中读取数据。可读数据包括：beta（回归系数）rmse（剩余标准差）residuals （残差）本题只要键入
beta,rmse,residuals
得
beta =
-312.5871
7.2701
-1.7337
-0.0228
0.0037
rmse =
16.6436
residuals =
6.6846
-12.6703
-0.2013
6.4855
-19.6533
7.9989
-11.4737
5.4303
-4.9932
22.3926
大家不妨用此命令选其它模型作一个比较。
注：关于各种模型，可查阅有关理论
参考书:数学实验,高等教育出版社
参考书:科学计算技术与MATLAB,科学出版社
MATLAB中还包括神经网络工具箱，小波分析工具箱，在网上还可以下载遗传算法工具箱，有兴趣的同学可以借这次机会，结合学习MA TLAB，好好学习一下相关理论知识。最后，祝大家学习，竞赛都取得成功。谢谢大家。
第一讲Matlab基本数值计算(可编辑修改word版)
第一讲Matlab 基本数值计算 一、矩阵 在Matlab 中，一个矩阵可以使数学意义上的矩阵，也可以是标量或者向量。对于一个标量（一个数）可以将之作为1?1 的矩阵，而向量（一行或一列）则可以认为是1?n 或者n ?1的矩阵。另外，一个0 ?0 矩阵在Matlab 中被认为是空矩阵，用“[]”表示。 1、矩阵的创建 矩阵的创建可以有以下几种形式 ⑴直接输入 >> A=[1 2 3;4 3 7;2 4 1] 注意：每行间的元素用逗号或空格分开，行与行之间用分号或回车分开，矩阵标示是一对中括号[ ]。 也可以采用数组编辑器（Array Editor）像在Excel 电子表格中据那样输入数据。 ⑵通过语句和函数产生 常用的特殊矩阵：zeros:全零矩阵,ones：全1 矩阵，eye：单位矩阵，rand:随机矩阵，diag:对角阵等。 例：>> A=ones(3,4) >> E=eye(3) >> D=diag([3 5 2]) ⑶对矩阵进行裁剪或拼接 ⑷从外部文件装入数据
外部数据文件可以是以保存的Matlab 工作空间，也可以是文本（.txt）文件，或者是电子表格创建的文件（.xls）. 例：已知一个文本格式的数据文件E:\Mathmodel\data1.txt >> load e:\Mathmodel\data1.txt 得到一个变量名与文件名相同的矩阵（data1）。注意：文件的扩展名不能省略。 例：已知一个Excel 文件的路径为E:\Mathmodel\data2.xls a.缺省操作: >> NUMBER=xlsread('E:\Mathmodel\data2.xls') >>[NUMBER,TXT]=xlsread('E:\Mathmodel\data2.xls') 默认操作是从第一个工作表（sheet1）中提取数据。 b.从指定的工作表（而不是第一个）中提取数据： >> NUMBER=xlsread('E:\Mathmodel\data2.xls','S2') 或者 >> NUMBER=xlsread('E:\Mathmodel\data2.xls',2) c.从指定的工作表中读取指定区域的数据： >> NUMBER=xlsread('E:\Mathmodel\data2.xls',2,'g3:i8') 2、Matlab 的矩阵运算 ⑴基本运算 矩阵的加（+）、减（-）、乘（*）、乘方（^）运算法则与代数中的定义完全一致。例如： >> A=[1 2;3 4];B=[3 1;4 8];
数模培训1(MATLAB基础)
数模培训 一、曲线插值与拟合 二、数值微分与积分 三、微分方程数值解 四、优化问题 五、回归分析 1.一维插值 对表格给出的函数，求出没有给出的函数值。 在实际工作中，经常会遇到插值问题。 例1：表1是待加工零件下轮廓线的一组数据，现需要得到x坐标每改变0.1时所对应的y的坐标. 下面是关于插值的两条命令(专门用来解决这类问题)： y=interp1(x0,y0,x) 分段线性插值 y=spline(x0,y0,x) 三次样条插值 其中x0,y0是已知的节点坐标，是同维向量。y对应于x处的插值。y与x是同维向量。解决上述问题,我们可分两步: 一用原始数据绘图作为选用插值方法的参考. 二确定插值方法进行插值计算 对于上述问题,可键入以下的命令: x0=[0,3,5,7,9,11,12,13,14,15]'; y0=[0,1.2,1.7,2.0,2.1,2.0,1.8,1.2,1.0,1.6]' plot(x0,y0) %完成第一步工作 x=0:0.1:15; y=interp1(x0,y0,x'); %用分段线性插值完成第二步工作 plot(x,y) y=spline(x0,y0,x'); plot(x,y) %用三次样条插值完成第二步工作 练习：对y=1/(1+x2),-5≤x≤5，用n（=11）个节点（等分）作上述两种插值，用m（=21）个插值点（等分）作图，比较结果。 解：键入并运行如下命令 n=11;m=21;x=-5:10/(m-1):5;y=1./(1+x.^2); xo=-5:10/(n-1):5;yo=1./(1+xo.^2); y1=interp1(xo,yo,x); y2=spline(xo,yo,x); plot(x,y,'r',x,y1,'b',x,y2,'k')
数学建模(Matlab)
数学规划作业(MatLab) 1、某厂向用户提供发动机，合同规定，第一、二、三季度末分别交货40台、60台、80台．每季度的生产费用为()2 =+ f x ax bx （单位:元）, 其中x是该季度生产的台数．若交货后有剩余，可用于下季度交货，但需支付存储费，每台每季度c元．已知工厂每季度最大生产能力为100台，第一季度开始时无存货，设a=50、b=0.2、c=4，问:工厂应如何安排生产计划，才能既满足合同又使总费用最低．讨论a、b、c变化对计划的影响，并作出合理的解释． 解： 问题的分析和假设： 分析： 问题的关键在于由于工厂的生产能力足以满足每个季度用户的需求，但是为了使总费用最少，那么利用每个季度生产费用的不同，可用利用上个生产费用低的季度多生产来为下个季度进行准备，前提是本月节省下的费用减去总的发动机存储费用还有剩余，这样生产才有价值，才可能满足合同的同时又能使总费用最低。基本假设：1工厂的生产能力不受外界环境因素影响。2为使总费用最低，又能满足合同要求，各个季度之间的生产数量之间是有联系的。3第一季度开始时无存货。4工厂每季度的生关费用与本季度生产的发动机台数有关。5生产要按定单的数量来进行，生产的数量应和订单的数量相同，以避免生产出无用的机器。 符号规定：X1―――第一季度生产发动机的数量 X2―――第二季度生产发动机的数量
X3―――第三季度生产发动机的数量 建模： 1.三个季度发动机的总的生产量为180台。 2.每个季度的生产量和库存机器的数量之和要大于等于本季度的交货数量。 3.每个月的生产数量要符合工厂的生产能力。 4.将实际问题转化为非线性规划问题，建立非线性规划模型 目标函数 min f(x)=50(x1+x2+x3)+0.2(x12+x22+x32)+4(x1-40)+4(x1+x2-100) 整理，得 min f(x)=50(x1+x2+x3)+0.2(x12+x22+x32)＋4(2x1+x2-140) 约束函数s.t x1+x2≥100; x1+x2+x3=180; 40≤x1≤100; 0≤x2≤100; 0≤x3≤100; 求解的Matlab程序代码： M-文件 fun.m: function f=fun (x); f=50*(x(1)+x(2)+x(3))+0.2*(x(1)^2+x(2)^2+x(3)^2)+4*(2*x(1) +x(2)-140)主程序fxxgh.m：
【数学建模学习】matlab作图
基本形式 >> y=[0 0.58 0.70 0.95 0.83 0.25]; >> plot(y) 生成的图形是以序号为横坐标、数组y的数值为纵坐标画出的折线。 >> x=linspace(0,2*pi,30); % 生成一组线性等距的数值 >> y=sin(x); >> plot(x,y) 生成的图形是上30个点连成的光滑的正弦曲线。 多重线在同一个画面上可以画许多条曲线，只需多给出几个数组，例如>> x=0:pi/15:2*pi; >> y1=sin(x)；>> y2=cos(x)；>> plot(x,y1,x,y2) 则可以画出多重线。 另一种画法是利用hold命令。在已经画好的图形上，若设置hold on，MATLA将把新的plot命令产生的图形画在原来的图形上。而命令hold off 将结束这个过程。例如： >> x=linspace(0,2*pi,30); y=sin(x); plot(x,y) >> hold on >> z=cos(x); plot(x,z) >> hold off 线型和颜色MATLAB对曲线的线型和颜色有许多选择，标注的方法是在每一对数组后加一个字符串参数，说明如下：线型线方式：- 实线:点线-. 虚点线- - 波折线。线型 点方式： . 圆点+加号* 星号x x形o 小圆颜色：y黄；r红；g绿；b蓝；w白；k黑；m 紫；c青. 以下面的例子说明用法：>> x=0:pi/15:2*pi; >> y1=sin(x); y2=cos(x); >> plot(x,y1,’b:+’,x,y2,’g-.*’) 网格和标记在一个图形上可以加网格、标题、x轴标记、y轴标记，用下列命令完成这些工作。>> x=linspace(0,2*pi,30); y=sin(x); z=cos(x); >> plot(x,y,x,z) >> grid >> xlabel(‘Independent Variable X’) >> ylabel(‘Dependent Variables Y and Z’) >> title(‘Sine and Cosine Curves’) 也可以在图形的任何位置加上一个字符串，如用：>> text(2.5,0.7,’sinx’) 表示在坐标x=2.5, y=0.7处加上字符串sinx。 更方便的是用鼠标来确定字符串的位置，方法是输入命令：>> gtext(‘sinx’) 在图形窗口十字线的交点是字符串的位置，用鼠标点一下就可以将字符串放在那里。 坐标系的控制在缺省情况下MATLAB自动选择图形的横、纵坐标的比例，如果你对这个比例不满意，可以用axis命令控制，常用的有：axis([xmin xmax ymin ymax]) [ ]中分别给出x轴和y轴的最大值、最小值axis equal 或axis(‘equal’) x轴和y轴的单位长度相同axis square 或axis(‘square’) 图框呈方形axis off 或axis(‘off’) 清除坐标刻度还有axis auto axis image axis xy axis ij axis normal axis on axis(axis) 用法可参考在线帮助系统。 多幅图形可以在同一个画面上建立几个坐标系, 用subplot(m,n,p)命令；把一个画面分成m×n个图形区域, p代表当前的区域号，在每个区域中分别画一个图,如>> x=linspace(0,2*pi,30); y=sin(x); z=cos(x); >> u=2*sin(x).*cos(x); v=sin(x)./cos(x); >> subplot(2,2,1),plot(x,y),axis([0 2*pi –1 1]),title(‘sin(x)’) >> subplot(2,2,2),plot(x,z),axis([0 2*pi –1 1]),title(‘cos(x)’) >> subplot(2,2,3),plot(x,u),axis([0 2*pi –1 1]),title(‘2sin(x)cos(x)’) >> subplot(2,2,4),plot(x,v),axis([0 2*pi –20 20]),title(‘sin(x)/cos(x)’) 图形的输出在数学建模中，往往需要将产生的图形输出到Word文档中。通常可采用下述方法：首先，在MATLAB图形窗口中选择【File】菜单中的【Export】选项，将打开图形输出对话框，在该对话框中可以把图形以emf、bmp、jpg、pgm等格式保存。然后，再打开相应的文档，并在该文档中选择【插入】菜单中的【图片】选项插入相应的图片即可。
【数学建模学习】MatLab2画图
MatLab & 数学建模 授课：唐静波（九江学院理学院） 第二讲MatLab图形绘制功能 一、二维平面图形 基本绘图函数
hold on 命令用于在已画好的图形上添加新的图形 plot是绘制一维曲线的基本函数，但在使用此函数之前，我们需先定义曲线上每一点的x及y座标。下例可画出一条正弦曲线： x=0:0.001:10; % 0到10的1000个点的x座标 y=sin(x); % 对应的y座标 plot(x,y); % 绘图 Y=sin(10*x); plot(x,y,'r:',x,Y,'b') % 同时画两个函数
?若要改变颜色，在座标对後面加上相关字串即可： x=0:0.01:10; plot(x,sin(x),'r') 若要同时改变颜色及图线型态（Line style），也是在坐标对後面加上相关字串即可： plot(x,sin(x),'r*')
用axis([xmin,xmax,ymin,ymax])函数来调整图轴的范围axis([0,6,-1.5,1]) MATLAB也可对图形加上各种注解与处理：xlabel('x轴'); % x轴注解 ylabel('y轴'); % y轴注解 title('余弦函数'); % 图形标题
legend('y = cos(x)'); % 图形注解 gtext('y = cos(x)'); % 图形注解 ,用鼠标定位注解位置 grid on; % 显示格线 fplot 的指令可以用来自动的画一个已定义的函数分布图，而无须产生绘图所须 要的一组数据做为变数。其语法为fplot('fun',[xmin xmax ymin ymax])，其中 fun 为一已定义的函数名称，例如 sin , cos 等等；而 xmin , xmax , ymin , ymax 则是设定绘图横轴及纵轴的下限及上限。 以下的例子是将一函数 f(x)=sin(x)/x 在-20
MATLAB数学建模算法学习笔记
关于MATLAB的数学建模算法学习笔记 目录 线性规划中应用： (3) 非线性规划： (3) 指派问题;投资问题：（0-1问题） (3) 1）应用fmincon命令语句 (3) 2）应用指令函数：bintprog (5) 重新整理矩阵类型 (6) 1）应用reshape (6) 2）应用命令：nonzeros (7) 非线性的最小值得求法：含有一个变量时，应用命令:fminsearch(@fun,x0) (7) 含有多个变量时用：fminunc() (7) 求解非线性多变量等式应用命令fsolve (8) 二次规划问题应用：quadprog (8) 把有条件的问题转化成无条件问题。罚函数法：fminunc (9) 在Matlab中求解极值问题函数有： (9) 1）fminbnd (9) 1：在Matlab中求解距离的函数为：dist (9) 最小生成树 (9) prim算法 (10) Find函数的应用 (10) 关于图论的Matlab工具箱相关命令 (10) 这些命令基本上都用到稀疏阵,产生稀疏阵用sparse命令 (10) 查看网图用view (11) 积分命令quadl (11) Matlab插值工具箱 (11) 一维插值：interp1 (11) 二维插值： (11) 插值接点为网格节点：interp2 (11) 插值节点为散乱节点：griddata (11) 最小二乘法 (11) 2）应用lsqlin命令语句 (12) 三次样条差 (12) 积分函数命令：quadl (13) 同一组数据用不同插值方法效果比较线性插值、三次样条插值 (13) 参数估计 (14) 1）非线性最小拟合 (14) 命令：lsqcurvefit解决非线性拟合问题。 (14) 2）线性最小二乘法 (15) 解微分方程 (16)
数学建模培训作业(MATLAB编程部分)
数学建模培训作业 （MATLAB 编程部分） 1． 请使用switch 语句将百分制的学生成绩转换为五级制的成绩输出。 2． 猜数游戏：首先由计算机随机产生一个 [1,100] 之间的一个整数，然后由用户猜测所产生的这个数。根据用户猜测的情况给出不同的提示，如果猜测的数大于产生的数，则显示 “High” ，小于则显示 “ Low ” ，等于则显示 “You won ！”，同时退出游戏。用户最多有 7 次机会。 3． 使用for 结构计算1+2+3+…+100。 4． 设计一个九九乘法表。 5． 使用while 结构计算1+2+3+…+100。 6． 求1!+2!+ …+10!的值。 7． 编程生成三对角矩阵。 1 10000011100000001110000000111000000011100 000001110000000111000000011100 1 1 轾犏犏犏犏犏犏犏 犏犏 犏犏犏犏犏犏犏犏臌 8． 计算分段函数的值，要求根据不同的x 输入，给出相应的结果。 223135 x x y x x ì?-???+??=í?-????+?? 110011x x x x ?-< ?> 9． 已知1 1111 1(1)435721 n n p -?+-++-- ，编程求 的近似值。 10． 输入下面的矩阵
12345678 910111213141516A 轾犏犏犏=犏犏 犏臌 编程求该矩阵的对角线元素之和，并找出最大和最小元素的值以及其所在的行、列号。 11． 求水仙花数。如果一个三位数的个位数、十位数和百位数的立方和等于该数自身， 则称该数为水仙花数。编一程序求出所有的水仙花数。 12． 给定两个实数a 、b 和一个正整数n ，计算()k a b +和()k a b -，其中n k ,,2,1 。 13． 编写函数，生成一个1!，2!，…，n!的阶乘表。 14． 编一个函数统计字符串中单词的个数。 15． 求n 阶勒让德多项式的值，其递归公式为： (,)((2* 1)**(1,) (p n x n x p n x n p n n =----- (0,) 1;(1,p x p x x == 16． 编写一个判断任意输入正整数是否为素数的函数文件，并在命令窗口调用。 17． 编写一个万年历计算程序，当输入年月日后，能够计算出该日是星期几。
MATLAB培训策划书
MATLAB培训策划书 1、活动背景： 我们都知道在数学建模中，不仅仅是数学模型的建立好就行，还要面对大量数据的各种处理，要涉及许多多元统计分析知识，是数学建模的难点。这就需要数学软件来帮助解决这个问题，数学软件的引入，大大提高了同学们的数据处理能力，所以在此特举办数学软件MATLAB的应用培训。 2、活动目的： 目前数学建模逐渐和计算机技术一样成为任何一门科学发展过程中的必备工具，随着数学软件在数学建模中的应用不断加强，计算机已成为一种强大的建模工具，在我们的数学建模过程中起着越来越重要的作用。因此，本次培训活动的目的在于使大一新生能熟悉并掌握常用的数学软件，为他们在今后的数学建模中打下良好的技术基础。 3、活动主题： 东华理工大学信息工程学院计算机与通信协会MATLAB 应用软件培训 4、活动对象： 全体计算机与通信协会的会员 5、活动时间： 星期五晚上（如遇特殊情况，另行通知），共五周课。 6、活动地点： 信工楼多媒体阶梯教室 7、活动说明： 1、在经过大一前段时间的学习以及高中的数学基础知识积累，同学们对数学建模有了或多或少的了解。因此，此次MATLAB 应用软件的培训能够促使会员们更好地运用数学知识解决问题; 2、本次培训活动以老师及学长为会员们讲解数学软件MATLAB 的相关应用知识。 8、活动准备： 1、由计通协会长邀请能够熟练使用数学软件的学长; 2、提前做好MATLAB 培训的相关PPT，为老师及学长的讲 解做准备; 3、提前向学院申请培训地点（多媒体阶梯教室）;
4、以上时间、地点只是初步拟定的，确切时间由本协会干事 与会员们协商解决; 5、由干事以发送短信的方式通知会员培训。 9、活动流程： 1、在每次培训课之前干事做好会员签到工作; 2、由邀请的老师及学长介绍本学期的培训安排并开始培训; 3、讲课时，做好讲课记录; 4、讲课结束后，写出工作成果总结. 10、活动经费： 活动宣传 25 元 培训材料 80 元 授课老师培训费 200元 总计 305元 计算机与通信协会 2012年11月8日
数模中常用的matlab语言和方法
数模培训 一、矩阵的输入 二、循环语句 三、M文件及M函数 四、优化问题 五、回归分析 一矩阵的输入 常用的矩阵输入方法有三种: 1.直接输入法 2 下载法 3. 编程输入法 下面以例解释下载法 例1．现要根据瑞士地图计算其国土面积。于是对地图作如下的测量：以西东方向为横轴，以南北方向为纵轴。（选适当的点为原点）将国土最西到最东边界在x轴上的区间划取足够多的分点x i，在每个分点处可测出南北边界点的对应坐标y1，y2。用这样的方法得到下表 41288km2）。 解题思路：数据实际上表示了两条曲线，实际上我们要求由两曲线所围成的图形的面积。解此问题的方法是数值积分的方法。具体解时我们遇到两个问题：1。数据如何输入；2。没有现成的命令可用。 解：对于第一个问题，我们可把数据考备成M文件（或纯文本文件）。 然后，利用数据绘制平面图形。键入 load mianji.txt a=mianji; A=[a(1:3,1:9),a(4:6,1:9),a(7:9,1:9)]'; plot(A(:,1),A(:,2),'r',A(:,1),A(:,3),'g')
接下来可以计算面积。键入： a1=trapz(A(:,1)*40/18,A(:,2)*40/18); a2=trapz(A(:,1)*40/18,A(:,3)*40/18); d=a2-a1 d = 4.2414e+004 至此，问题可以说得到了解决。之所以说还有问题，是我们觉得误差较大。但计算方法的理论给了我们更精确计算方法。只是MA TLAB 没有相应的命令。想得到更理想的结果，我们可以自己设计解决问题的方法。（可以编写辛普森数值计算公式的程序，或用拟合的方法求出被积函数，再利用MA TLAB 的命令quad,quad8） 2. 下面以例解释编程输入法 若矩阵的各元素之间存有某种关系，那么可利用存有的关系，通过编写程序来建立矩阵；还可以利用MATLAB 的矩阵命令来建立矩阵。 例：输入矩阵)1 1( 5 5-+=?j i H for i=1:5 % 变量i 的步长为1，从1变到5 for j=1:5 % 变量j 的步长为1,从1变到5 H(i,j)=1/(i+j-1); % 将值1/(i+j-1)赋予矩阵H 的第i 行第j 列的元 素 end end 循环语句 MATLAB 也有控制流语句,用于控制程序的流程.主要有for 循环、while 循环、if 和break 三种控制语句.虽然语句很少,但功能很强.
数学建模matlab例题参考及练习讲课稿
数学建模m a t l a b例题参考及练习
数学实验与数学建模 实验报告 学院： 专业班级： 姓名： 学号： 完成时间：年月日
承 诺 书 本人承诺所呈交的数学实验与数学建模作业都是本人 通过学习自行进行编程独立完成，所有结果都通过上机验 证，无转载或抄袭他人，也未经他人转载或抄袭。若承诺 不实，本人愿意承担一切责任。 承诺人： 年 月 日 数学实验学习体会 （每个人必须要写字数1200字以上，占总成绩的20%） 练习1 一元函数的图形 1． 画出x y arcsin =的图象. 2． 画出x y sec =在],0[π之间的图象. 3． 在同一坐标系中画出x y = ，2x y =，3x y =，3x y =，x y =的图象. 4． 画出3232)1()1()(x x x f ++-=的图象，并根据图象特点指出函数)(x f 的 奇偶性. 5． 画出)2ln(1++=x y 及其反函数的图象.
6． 画出321+=x y 及其反函数的图象. 练习2 函数极限 1． 计算下列函数的极限. (1)x x x 4cos 12sin 1lim 4-+π→. 程序： sym x ; f=(1+sin(2*x))/(1-cos(4*x)); limit(f,x,pi/4) 运行结果： lx21 ans = 1 (2). 程序： sym x ; f=(1+cos(x))^(3*sec(x)); limit(f,x,pi/2) 运行结果： lx22 ans = exp(3) (3)22)2(sin ln lim x x x -ππ→. 程序： sym x ; f=log(sin(x))/(pi-2*x)^2; limit(f,x,pi/2) 运行结果： lx23 ans = x x x sec 3 2 ) cos 1 ( lim + π →
MATLAB数学建模14个范例
1.整数规划的蒙特卡洛解法2015-06-10 (2) 2. 罚函数法 2015-06-11 (3) 3. 层次分析 2015-06-12 (4) 4. 粒子群优化算法的寻优算法--非线性函数极值寻优 2015-06-13 (5) 5有约束函数极值APSO寻优 2015-06-14 (12) 6.模拟退火算法 TSP问题2015-06-15 (17) 7. 右端步连续微分方程求解2015-06-16 (19) 8. 多元方差分析 2015-06-17 (22) 9. 基于MIV的神经网络变量筛选 2015-06-18 (25) 10. RBF网络的回归--非线性函数回归的实现 2015-06-19 (29) 11. 极限学习机在回归拟合中的应用 2015-06-20 (32) 12. 极限学习机在分类中的应用 2015-06-21 (34) 13. 基于PSO改进策略 2015-06-22 (37) 14. 神经网络遗传算法函数极值寻优 2015-06-23 (46) 1.
1.整数规划的蒙特卡洛解法2015-06-10 已知非线性整数规划为： ????? ?? ??≤++≤++≤++++≤++++=≤≤-----++++=200 5200 62800622400)5,....,1(9902328243max 54233 21 64321543215 43212 524232221x x x x x x x x x x x x x x x x i x x x x x x x x x x x z i 如果用显枚举试探，共需要计算100^5=10^10个点，其计算量非常大。然而应用蒙特卡洛 去随机模拟计算10^6个点，便可以找到满意解，那么这种方法的可信度究竟怎么样呢？ 下面就分析随机采样10^6个点计算时，应用概率理论估计下可信度。 不是一般性，假设一个整数规划的最优点不是孤立的奇点。 假设目标函数落在高值区的概率分别为0.01，0.00001，则当计算10^6个点后，有任一个点落在高值区的概率分别为： 1-0.99^1000000=0.99...99(100多位) 1-0.99999^1000000=0.999954602 解 （1）首先编写M 文件 mengte.m 定义目标函数f 和约束向量g,程序如下： function [f,g]=mengte(x); f=x(1)^2+x(2)^2+3*x(3)^2+4*x(4)^2+2*x(5)-8*x(1)-2*x(2)-3*x(3)-... x(4)-2*x(5); g=[sum(x)-400 x(1)+2*x(2)+2*x(3)+x(4)+6*x(5)-800 2*x(1)+x(2)+6*x(3)-200 x(3)*x(3)+x(4)+5*x(5)-200]; （2）编写M 文件mainint.m 如下求问题的解： rand('state',sum(clock)); p0=0; tic for i=1:10^5 x=99*rand(5,1); x1=floor(x);%向下取整 x2=ceil(x);%向上取整 [f,g]=mengte(x1); if sum(g<=0)==4 if p0<=f x0=x1; p0=f; end end [f,g]=mengte(x2); if sum(g<=0)==4 if p0<=f
相关文档
- MATLAB数学建模算法学习笔记
- 建模培训 matlab教程之二图像及其补充
- 数学建模培训-matlab基础讲义剖析
- 2019数学建模培训matlab基本使用
- 第一讲Matlab基本数值计算(最新整理)
- 2012暑假建模Matlab培训
- Matlab建模仿真
- 数学建模matlab课件学习
- MATLAB培训策划书
- 第一讲Matlab基本数值计算(可编辑修改word版)
- 2016全国大学生数学建模竞赛培训(MATLAB)
- 2019精品福建农林大学数学建模培训演示数学
- 数学建模MATLAB培训2
- 国赛培训：MatLab基础及在数学建模中的应用
- Matlab建模培训讲义
- Matlab建模仿真
- 应用MATLAB建模实例
- 数学建模培训m文件
- 数学建模培训作业(MATLAB编程部分)
- 数学建模培训matlab基础讲义