数学建模实验作业参考例子

例1.1 求解线性方程组

命令如下：

>> a=[2,3,-1;8,2,3;45,3,9];

>> b=[2;4;23];

>> x=inv(a)*b

x =

0.5531

0.2051

-0.2784

例1.2 绘制正弦曲线和余弦曲线。

命令如下：

>> x=[0:0.5:360]*pi/180;

>> plot(x,sin(x),x,cos(x)) ;

例1.3 输入10个学生的成绩并对成绩按升序排序。

＞＞g=[45 56 34 24 76 89 32];

＞＞f=sort(g)

f =

24 32 34 45 56 76 89

例1.4 设有常微分方程初值问题，试求其数值解，并与精确解相比较。

(1)建立函数文件funt.m：

function yp=funt(t,y)

yp=(y^2-t-2)/4/(t+1);

(2)求解微分方程：

t0=0;tf=10;y0=2;

[t,y]=ode23('funt',[t0,tf],y0);

y1=sqrt(t+1)+1;

t'

例2.1 计算表达式的值，并将结果赋给变量x，然后显示出结果。

在MA TLAB命令窗口输入命令：

>>x=(5+cos(47*pi/180))/(1+sqrt(7)-2*i) %计算表达式的值

2.2 A=[1 2 3;4 5 6;7 8 9]

A =

1 2 3

4 5 6

7 8 9

例2.3 分别建立3×3、3×2和与矩阵A同样大小的零矩阵。

(1)建立一个3×3零矩阵：zeros(3)

>> zeros(3)

ans =

0 0 0

0 0 0

0 0 0

)建立一个3×2零矩阵：zeros(3,2)

>> zeros(3,2)

ans =

0 0

0 0

0 0

(3)建立与矩阵A同样大小零矩阵：zeros(size(A))

>> A=[2 4 7 5;7 8 9 10];zeros(size(A))

ans =

0 0 0 0

0 0 0 0

例2.4 产生5阶随机方阵A，其元素为[10,90]区间的随机整数，然后判断A的元素是否能被3整除。

(1) 生成5阶随机方阵A。

>> A=fix((90-10+1)*rand(5)+10)

A =

86 71 59 42 14

28 46 74 85 38

59 11 84 84 75

49 76 69 43 10

82 46 24 82 21

(2) 判断A的元素是否可以被3整除。

>> P=rem(A,3)==0

P =

0 0 0 1 0

0 0 0 0 0

0 0 1 1 1

0 0 1 0 0

0 0 1 0 1

例3.1 建立一个命令文件将变量a,b的值互换，然后运行该命令文件。

首先建立命令文件并以文件名exch.m存盘：

clear;

a=1:10; b=[11,12,13,14;15,16,17,18];

c=a;a=b;b=c;

a

b

ant

a =

11 12 13 14

15 16 17 18

b =

1 2 3 4 5 6 7 8 9 10

例3.2 建立一个函数文件将变量a,b的值互换，然后在命令窗口调用该函数文件。

首先建立函数文件fexch.m：

function [a,b]=fexch(a,b) % exch

c=a;a=b;b=c;

clear;

x=1:10; y=[11,12,13,14;15,16,17,18];

[x,y]=fexch(x,y)

例3.5 输入一个字符，若为大写字母，则输出其后继字符，若为小写字母，则输出其前导字符，若为数字字符则输出其对应的数值，若为其他字符则原样输出。

c=input('请输入一个字符','s');

if c>='A' & c<='Z'

disp(setstr(abs(c)+1));

elseif c>='a'& c<='z'

disp(setstr(abs(c)-1));

elseif c>='0'& c<='9'

disp(abs(c)-abs('0'));

else

disp(c);

end

例3.6 某商场对顾客所购买的商品实行打折销售，

已知打折标准为：

商品价格小于200元时，折扣率为0;

商品价格大于或等于200元且小于500时，折扣率为0.03;

商品价格大于或等于500且小于1000元时，折扣率为0.05;

商品价格大于或等于1000且小于2500元时，折扣率为0.08;

商品价格大于或等于2500且小于5000元时，折扣率为0.1;

商品价格大于或等于5000元时，折扣率为0.14;

求所售商品的实际销售价格。

＞＞price=input('请输入商品价格');

switch fix(price/100)

case {0,1}

rate=0;

case {2,3,4}

rate=3/100;

case num2cell(5:9)

rate=5/100;

case num2cell(10:24)

