常州大学数值分析作业

1.解：

（1）x = [ 3*π/8 π/2];

Y = cos(x);

x0 = π/3;

[A,Y] = lagrange(x,y,x0);

P1 = vpa(poly2sym(A),3)

结果如下：

P1 = 1.53*x - 0.974

Y = 0.5102

（2）x = [π/4 3*π/8π/2];

Y = cos(x);

[A,Y] = lagrange(x,y,x0);

P2=vpa(poly2sym(A),3)

结果如下：

P2 = 1.18*x^2 - 0.455*x - 0.189

Y = 0.4973

（3）x = [0 π/8 π/4 3*π/8 π/2];

Y = cos(x);

[A,Y]=lagrange(x,y,x0);

结果如下：

P3 = x^4 + 0.00282*x^3 - 0.514*x^2 + 0.0232*x + 0.0287

Y = 0.5001

7.function [T]=aitken(x,y,x0,T0)

If nargin == 3

T0=[];

end

n0=size(T0,1);

m=max(size(x));

n=n0+m;T=zeros(n,n+1);

T(1:n0,1:n0+1)=T0;T(n0+1:n,1)=x;T(n0+1:n, 2)=y;

if

n0==0

i0=2;

else

i0=n0+1;

End

x=[0 1];

y=[0.5 1.25];

x0=2.8;

T0=aitken(x,y,x0);

T=T0;

x=[3.0,4.0]';

y=[3.5,2.75]';x0=2.8;

T=aitken(x,y,x0,T0);

n=max(size(x))+size(T0,1);

for i=1:n

for j=1:i+1

fprintf('%10.4f',T(i,j));

end

fprintf('\n');

End

Return

0.0000 0.5000 0 0 0

1.0000 1.2500

2.6000 0 0

3.0000 3.5000 3.3000 3.2300 0

4.0000 2.7500 2.0750 2.2850 3.4190

16.function [C,D,Y]=newpoly(x0,y0,x)

if nargin < 2 | nargin> 3

error( 'Incorrect Number of Inputs'); end

if length(x0)~=length(y0)

error('The length of x0 must be equal to it of y0');

end n=length(x0); D=zeros(n,n); D(:,1)=y0'; for j=2:n%计算差商表

for k=j:n

If abs(x0(k)-x0(k-j+1))

error('DividedbyZero,therearetwonodesarethes ame');

end

D(k,j)=(D(k,j-1)-D(k-1,j-1))/(x0(k)-x0(k-j+1));

End

End

C=D(n,n);

For k=(n-1):-1:1

C=conv(C,poly(x0(k)));

m=length(C);

C(m)=C(m)+D(k,k);

end

If nargin==3

Y=polyval(C,x);

endC=fliplr(C);

Return

x = [0 1 2 3 4 ];

y = [0.5,1.25,2.75,3.5,2.75];

[A,Y]=lagrange(x,y,x0)

x0 = [0 1 2 3 4 ];

y0 = [0.5,1.25,2.75,3.5,2.75];

[C,D,X]=newpoly(x0,y0,x)

plot(x,Y,'b-',x0,X,'r:')

A = 0.5000 -0.3125 1.4687 -0.4375 0.0313

Y = 0.5000 1.2500 2.7500

3.5000 2.7500

C = 0.0313 -0.4375 1.4688 -0.3125 0.5000

D = 0.5000 0 0

0 0

1.2500 0.7500 0 0 0

2.7500 1.5000 0.3750 0 0

3.5000 0.7500 -0.3750 -0.2500 0

2.7500 -0.7500 -0.7500 -0.1250 0.0313

X = 0.5000 1.2500 2.7500 3.5000 2.7500

fl(x)=0.5*x^4 - 0.312*x^3 + 1.47*x^2 - 0.438*x + 0.0312

fn(x)=0.5*x^4 - 0.312*x^3 + 1.47*x^2 - 0.438*x + 0.0312

6. 解：

对题中函数进行变形：

原式→y/x = a* exp(b*x) →ln(y/x) = ln(a) + b*exp(x) 化为线性形式

计算：

>> a = [1 2 3 4 5];

>> b = [1.222 2.984 5.466 8.902 13.592];

>> x = exp(a);

>> y = log(b)-log(a);

>> n = 1; >> [C]=lspoly(x,y,n);

>> y = vpa(poly2sym(C),3)

结果如下：

y = 0.00464*x + 0.384

写成题中拟合函数的形式即为：y = 1.4679*x*exp(0.00464*x)

7.function [a0,a1,a2]=ec2(h,w)

S=log(s)';N=length(h);

A=zeros(N,3);

for i=1:5

A(i,1)=1;

A(i,2)=log(h(i));

A(i,3)=log(w(i));

end

c=inv(A'*A)*(A'*S); a0=exp(c(1)); a1=c(2); a2=c(3);

return

%给出数据

h=[175 172 183 164 156]; w=[80 90 80 70 65];

s=[1000 900 1200 750 800]; [a0, a1,a2]=ec2(h,w,s)

输出结果为：

a0 =1.614815742043648e-04

a1 =3.383163094165866

a2 =-0.419165011582663

8.x= lsqcurvefit(fun,x0,xdata,ydata)

x=lsqcurvefit(fun,x0,xdata,ydata,lb,

ub)

x=lsqcurvefit(fun,x0,xdata,ydata,lb, ub,options)

[x,resnorm] = lsqcurvefit(…)

[x,resnorm,residual]=lsqcurvefit(…)

[x,resnorm,residual,exitflag]= lsqcu rvefit(…)

[x,resnorm,residual,exitflag,output] = lsqcurvefit(…)

[x,resnorm,residual,exitflag,output,l ambda] = lsqcurvefit(…)

[x,resnorm,residual,exitflag,output,l ambda,jacobian] =lsqcurvefit(…)

function F = myfun(x,xdata)

F=(x(1).*xdata).*(exp(x(2).*xdat a));

xdata=[1,2,3,4,5];

ydata=[1.222,2.984,5.466,8.902,13. 592]; x0=[0,0];

[x,resnorm]=lsqcurvefit(@myfun,x0, xdata,ydata)

输出结果为：

Local minimum found.

Optimization completed because t he size of the gradient is less than t he default value of the function toler ance.

x =

0.999958348976391 0.2000141328 12834

aaresnorm = 8.067930437509675e -7