2023南京邮电大学通达学院《数学实验》MATLAB实验答案

四月维夏,六月徂暑。
勤将励勉,勿望再晨。

——赠nmy

一 声明

南京邮电大学通达学院《数学实验》MATLAB实验答案
答案更新时间:2023.04.28,修改了4.2的存疑部分。已更新完成,如无错误不在更新

为了方便核算,我在代码中单独将m定义为自变量运算或者直接以m=117代入,作业中可以直接代入,即代码中不出现m。本机版本为 MATLAB R2020b在这里插入图片描述


由于作者解答能力有限,难免有瑕疵错误之处,还请多多海涵!本答案仅供学习参考之用,请勿直接抄袭。有错漏之处,烦请指正。联系QQ:1415520898,如有问题可通过qq或者评论区留言方式交流。

二 MATLAB下载

这里引用@dew_142857博主的相关文章最新MATLAB R2020b超详细安装教程(附完整安装文件)实测有效,按照步骤一步步来即可,为方便同学下载,这里将文中所提向公众号索要的百度网盘链接放在下方
另外安装好的MATLAB约为96.6 GB ,请提前规划好磁盘空间。

链接:https://pan.baidu.com/s/1NExZ_v-QN4Xbu4Jk1C0dEA
提取码:7won

也可以在https://matlab.mathworks.com/注册一个账户,直接在线使用
在这里插入图片描述

《数学实验》练习一

1.1

在这里插入图片描述
在这里插入图片描述在这里插入图片描述

log(x)——>lnx;inf——>无穷

1.2

请添加图片描述
在这里插入图片描述

exp(x)——>eˣ;diff(y,x,n)——>y对x的n阶导函数

1.3

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

第一小问答案不要忘记+C;int——>处理定积分、不定积分

1.4

在这里插入图片描述
2020版本
在这里插入图片描述

写全应该是taylor((117/200+sin(x))*cos(x),x,‘Order’,5,‘ExpansionPoint’,0),在x=0处可省略。

2010版
在这里插入图片描述

1.5

在这里插入图片描述

在这里插入图片描述

本次随机的中间数据为:

