南京邮电大学数学实验最新版

最新推荐文章于 2024-03-18 16:58:46 发布

Toas7

最新推荐文章于 2024-03-18 16:58:46 发布

阅读量7.6w

收藏 1k

点赞数 402

文章标签：算法学习

本文链接：https://blog.csdn.net/Toas7/article/details/128008737

版权

1.1(1)

syms x;
m=225;
y=(log(1+x-m*x^2)-x)/(1-cos(x));
limit(y,x,0)

ans=-1

1.1(2)

syms x;
m=225;
y=(sqrt(2+m*x^2)-atan(m*x))/x;
limit(y,x,inf)

ans=15

1.2

syms x;
m=225;
y=exp(m*x)*sin(x);
s1=diff(y,x,2)
s2=diff(y,x,6);
subs(s2,x,0)

y =
 
exp(225*x)*sin(x)
 
 
s1 =
 
450*exp(225*x)*cos(x) + 50624*exp(225*x)*sin(x)
 
 
ans =
 
3459674532600

1.3

syms x;
m=225;
y=(x+sin(x))/(1+cos(x));
int(y,x)

f=log(1+m*x)-m*x;
int(f,0,1)

 
ans =
 
x*tan(x/2)
 
 
ans =
 
(226*log(226))/225 - 227/2

1.4

syms x;
m=225;
y=cos(x)*(m/200+sin(x));
taylor(y,x,0,'order',5)

ans =
 
(3*x^4)/64 - (2*x^3)/3 - (9*x^2)/16 + x + 9/8

1.5

x(1)=rand(1);
m=225;
for n=2:10
x(n)=sqrt(m/100+x(n-1));
end;
x

x =

    0.1270    1.5417    1.9472    2.0487    2.0733    2.0793    2.0807    2.0810    2.0811    2.0811

1.6

m=225;
A=[4,2,m-2;-3,0,5;1,5,2*m];
B=[3,4,0;2,0,-3;-2,1,1];
det(A)
inv(A)
eig(A)
[P,D]=eig(A)
A/B
A\B
C=[A B]
[RC,P]=rref(C)

ans =

 -735.0000


ans =

    0.0340   -0.2925   -0.0136
   -1.8435   -2.1456    0.9374
    0.0204    0.0245   -0.0082


ans =

    3.8820
   -0.4202
  450.5383


P =

    0.7944   -0.1202    0.4468
   -0.6074   -0.9927    0.0070
    0.0050    0.0113    0.8946


D =

    3.8820         0         0
         0   -0.4202         0
         0         0  450.5383


ans =

   18.6400  -98.5200  -72.5600
    0.0400   -1.7200   -0.1600
   36.9200 -197.5600 -142.6800


ans =

   -0.4558    0.1224    0.8639
  -11.6966   -6.4367    7.3741
    0.1265    0.0735   -0.0816


C =

     4     2   223     3     4     0
    -3     0     5     2     0    -3
     1     5   450    -2     1     1


RC =

    1.0000         0         0   -0.4558    0.1224    0.8639
         0    1.0000         0  -11.6966   -6.4367    7.3741
         0         0    1.0000    0.1265    0.0735   -0.0816


P =

     1     2     3

1.7(1)

f.m

function y=f(x)
if x>=0&&x<=1/2
    y=2*x;
elseif x>1/2&&x<=1
    y=2*(1-x);
end

g.m

function y=g(x,f)
n=length(x);
for i=1:n
    y(i)=f(x(i));
end
end

main

fplot(@(x)g(x,@f),[0,1])

1.7(2)

f1.m

function y1=f1(x)
if -pi<=x&&x<0
    y1=x-1;
elseif 0<=x&&x<=pi
    y1=x+1;
end

g1.m

function y1=g1(x,f1)
n=length(x);
for i=1:n
    y1(i)=f1(x(i));
end
end

main

syms x;
fplot(@(x)g1(x,@f1),[-pi,pi])

1.8(1)

syms x;
m=225;
t=-m/25:0.1:m/25;
x=m/20*cos(t);
y=m/20*sin(t);
z=t;
plot3(x,y,z);
grid on

