成績:

Matlab基礎及科學計算

姓名：黃曉陽

班號：HC001501

學號：2()1530()005

2345

4321

1、用Uatlab語句創建矩陣A,A=，如果給出“A(5,6)=你的學號”命令

6437

8946

將得到什么結果？

答:A方陣將擴充成5X6的矩陣,在第5行第6列處填充2015300005,多出來的其他未賦值處

填充0,原位置處為初始值.

源代碼：

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

A(5,6)=2015300005

2、川magic命令生成一個六階的魔方矩陣，并將其全部偶數行提取出來。

答:提取的結果為[8,30,4,32,28,36,6,2,34,26,22,12,10,14,18,24,20,16]

源代碼：

Magic=reshape(magic⑹，[1,36]);

result=(];

fori=l:length(Magic)

ifrem(Magic(i),2)==0

result(length(result)+1)=Magic(i);

end

end

result

3、輸入A為3x3的魔方陣，B為3x3的單位陣，由小矩陣組成3x6的大矩陣C和6x3

的大矩陣D,將D矩陣的最后1列構成向量E。

源代碼：

A=magic(3)

B=ones(3,3)

C=zeros(3,6);

D=zeros(6,3);

fori=l:3

C(iz:)=[A(iz:),B(iz:)];

end

C

D=[A;B]

E=D(:,3)

123

4、輸入矩陣A=456，使用全下標方式取出元素“3”，使用單下標方式取出元素

789

13

“8”，取出后2行子矩陣塊,使用邏輯矩陣方式取出

79

答:詳細操作見代碼.

源代碼：

A=[l,2,3;4,5,6;7,8,9]

three=A(1,3)

three_=A(7)

last_two_line=A(2:3,:)

A.*(rem(A,2)~=0).*(A~=5)

50

5、用數值方法求出5=,2'=1+2+4+…+2”+25。，試不采用循環的形式給出和式的數

r=0

值解。

答:5=2251799813685247

源代碼：

symsk

symsum(2Ak,k,0,50)

6、試用循環結構找出“學號后三位”以下的所有質數。

答:[2,3,5]

源代碼：

fori=l:005

if(isprime(i))

i

end

end

7、Fibonacci數列定義為4=4-+%一2/=3,4…為一個數值增長很快的數歹力

試編寫程序尋找Fibonacci數列中第一個大于10000的元素

答:10946

源代碼：

i=l;

whilefibonacci(i)<10000

i=i+l;

and

fibonacci(i)

8、使用while循環計算滿足。(2燈-1)>學號后五位的最小，〃值及相應的積。

答:mmin=3乘積為：15

源代碼：

clear

clc

i=l;

c=l;

whilec*(2*i-l)<005

c=c*(2*i-l);

i=i+l;

end

i

c*(2*1-1)

9、下面給出一個迭代模型

%|=1+券一14%

=0.3玉

寫出求解該模型的M函數，如果取迭代初值為毛=0,%=0,那么請進行30000次迭代

求出一組x和y向量，然后在所有的'和以坐標處點亮一個點(注意不要連線)，最后繪制

出所需的圖形。(注：繪制出的圖形又稱為，它將迭代出來的隨機點吸引到一起，最后得出

貌似連貫的引力線圖)

T9(30000)

functionT9(n)

x=zeros(size(1,n))；

y=zeros(size(1,n));

for1=1:n-l

x(i+l)=l+y(i)-1.4*x(i)A2;

y(i+l)=0.3*x(i);

end

plot(x,y,1.1)

