﻿ 17074191 - 计算数学达人 - 专，学者，数值代数，微分方程数值解

### 17074191

1. Decide whether the equations form a cubic spline.

(b) S(x)=2x3 + x2 + 4x + 5 on [0,1]   (x − 1)3 + 7(x − 1)2 + 12(x − 1) + 12 on [1,2]

Solution:

Property 1

There are n=3 data points. (0,5) (1,12) (2,32)

Property 2

S1(x)=6X^2+2X+4 ; S2(x)=3(x-1)^2+14(x-1)+12

S1(1)=6+2+4=12 ; S2(1)= 12

So, S1(1)= S2(1)

Property 3

S1’’(x)=12X+2 ; S2’’(x)=6(x-1) +14

S1’’(1)=14 ; S2’’(1)= 14

So, S1’’(1)= S2’’(1)

Therefore,S(x) is a cubic spline.

Computer Problems

1. Find the equations and plot the natural cubic spline that interpolates the data points (a) (0,3), (1,5),(2,4),(3,1) (b) (−1,3),(0,5),(3,1),(4,1),(5,1).

Solution:

Code:

%Program 3.5 Calculation of spline coefficients

%Calculates coefficients of cubic spline

%Input: x,y vectors of data points

% plus two optional extra data v1, vn

%Output: matrix of coefficients b1,c1,d1;b2,c2,d2;...

function coeff=splinecoeff(x,y)

n=length(x);v1=0;vn=0;

A=zeros(n,n); % matrix A is nxn

r=zeros(n,1);

for i=1:n-1 % define the deltas

dx(i)= x(i+1)-x(i); dy(i)=y(i+1)-y(i);

end

for i=2:n-1 % load the A matrix

A(i,i-1:i+1)=[dx(i-1) 2*(dx(i-1)+dx(i)) dx(i)];

r(i)=3*(dy(i)/dx(i)-dy(i-1)/dx(i-1)); % right-hand side

end

% Set endpoint conditions

% Use only one of following 5 pairs:

A(1,1) = 1; % natural spline conditions

A(n,n) = 1;

%A(1,1)=2;r(1)=v1;

%A(n,n)=2;r(n)=vn;

%A(1,1:2)=[2*dx(1) dx(1)];r(1)=3*(dy(1)/dx(1)-v1);

%clamped

%A(n,n-1:n)=[dx(n-1) 2*dx(n-1)];r(n)=3*(vn-dy(n-1)/dx(n-1));

%A(1,1:2)=[1 -1];

% parabol-term conditions, for n>=3

%A(n,n-1:n)=[1 -1];

%A(1,1:3)=[dx(2) -(dx(1)+dx(2)) dx(1)];

% not-a-knot, for n>=4

%A(n,n-2:n)=[dx(n-1) -(dx(n-2)+dx(n-1)) dx(n-2)];

coeff=zeros(n,3);

coeff(:,2)=Ar; % solve for c coefficients

for i=1:n-1 % solve for b and d

coeff(i,3)=(coeff(i+1,2)-coeff(i,2))/(3*dx(i));

coeff(i,1)=dy(i)/dx(i)-dx(i)*(2*coeff(i,2)+coeff(i+1,2))/3;

end

coeff=coeff(1:n-1,1:3);

(a) splinecoeff([0;1;2;3],[3;5;4;1])

ans =

2.6667         0   -0.6667

0.6667   -2.0000    0.3333

-2.3333   -1.0000    0.3333

S(x)=

(b) splinecoeff([-1;0;3;4;5],[3;5;1;1;1])

ans =

2.5629         0   -0.5629

0.8742   -1.6887    0.3176

-0.6824    1.1698   -0.4874

0.1950   -0.2925    0.0975

S(x)=

Chap 4.2

Exercises

1. Fit data to the periodic model y =F3(t)=c1 + c2cos2πt + c3sin2πt. Find the 2-norm error and the RMSE.

Solution:

(a)

1=c1+c2+0*c3

3=c1+0*c2+c3

2=c1-c2+0*c3

0=c1+0*c2-c3

A=[1,1,0;1,0,1;1,-1,0;1,0,-1]

b=[1;3;2;0]

c =

1.5000

-0.5000

1.5000

r =

0

0

0

0

n=norm(r,2)

n =

0

RMSE=0

(b)

1=c1+c2+0*c3

3=c1+0*c2+c3

2=c1-c2+0*c3

1=c1+0*c2-c3

A=[1,1,0;1,0,1;1,-1,0;1,0,-1]

b=[1;3;2;1]

c =

1.7500

-0.5000

1.0000

r =

-0.2500

0.2500

-0.2500

0.2500

n=norm(r,2)

n =

0.5000

RMSE=0.2500

(c)

3=c1+c2+0*c3

1=c1-c2+0*c3

3=c1+c2+0*c3

2=c1-c2+0*c3

A=[1,1,0;1,-1,0;1,1,0;1,-1,0]

b=[3;1;3;2]

A为奇异矩阵，无解

2. Fit data to the exponential model by using linearization. Find the 2-norm of the difference between the data points yi and the best model c1ec2ti.

Solution:

lnyi=lnc1+c2*ti

(a)

0=lnc1-2c2

ln2=lnc1

ln2=lnc1+c2

ln5=lnc1+2c2

A=[1,-2;1,0;1,1;1,2]

b=[0;ln2;ln2;ln5]

c =

0.6586

0.3615

exp(0.658559070187365)

ans =

1.9320

So,yi=1.9320*e0.3615*ti

(b)

0=lnc1

0=lnc1+c2

ln2=lnc1+c2

ln4=lnc1+2c2

A=[1,0;1,1;1,1;1,2]

b=[0;0;ln2;ln4]

c =

-0.1733

0.6931

exp(-0.173286795139986)

ans =

0.8409

So,yi=0.8409*e0.6931ti

Computer Problems

1. Fit the monthly data for Japan 2003 oil consumption, shown in the following table, with the periodic model (4.9), and calculate the RMSE:

Solution:

yi=c1ec2ti.

A=[1,6.224;1,6.665;1,6.241;1,5.302;1,5.073;1,5.127;1,4.994;1,5.012;1,5.108;1,5.377;1, 5.510 ;1,6.372]

b=[ln1; ln2; ln3; ln4; ln5; ln6;ln7;ln8;ln9;ln10; ln11;ln12]

c=(A'*A)(A'*b)

c =

5.2284

-0.6381

exp(5.228409900541486)

ans =

186.4960

So,yi=186.4960*e-0.6381ti

r=b-A*c

r =

-1.2571

-0.2825

-0.1476

-0.4591

-0.3821

-0.1653

-0.0960

0.0490

0.2281

0.5051

0.6852

1.3223

>> norm(r,2)

ans =

2.1454

RMSE=(norm(r,2)^2/12)^(1/2)

RMSE =

0.6193

3. Consider the world population data of Computer Problem 3.1.1. Find the best exponential ﬁt of the data points by using linearization. Estimate the 1980 population, and ﬁnd the estimation error.

Solution:

yi=c1ec2ti.

A=[1,3039585530;1,3707475887;1,5281653820;1,6079603571]

b=[ln1960; ln1970; ln1990; ln2000]

` `