1.8(2)

syms x;
m=225;
t=-m/25:0.1:m/25;
x=cos(t)+t.*sin(t);
y=sin(t)-t.*cos(t);
z=-t;
plot3(x,y,z);
grid on;

1.9

syms x;
m=225;
a=[1000/m,500/m,100/m];
col=['r','b','g'];
x=linspace(-10,10);
for i=1:3
y=(1/(sqrt(2*pi).*a(i)))*exp(-x.^2/(2.*a(i).^2));
plot(x,y,'color',col(i))
hold on;
end
legend('1000','500','100')

1.10

syms x y;
ezplot('log(x^2+225*y)-x^3*y-sin(x)',[0,5,-5,5])

1.11

x=-5:0.1:5;y=x;
m=225;
[X Y]=meshgrid(x,y);
Z=m.*X.*exp(-X.^2-Y.^2);
mesh(X,Y,Z)

1.12(1)

syms x;
m=225;
ezplot('x^3+sqrt(225)*x^2+(225/3-3)*x-sqrt(225)*(1-225/27)',[-sqrt(m)/3-2,-sqrt(m)/3+2])

1.12(2)

fun.m

function y=fun(x);
y=x^3+sqrt(225)*x^2+(225/3-3)*x-sqrt(225)*(1-225/27);

Main.m

fzero(@fun,-6)
fzero(@fun,-5)
fzero(@fun,-4)

2.1(1)

dd.m

function y=dd(f2,x,n)
p=[x];
for i=2:n
    p(i)=f2(p(i-1))
end
end

main

syms x;
m=225;
dd(@(x)(2*x+1)/(x-m),-10,20)
dd(@(x)(2*x+1)/(x-m),0,20)
dd(@(x)(2*x+1)/(x-m),10,20)
dd(@(x)(2*x+1)/(x-m),20,20)

2.2

syms x;
m=225;
f=inline('1-2*abs(x-1/2)');
x0=1/m;
for i=1:1:50
	plot(i,f(x0),'.');
	x0=f(x0);
	hold on
end
hold off

2.3

Martin.m

function Martin(a,b,c,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))];
end
plot(m(1,:),m(2,:),'kx');
axis equal

m=225;
Martin(m,m,m,5000)
Martin(m,m,m,10000)
Martin(m,m,m,15000)
Martin(m,m,m,20000)

2.4

f=inline('(100*x+225)/(x^2+100)');
x0=5;
for i=1:1:20
     x0=f(x0);
     fprintf('%g,%g\n',i,x0)
end

2.5(1)

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.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.1

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

结果

xn =
 
(661*2^(1/2)*4^n)/4 - (661*0^n*5^(1/2))/10 - (657*0^n*2^(1/2))/4 + (657*5^(1/2)*6^n)/10
  (661*0^n*5^(1/2))/5 - (657*0^n*2^(1/2))/4 + (661*2^(1/2)*4^n)/4 - (657*5^(1/2)*6^n)/5

3.2

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

结果

xn =
 
(661*(2/5)^n*2^(1/2))/4 - (661*0^n*5^(1/2))/10 - (657*0^n*2^(1/2))/4 + (657*(3/5)^n*5^(1/2))/10
  (661*0^n*5^(1/2))/5 - (657*0^n*2^(1/2))/4 + (661*(2/5)^n*2^(1/2))/4 - (657*(3/5)^n*5^(1/2))/5

3.3(1)

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')

3.3(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;
结论：在迭代14 次后，发现数列存在极限为0.4

3.4

m=225;
A1=[m-1,m;1-m,-m]
p=[0.5;0.5];
[P,D]=eig(A1)
for i=1:20
    p(:,i+1)=A1*p(:,i);
end
fprintf('%2f,%2f\n',p)

4.1

for b=1:998
a=sqrt((b+2)^2-b^2);
if(a==floor(a))
fprintf('a=%i,b=%i,c=%i\n',a,b,b+2)
end
end

4.3

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

4.4

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)

4.5（1）

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,226.5];
t1=1790;x1=3.9;
t2=1900;x2=76.0;
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)