title('Henon引力線圖，)

10、繪制出下面的圖形sin（；）,其中由于該函數在1=0附近震蕩劇烈，因此,

嘗試在在0附近進行調整步長。繪制并對比兩種時間序列下的結果。第一個是采用0.02的

時間間隔，另外一個在/w[-0.2,0.2]附近將時間間隔設置為0.001。

*p，0.02step?0.001

源代碼：

t=[-l:0.02:l];

subplot(1,2,1)

plot(t,oin(1./t))

title(*step=0.021)

t=-0.2:0.001:0.2;

subplot(1,2,2)

plot(t,sin(1./t))

title(*step=0.0011)

2*X2+1x>=\

11、繪出函數》二,0-Icxvl的圖像。

-xx<=-1

源代碼：

a=[-10:0.01:10];

plot(a(a>=l),2.*a(a>=l).A2+l,a(a<l&a>-l),0,a(a<=-l),-a(a<=-l).A3)

sin(86/3)

12、試繪制出p-的極坐標曲線，選擇，變量的范圍為(0,4%)和

學號尾號-cos?(3。/2)

(0,6幻，觀察極坐標曲線的情況

源代碼：

subplot(1,2,1)

theta=[0:0.01:4*pi];

rho=sin(8*theta/3)./(005-cos(3*theta/2).A2);

polar(theta,rho)

subplot(1,2,2)

theta=[0:0.01:6*pi];

rho=sin(8*theta/3)./(005-cos(3*theta/2).A2);

polar(theta,rho)

13、請分別繪制出p和sin。〉，）的三維圖和等高線

10

源代碼：

[x,y]=meshgrid(-10:0.1:10);

subplot(2,2,1)

meoh(x,y,x.*y)

subplot(2,2,2)

contour(x,y,x.*y,80)

subplot(2,2,3)

mesh(x,y,sin(x.*y))

subplot(2,2,4)

contour(x,y,sin(x.*y))

14、假設A為四階希爾伯特矩陣，分別找出由矩陣A的第1、2、3、4條對角線上的元

素組成的列矢量。

-1

'1/2-

1/3~l/3~

答:diag_l=diag_2=1/4diag_3=diag_4=[l/4]

1/5_l/5_

1/6

1/7_

源代碼：

A=hilb(4)

diag(A,0)

diag(A,1)

diag(A,2)

diag(A,3)

123

15、求矩陣A=3I2的轉置矩陣、逆矩陣、矩陣的秩、矩陣的行列式值、矩陣的

_234

三次暴、矩陣的無窮范數，矩陣的特征值和特征向量。

-123~-132~909614「

答:A=312A'=2I3A3=8790132

_234_324132141207

-

2/31/31/3

inv(A)=8/3-2/1§7/3

7/31/3-5/3

7.290800

特征值矩陣為：（）-0.71680

00-0.5739

0.49720.17830.1051

對應的特征向量矩陣為：0.4691-0.8628-0.8520

0.72990.47310.5129

rank(A)=3|A|=3||A||-二9

源代碼：

A=[l,2,3;3,l,2;2,3,4]

zhuanzhi=A'

ni=inv(A)

zhi=rank(A)

hanglieshi=det(A)

three_mi=AA3

wuqiongfanshu=norm(A,inf)

[tezhengxiangliang,tezhengzhi]=eig(A);

tezhengzhi

tezhengxiangliang

x+y+z=10

16、求方程組＜)，+z=0的解

2x-y-z=-4

x=2

答：，y=5

z=3

源代碼：

[x,y,z]=solve(*x+y+z=10,x-y+z=O,2*x-y-z=-4','x,yzz*)

x}+x2-3X3-X4=1

17、求非齊次線性方程組卜%-3&+4七=4的通解

N+5X2-9X3-8X4=0

0"3/2--3/4

03/27/4

答：x=+G+C2

-8/1510

3/5_01

源代碼：

A=[1,1,-3,-1;3,-1,-3,4;1,5,-9,-8];

b=[1;4;0];

B=[A,b];

teShuJie=A\b

tongjie=null(A,'r1)

18、求解下面的Lyapunov方程，并檢驗所得解的精度

0.31960.2360-0.2096

答：X:0.2624-0.18200.2844，將結果帶入到所給矩陣運算式左端,可得右端為

-0.11780.5393-0.1797

1.000000000000003.000000000000002.00000000000000

3.000000000000004.00000000000000l.OOOOODOOOOOOOO，與上式右端相等，每

5.000000000000002.000000000000001.00000000000000

個元素的誤差在10“4的量級.

源代碼:

A=[l,2,3;4,5,6;7,8,0];

B=[2Z3,6;3,5Z2;4Z3Z1];

C=[l/3,2;3,4/l;5,2/L];

X=lyap(A,B,-C)

C_=A*X+X*B

C

A=8,.24k=0fl

19、計算zL，y——

.k\占Ax(A+l)

x

e'I

2

G:S,=A:(---)S,=l

X~X

源代碼：

symskx

S2=symsum(1/(k*(k-1)),k,1,Inf)

A

SI=symsum(x(2*k)/factorial(k)zk,1,Inf)

20、對函數/=求函數的的一階導和二階導。

夕”

xcos(x)sin(x)xsin(x)

答：〈

xsin(x)cos(x)xcos(x)sin(x)

4

源代碼：

symsx

yijiedao=diff(x*sin(x)/exp(2*x),1)

erjiedao=diff(x*sin(x)/exp(2*x),2)

21、人造地球衛星的軌跡可視為平面上的橢圓。我國第一顆人造地球衛星近地點距離地

球表面439km,遠地點距離地球平面2384km,地球半徑為6371km,求該衛星的軌跡長度。

人造地球衛星的軌跡可用橢圓的參數方程來表示，即

.r=?cosr

(0<t<2兀;a,b>0)

y=bsh\t

其中，a=6371+2384=8755(k⑼。=6371+439=6810(癡)。

由曲線積分可知，衛星估計長度L為積分L=4亢°2si/r+b28s2dt的值。

答:該衛星軌跡長度為4.9095e+04(km).

源代碼：

%T21(2384+6371,439+6371,10000)

functionT21(a,b,n)

len=0;

step=pi/2/n;

fortheta=0:step:pi/2

len=len+Len(a*cos(theta),a*cos(theta+stap),...

b*sin(theta),b*sin(theta+step));

and

perimeter=len*4

end

functiony=Len(xlzx2,yl,y2)

y=((xl-x2)A2+(yl-y2)A2)A0.5;

end

22、求二重積分J'JJ(ysinx+xsi”)dMv，]了,—+7)，+6"dt

答:S]=-98.6960S2=14.5

源代碼：

symsxy

f1=@(xzy)y.*sin(x)-x.*sin(y);

f2=@(x,y)xJ2+7.*yJ2+6;

integral2(fl,3*pi,4Tpi,pi,2*pi)

integral2(￡2,0,1,@(x)x,@(x)x+1)

23、求解微分方程的通解和初始條件為M0=0，W0=l時的特解，并采用ode45()的數

值求解方法解該方程組，繪圖并比較兩種求解方法。

公

-

小

電

力

ode<45)求桿x(t)用像