rate=8/100;

case num2cell(25:49)

rate=10/100;

otherwise

rate=14/100;

end

price=price*(1-rate)

请输入商品价格890

price =

845.5000

例3.7 矩阵乘法运算要求两矩阵的维数相容，否则会出错。先求两矩阵的乘积，若出错，则自动转去求两矩阵的点乘。

程序如下：

A=[1,2,3;4,5,6]; B=[7,8,9;10,11,12];

try

C=A*B;

catch

C=A.*B;

end

C

lasterr %显示出错原因

C =

7 16 27

40 55 72

ans =

Error using ==> *

Inner matrix dimensions must agree.

例3.8 已知求y 的表达式，

当n=100时，求y 的值。

程序如下

y=0;n=100;

for i=1:n

y=y+1/i/i;

end

y

补充知识： 在实际MATLAB 编程中，为提高程序的执行速度，常用向量运算来代替循环操作：

n=100; i=1:n;

f=1./i.^2;

y=sum(f)

例3.9 求定积分

程序如下：

a=0;b=3*pi;n=1000; h=(b-a)/n;

x=a:h:b; f=exp(-0.5*x).*sin(x+pi/6);

for i=1:n

s(i)= (f(i)+f(i+1))*h/2;

end

s=sum(s) 补充知识：

事实上，MATLAB 提供了有关数值积分的标准函数，实际应用中可直接调用这些函数求数值积分

例3.10 已知5个学生4门功课的成绩，求每名学生的总成绩。

程序如下：

1320()sin()6

x e x dx ππ-+?

s=0;

a=[65,76,56,78;98,83,74,85;76,67,78,79;98,58,42,73;67,89,76,87]

for k=a

s=s+k;

end

disp(s');

结果为：

a =

65 76 56 78

98 83 74 85

76 67 78 79

98 58 42 73

67 89 76 87

275 340 300 271 319 例3.11 产生一个长度为６的调和数列。

程序如下：

X=zeros(1,6);

i=1;

while i<=6

X(i)=1/i;

i=i+1;

end

format rat %使用分数表示数值

disp(X

例3.12 用while语句实现例3.8。

程序如下：

y=0; i=1;

while 1

f=1/i/i;

y=y+f;

if i==100

break;

end

i=i+1;

end

y

例3.13 用筛选法求某自然数范围内的全部素数。

程序如下：

m=input('m=');

p=2:m;

for i=2:sqrt(m)

n=find(rem(p,i)==0&p~=i);

p(n)=[];

end

p

输入一个值为：

m=100

p =

Columns 1 through 5

2 3 5 7 11

Columns 6 through 10

13 17 19 23 29

Columns 11 through 15

31 37 41 43 47

Columns 16 through 20

53 59 61 67 71

Columns 21 through 25

73 79 83 89 97

例3.14 编写函数文件求半径为r的圆的面积和周长。函数文件如下：function [s,p]=fcircle(r)

%CIRCLE calculate the area and perimeter of a circle of radii r

%r 圆半径

%s 圆面积

%p 圆周长

%2001年7月30日编

s=pi*r*r;

p=2*pi*r;

将以上函数文件以文件名fcircle.m存入c:\matlabr11\work下，然后在MA TLAB命令窗口调用该函数：

[s,p]=fcircle(10)

s =

13823/44

p =

7100/113

例3.15 利用函数文件，实现直角坐标(x,y)与极坐标(ρ,θ)之间的转换。

函数文件tran.m：

function [rho,theta]=tran(x,y)

rho=sqrt(x*x+y*y);

theta=atan(y/x);

调用tran.m的命令文件main1.m：

x=input('Please input x=:');

y=input('Please input y=:');

[rho,the]=tran(x,y);

rho

the

例3.17 nargin用法示例。

函数文件examp.m：

function fout=examp(a,b,c)

if nargin==1

fout=a;

elseif nargin==2

fout=a+b;

elseif nargin==3

fout=(a*b*c)/2;

end

命令文件mydemo.m：

x=[1:3];y=[1;2;3];

examp(x)

examp(x,y')

examp(x,y,3)

例3.18 全局变量应用示例。

先建立函数文件wadd.m，该函数将输入的参数加权相加。

function f=wadd(x,y)

global ALPHA BETA

f=ALPHA*x+BETA*y;

在命令窗口中输入：

global ALPHA BETA

ALPHA=1;

BETA=2;

s=wadd(1,2)