4.5(2)

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,226.5];
d0=length(t);d1=sum(t);d2=sum(t.^2);
A=[d0,d1;d1,d2];b1=sum(log(x));
b2=sum(log(x).*t);b=[b1;b2];
u=A\b;x0=exp(u(1))
k=u(2)
error=sum((x0*exp(k*t)-x).^2)

4.6

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];
f=@(a,x)a(1)./(1+a(2)*exp(-a(3)*x));
[A,resnorm]=lsqcurvefit(f,a,x,y)
t=26;
while    abs(A(1)/(1+A(2)*exp(-A(3)*t))-A(1)/(1+A(2)*exp(-A(3)*(t+1))))>1
    t=t+1;
end
t
A(1)/(1+A(2)*exp(-A(3)*(t+1)))
t=1:50;
plot(x,y,'*')
hold on;
 plot(t,A(1)./(1+A(2)*exp(-A(3)*t)));

4.7

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=26;
while    abs(A(1)/(1+A(2)*exp(-A(3)*t-A(4)*t.^2))-A(1)/(1+A(2)*exp(-A(3)*(t+1)-A(4)*(t+1).^2)))>1
    t=t+1;
end
t
A(1)/(1+A(2)*exp(-A(3)*(t+1)-A(4)*(t+1).^2))
t=1:50;
plot(x,y,'*')
hold on;
 plot(t,A(1)./(1+A(2)*exp(-A(3)*t-A(4)*(t+1).^2)));

4.8

x=1:27;
y=[21,65,127,281,558,923,1321,1801,2420,3125,3886,4638,5306,6028,6916,7646,8353,9049,9503,10098,10540,10910,11287,11598,11865,12086,12251];
a=[13000,2,0.1];
f=@(a,x)a(1)./(1+a(2)*exp(-a(3)*x));
[A,resnorm]=lsqcurvefit(f,a,x,y)
t=27;
while    abs(A(1)/(1+A(2)*exp(-A(3)*t))-A(1)/(1+A(2)*exp(-A(3)*(t+1))))>1
    t=t+1;
end
t
A(1)/(1+A(2)*exp(-A(3)*(t+1)))
t=1:50;
plot(x,y,'*')
hold on;
plot(t,A(1)./(1+A(2)*exp(-A(3)*t)));

4.9

x=1:27;
y=[21,65,127,281,558,923,1321,1801,2420,3125,3886,4638,5306,6028,6916,7646,8353,9049,9503,10098,10540,10910,11287,11598,11865,12086,12251];
a=[4,0.1,0.1];
f=@(a,x)a(1)*exp(exp(a(2)*x+a(3)));
[A,resnorm]=lsqcurvefit(f,a,x,y)
t=27;
while  abs(f(A,t)-f(A,t+1))>1
	t=t+1;
end
t
f(A,t+1)
t=1:50;
plot(x,y,'*')
hold on;
plot(t,f(A,t));

实验一

//散点图

plot(x,y,'.');
p=polyfit(x,y,4);
polyval(p,65)
plot(x,y,'*',x,polyval(p,x));

//模型建立

d0=length(x);d1=sum(x);d2=sum(x.^2);
A=[d0,d1;d1,d2];
b1=sum(log(y));b2=sum(log(y).^x);
B=[b1;b2];
u=A\B,y0=exp(u(1));
k=u(2);
error=sum((y0*exp(k*x)-y).^2)

x1=1;y1=178;
x2=50;y2=50;
A=[1,x1;1,x2];
B=[log(y1);log(y2)];
u=A\B;
y0=exp(u(1))
k=u(2)
error=sum((y0*exp(k*x)-y).^2)

//数据预测

a=[182.6729,-0.0259];
f=@(a,x)a(1)*exp(a(2)*x);
[A,resnorm]=lsqcurvefit(f,a,x,y)
t=61;
while  f(A,t)>1
	t=t+1;
end
t
t=1:100;
plot(x,y,'*')
hold on;
plot(t,f(A,t));
while j<=t
j=j+1;
S=S+f(A,j);
end
S