源代碼：

functionT23()

[x,y]=dsolve('Dx=y+x','Dy=2*x1,'t*);

disp(x)

disp(y)

[x,y]=dsolve('Dx=y+x*,'Dy=2*x1,'x(0)=0*,'y(0)=1*,'t*);

disp(x)

disp(y)

subplot(2,2,1)

ezplot(x)

title(*dsolve求解x(二)圖像，)

gridon

axis([0,4,0,1000])

subplot(2,2,2)

ezplot(y)

title('dsolve求解y(二)圖像，)

gridon

axis([0z4,0z1000])

g=@(t,y)Wf(t,y);

[t,x]=ode45(g,[0,4.2],[0,1]);

subplot(2,2,3)

plot(t,x(:,1));

title(?ode(45)求解x(t)圖像,)

gridon

axis([0z4,0z1000])

subplot(2,2,4)

plot(t,x(:,2))

title('ode(45)求解y(t)圖像，)

gridon

axis([0z4z0,1000])

end

functiondy=Wf(t,y)

dy=zeros(2,1);

dy(l)=y(2)+y(l);

dy(2)=2*y(l);

end

24、求解下面的最優化問題

,，牒“(xTx+0

兌％=1

答:《

x2=0

源代碼：

functionT24()

fmincon(@(x)x(1)A2-2*x(1)+x(2),[50,300]，匚,[],[],[],[],[],@nonl

con)

end

function[C,Ceq]=nonlcon(x)

C(l)=4*x(l)A2+x(2)A2-4;

C(2)=-x(l);

C(3)=-x(2);

Ceq=[];

end

25、求解/(1)=10012-4『+(1-西)2在1121]上的最小值

答:當未知量都取1時，有最小值f(X)min=0

源代碼：

functionT25()

fmincon(@(x)100*(x(2)-x(l)A2)A2+(l-x(l))A2,[0,0],[],[],[],[],[J,[：/

un)

end

function[c,ceq]=nun(x)

c(l)=-x(l)-1.2;

c(2)=x(l)-l;

c(3)=-x(2)-1.2;

c(4)=x(2)-1;

ceq=[];

end

min(X)+3/-5x3)

的線性規劃

26、計算4

3%+鼻W4

xpx2,x3>0

尤i=0

答:當未知量取值為,占=0時，有最小值-20

￡=4

源代碼：

functionT26()

fmincon((a(x)x(l)+3*x(2)-5*x(3)/[0,0,0],[],：],[],[],[],[],@nun)

end

function[c,ceq]=nun(x)

ceq=[];

c(1)=x(1)-x(2)+x(3)-8;

c(2)=3*x(l)+x(3)-4;

c(3)=-x(1);

c(4)=-x⑵;

c(5)=-x(3);

end

27、已知x=0:0.2:2,函數y=(x?-3x+學號尾號)，"sinx。對外=0:0.03:2,分別采

用最近鄰插值、后臨插值、線性插值、三次樣條插值方法進行插值，并畫出響應的圖形。

近鄰插僮E鄰插俶

源代碼：

clea