[8226958330713791/9007199254740992, (2^(1/2)*469134536469018791^(1/2))/671088640, ((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2), (((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2), ((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2), (((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2), ((((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2), (((((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2), ((((((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2), (((((((((2^(1/2)*469134536469018791^(1/2))/671088640 + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2) + 117/100)^(1/2)]

1.6

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

本题用到的符号较多,进行下一题时使用clear清除变量

1.7

在这里插入图片描述
在这里插入图片描述
1.7.1
请添加图片描述
在这里插入图片描述
在这里插入图片描述
1.7.2
在这里插入图片描述

在这里插入图片描述
在这里插入图片描述

1.8

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

1.8.2
在这里插入图片描述
也可以使用下方代码,效果一样
在这里插入图片描述

在这里插入图片描述

1.9

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

1.10

plot是绘制二维图形,并且是x,y的表达式是已知的或者是形如y=f(x)这样确切的表达式plot函数的基本调用格式为:plot(x,y) 其中x和y为长度相同的向量,分别用于存储x坐标和y坐标数据。
ezplot是画出隐函数图形,是形如f(x,y)=0这种不能写出像y=f(x)这种函数的图形ezplot一元函数绘图函数ezplot(fun) ezplot(fun,[min,max])
fplot(y,[a,b])精确绘图

请添加图片描述

在这里插入图片描述
请添加图片描述

1.11

请添加图片描述
在这里插入图片描述
[X,Y] = meshgrid(-5:0.1:5);可以换成书上形式:
x=-5:0.1:5;y=x;
[X Y]=meshgrid(x,y);
在这里插入图片描述

《数学实验》练习二

2.1

请添加图片描述
在这里插入图片描述

第一个不动点为-0.0084
第二个不动点为119.0084

(2)先定义一个普世性的迭代方法,用M文件保存
在这里插入图片描述
在这里插入图片描述
函数收敛,只要初值不取14165^(1/2)/2+119/2 即第二个不动点,收敛值与初值的选取关系不大,总是收敛于-0.0084, 只有初值取 14165 ^(1/2)/2+119/2,迭代函数才以它为极限;
收敛值一定是不动点其一;

2.2

请添加图片描述

m=117;
syms x;
f=inline('1-2*abs(x-1/2)');%设定函数
x0=1/4;%设定初值
for i=1:1:10
plot(i,f(x0),'*');%用*作图,可以在括号内添加'MarkerSize',20放大点
x0=f(x0); %更新x0的值,x0类似于C语言的static类型变量
hold on %将各个点划在一张图上
end
hold off

在这里插入图片描述

几个图像最后都是趋于0,如果没有的话要将i的终值调大,我的后三个图的i=1:1:100;

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

2.3

该题是P76页例二
在这里插入图片描述
在这里插入图片描述

%MARTIN函数代码
function Martin(a,b,c,N) %N为迭代次数
f=@(x,y)(y-sign(x)*sqrt(abs(b*x-c)));
g=@(x)(a-x);
m=[0;0];
for n=1:N
    m(:,n+1)=[f(m(1,n),m(2,n)),g(m(1,n))];%表示矩阵m的第n+1列。冒号表示选择所有行
end
plot(m(1,:),m(2,:),'kx');
axis equal %横纵坐标采用相等单位长度
%循环迭代N次,N是预定义的数字。在循环内部,代码更新矩阵m中的值。 具体来说,该代码通过将其第一个元素设置为f(m(1,n),m(2,n)),将其第二个元素设置为g(m(1,n))来更新m的第n列。 第一行0后面的分号表示矩阵m初始化为两行N列的列向量。

m=117;
Martin(m,m,m,5000)
Martin(-m,-m,m,10000)
Martin(-m,m/1000,-m,15000)
Martin(m/1000,m/1000,0.5,20000)

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

2.4

在这里插入图片描述
在这里插入图片描述
(1)

%此小问无需在卷面作答,且每人选的数不一样,仔细看题目!!!
m=117;
syms x;
diff(subs((100*x+117)/(x^2+100),x,117^(1/3))) %对默认的变量进行一次的求导
%我取的是a=100,c=1;最后结果的绝对值应小于1才可以,否则另取
ans =
 
0

(2)

syms x;
m=117;
f=inline('(100*x+117)/(x^2+100)');
x0=10;% 任取一个初值
for i=1:20;
x0=f(x0);
fprintf('%g,%g\n',i,x0);
end

%我的运行结果
1,5.585
2,5.14893
3,4.99475
4,4.93387
5,4.90889
6,4.89849
7,4.89413
8,4.8923
9,4.89153
10,4.89121
11,4.89107
12,4.89101
13,4.89099
14,4.89098
15,4.89098
16,4.89097
17,4.89097
18,4.89097
19,4.89097
20,4.89097

(3)
根据个人体会回答

函数迭代的收敛速度与初值的选取关系不大;
迭代初值对迭代的收敛性存在影响,但是这种影响存在不确定性,没有发现可循的规律;

用自己的话改一下即可

2.5

在这里插入图片描述

syms x;
y=sin(x);
y1=taylor(sin(x),x,'Order',2);
y2=taylor(sin(x),x,'Order',4);
y3=taylor(sin(x),x,'Order',6);
fplot([y y1 y2 y3])
xlim([-3/2*pi 3/2*pi])
grid on
legend('sin(x)','approximation of sin(x) up to O(x^1)','approximation of sin(x) up to O(x^3)','approximation of sin(x) up to O(x^5)')

在这里插入图片描述

(2)

syms x;
y=sin(x);
y1=taylor(sin(x),x,'Order',8);
y2=taylor(sin(x),x,'Order',10);
y3=taylor(sin(x),x,'Order',12);
fplot([y y1 y2 y3])
xlim([-3/2*pi 3/2*pi])
grid on
legend('sin(x)','approximation of sin(x) up to O(x^7)','approximation of sin(x) up to O(x^9)','approximation of sin(x) up to O(x^(11))')

在这里插入图片描述

(3)

在这里插入图片描述

《数学实验》练习三

3.1

在这里插入图片描述

A=str2sym('[117,117-4;6-117,10-117]');%表示符号表达式
[P,D]=eig(A);
Q=inv(P);
syms n;
x=[1;2];
xn=P*(D.^n)*Q*x 


xn =
 
(339*6^n)/2 - (337*4^n)/2 - (559*0^n)/111
        2*0^n + (337*4^n)/2 - (333*6^n)/2

3.2

在这里插入图片描述

A=str2sym('[117,117-4;6-117,10-117]');
B=1/10.*A;
[P,D]=eig(B);
Q=inv(P);
syms n;
x=[1;2];
xn=P*(D.^n)*Q*x



xn =
 
(339*(3/5)^n)/2 - (337*(2/5)^n)/2 - (559*0^n)/111
        2*0^n + (337*(2/5)^n)/2 - (333*(3/5)^n)/2

3.3

在这里插入图片描述

%教材P136页原题
A=[9,5;2,6];
t=[];
for i=1:20
    x=2*rand(2,1)-1;
    t(length(t)+1,1:2)=x;
    for j=1:40
        x=A*x;
        t(length(t)+1,1:2)=x;
    end
end
plot(t(:,1),t(:,2),'*')
grid('on') 

在这里插入图片描述

(2)可以看到,迭代阵列似乎在一条通过原点的直线上。
(3)

A=[9,5;2,6]; a=[];
x=2*rand(2,1)-1; 
for i=1:20
a(i,1:2)=x;
x=A*x;
end
for i=1:20
if a(i,1)==0
else t=a(i,2)/a(i,1);
fprintf('%g,%g\n',i,t);
end
end
%结果
1,0.911983
2,0.551028
3,0.451391
4,0.418261
5,0.406586
6,0.402388
7,0.400867
8,0.400315
9,0.400115
10,0.400042
11,0.400015
12,0.400006
13,0.400002
14,0.400001
15,0.4
16,0.4
17,0.4
18,0.4
19,0.4
20,0.4

(4)
极限值是图像直线的斜率

按照自己语言组织下面任意一条

  1. 最终稳定值为迭代矩阵的特征值之一。
  2. 如果迭代矩阵有多个线性无关的特征向量对应于同一个特征值,那么最终稳定值将是这些特征向量线性组合的结果。
  3. 稳定值是迭代矩阵的特征向量,对应的特征值为1。而迭代矩阵的特征值和特征向量则可以通过特征方程来求得。

3.4

在这里插入图片描述
书P141相似题

m=117;
A=[m-1,m;1-m,-m];
p=[0.4;0.6];%选择合适初始向量,要求和为1
[P,D]=eig(A)%P每列是特征向量,D主对角线元素是特征值
for i=1:20
    p(:,i+1)=A*p(:,i);
end
fprintf('%2f,%2f\n',p)

还可以使用下面的方法求稳定值

m=117;
A=[m,1/4-m;m-3/4,1-m];
x0=[0.4;0.6];
n=10000;
y = A^n * x0

结果

%A=[m,6-m;m-2,8-m]
%A=[m,1/4-m;m-3/4,1-m]
%A=[m-1,m;1-m,-m]

(4)ps:本题较难,可适当放弃

在线性映射迭代中,迭代矩阵的稳定性取决于其特征值的大小和分布。特征值是矩阵的一个重要性质,它描述了矩阵在线性变换下的变化情况。

如果迭代矩阵的所有特征值的绝对值都小于1,那么迭代矩阵就是稳定的,每次迭代后矩阵的元素值都会趋近于一个稳定值。

但是,如果迭代矩阵存在特征值的绝对值大于等于1,那么迭代矩阵就是不稳定的。这种情况下,每次迭代后矩阵的元素值都会趋近于无穷大或无穷小,从而导致迭代结果失效。

另外,如果迭代矩阵存在多个特征值相同的情况,那么迭代矩阵也可能不稳定。这种情况下,迭代矩阵的特征向量可能会出现非常大的幅度波动,从而导致迭代结果不可靠。

因此,对于二维矩阵的线性映射迭代,需要对迭代矩阵的特征值进行分析,以确定其稳定性。如果迭代矩阵不稳定,需要采取一些措施,如调整迭代步长或使用更稳定的迭代算法,以确保迭代结果的可靠性。

3.5

在这里插入图片描述

%如果默认b>a
>> I=0;
>> m=[];
>> n=1000;
>> for a=1:n
for c=a+1:n
b=sqrt(c^2-a^2);
if(b==floor(b))&(b>a)&(c==b+2)
I=I+1;m(:,I)=[a,b,c];
end
end
end
>> m

m =

  1176     8    10    12    14    16    18    20    22    24    26    28    30    32    34    36    38
     8    15    24    35    48    63    80    99   120   143   168   195   224   255   288   323   360
    10    17    26    37    50    65    82   101   122   145   170   197   226   257   290   325   362

  182940    42    44    46    48    50    52    54    56    58    60    62
   399   440   483   528   575   624   675   728   783   840   899   960
   401   442   485   530   577   626   677   730   785   842   901   962

>>

公式:a=2m b=m^2-1 c=m^2+1(m>2,m为整数);
即:
{a,b,c}={(2u)^2,(u^2-1)^2,(u^2+1)^2}

上课时默认b>a,下面给出a、b关系不确定是时的代码,无需写在试卷上

abc0=zeros(1000,3);
k=0;
for c=3:1000
b=c-2;
a=sqrt(c^2-b^2);
if(mod(a,1)==0)
k=k+1;
abc0(k,:)=[a b c];
end
end
abc=abc0(1:k,:);
fprintf('所有勾股数 a b c=\n')
disp(abc)

3.6

在这里插入图片描述

for k=1:200 
	for b=1:999
		a=sqrt((b+k)^2-b^2);
		if((a==floor(a))&gcd(gcd(a,b),(b+k))==1)fprintf('%i,',k);
			break;
		end
	end
end

1,2,8,9,18,25,32,49,50,72,81,98,121,128,162,169,200

k为完全平方数或者完全平方数的二倍
预测k在[200,300]之间有200,225,242,288,289

《数学实验》练习四

4.1

在这里插入图片描述

% 方法一:通过法方程组求解
d0=9;
x=[1.5,1.8,2.4,2.8,3.4,3.7,4.2,4.7,5.3];
y=[8.9,10.1,12.4,14.3,16.2,17.8,19.6,22.0,24.1];
d1=sum(x);d2=sum(x.^2);b1=sum(y);b2=sum(y.*x);
A=[d0,d1;d1,d2];B=[b1;b2];
u=A\B;
a0=u(1)
a1=u(2)
error=sum((y-(a0+a1.*x)).^2)

a0 =

    2.8304


a1 =

    4.0244


error =

    0.2409
%方法二:直接求解
x=[1.5,1.8,2.4,2.8,3.4,3.7,4.2,4.7,5.3];
y=[8.9,10.1,12.4,14.3,16.2,17.8,19.6,22.0,24.1];
P=polyfit(x,y,1)


P =

    4.0244    2.8304


error=sum((y-(2.8304+4.0244.*x)).^2)%误差


error =

    0.2409

在这里插入图片描述

4.2

在这里插入图片描述
(1)

%我的学号尾数是7;故数据是到1920年,对应人口是106.5,第14,
%则t2=t(14),x2=x(14)
%这段代码是用来进行数据拟合的,其中变量t和x分别代表时间和数据点。代码用log函数将数据点转换成线性形式,然后使用线性回归来拟合两个数据点的斜率和截距,最后用指数函数求出x0和k,从而得到新的函数曲线。代码中的error表示新的函数曲线与原数据点的误差平方和
t=1790:10:1980;
x=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76,92,106.5,123.2,131.7,150.7,179.3,204.0,226.5];
t1=t(1);x1=x(1);
t2=t(14);x2=x(14); %此步根据学号不同而不同
A=[1,t1;1,t2];
b=[log(x1);log(x2)];
u=A\b;
x0=exp(u(1))
k=u(2)
error=sum((x0*exp(k*t)-x).^2)


x0 =

   6.5242e-20


k =

    0.0254


error =

   1.2278e+05

(2)

t=1790:10:1920;
x=[3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76,92,106.5];
%我的数据是到1920,所以上面的数据是截到1920年对应的106.5
y=log(x); 
m=length(t);
A=[m,sum(t);sum(t),sum(t.^2)]; 
b=[sum(y);y*t'];
u=A\b;
x0=exp(u(1))
k=u(2)
error=sum((x0*exp(k*t)-x).^2)


x0 =

   2.7207e-20


k =

    0.0260


error =

  681.9588

4.3

在这里插入图片描述

x=1:26;
y=[1807,2001,2158,2305,2422,2601,2753,2914,3106,3303,3460,3638,3799,3971,4125,4280,4409,4560,4698,4805,4884,4948,5013,5086,5124,5163];
 a=[6000,2,0.01];
f=@(a,x)a(1)./(1+a(2)*exp(-a(3)*x));
[A,resnorm]=lsqcurvefit(f,a,x,y)
f(A,20)



A =

   1.0e+03 *

    5.7882    0.0025    0.0001


resnorm =

   3.3995e+04


ans =

   4.7438e+03

4.4

在这里插入图片描述

x=1:26;
y=[1807,2001,2158,2305,2422,2601,2753,2914,3106,3303,3460,3638,3799,3971,4125,4280,4409,4560,4698,4805,4884,4948,5013,5086,5124,5163];
a=[6000,2,0.1,0.1];
f=@(a,x)a(1)./(1+a(2)*exp(-a(3)*x-a(4)*x.^2));
[A,resnorm]=lsqcurvefit(f,a,x,y)
 t=27;
while    f(A,t+1)-f(A,t)>=10
 t=t+1;
end
f(A,t)




A =

   1.0e+03 *

    5.3860    0.0021    0.0001    0.0000


resnorm =

   9.1025e+03


ans =

   5.3409e+03

4.5

在这里插入图片描述

x=1:26;
y=[1807,2001,2158,2305,2422,2601,2753,2914,3106,3303,3460,3638,3799,3971,4125,4280,4409,4560,4698,4805,4884,4948,5013,5086,5124,5163];
a=[2,0.1,0.1];%r、k、a
f=@(a,x)a(1)*exp(a(2)*x+a(3));
[A,resnorm]=lsqcurvefit(f,a,x,y)
 t=27;
while    f(A,t+1)-f(A,t)>=10
 t=t+1;
end
f(A,t)

A =

    2.4511    0.3152   -0.0841


resnorm =

   2.3628e+08


ans =

   Inf
实验是一种重要的学习方法,可以帮助学生更好地掌握理论知识并提升实际操作能力。南京邮电大学通达学院开设了许多实验课程,其中一门是"熟悉MATLAB环境"。 MATLAB是一种常用的数学软件和编程语言,被广泛应用于不同领域的科学与工程计算中。学习MATLAB的环境,对于提高学生的编程和计算能力具有重要意义。 在实验一中,学生可以通过模拟实践的方式,深入学习MATLAB环境的使用。首先,学生会了解到MATLAB的特点和功能,掌握其基本操作和命令的使用。例如,学生可以学习如何创建变量、矩阵和数组,并进行数值计算、绘图和数据处理。 其次,在实验过程中,学生还可以学习如何使用MATLAB解决各种数学问题。通过编写MATLAB程序和脚本,学生可以实现自己的数学计算和模拟,并通过调试和修改程序,提高代码的效率和准确性。 此外,在实验一中,学生还可以学习MATLAB的应用案例,包括信号处理、图像处理、控制系统等方面的应用。通过这些案例,学生可以更好地理解MATLAB在科学与工程计算中的实际应用,并将其运用到自己的学习和研究中。 总之,实验一"熟悉MATLAB环境"是南京邮电大学通达学院为学生提供的一门重要实验课程,通过实践学习MATLAB的使用,能够帮助学生提高编程和计算能力,并为将来的学术研究和工程实践打下坚实基础。
评论 71
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值