并针对下面一组具体实验数据 0.25 0.5000 0.3 0.5477 0.39 0.6245 0.45 0.6708 .
0.53 0.7280 求解,其中边界条件为
解:Matlab计算程序为:
clear clc
x=[0.25 0.3 0.39 0.45 0.53]
y=[0.5000 0.5477 0.6245 0.6708 0.7280] n=length(x); for i=1:n-1
h(i)=x(i+1)-x(i); end
for i=1:n-2
k(i)=h(i+1)/(h(i)+h(i+1)); u(i)=h(i)/(h(i)+h(i+1)); end
for i=1:n-2
gl(i)=3*(u(i)*(y(i+2)-y(i+1))/h(i+1)+k(i)*(y(i+1)-y(i))/h(i)); end
g0=3*(y(2)-y(1))/h(1);
g00=3*(y(n)-y(n-1))/h(n-1); g=[g0 gl g00]; g=transpose(g) k1=[k 1]; u1=[1 u];
Q=2*eye(5)+diag(u1,1)+diag(k1,-1) m=transpose(Q\\g) symsX; for i=1:n-1
p1(i)=(1+2*(X-x(i))/h(i))*((X-x(i+1))/h(i))^2*y(i); p2(i)=(1-2*(X-x(i+1))/h(i))*((X-x(i))/h(i))^2*y(i+1); p3(i)=(X-x(i))*((X-x(i+1))/h(i))^2*m(i); p4(i)=(X-x(i+1))*((X-x(i))/h(i))^2*m(i+1); p(i)=p1(i)+p2(i)+p3(i)+p4(i); p(i)=simple(p(i)); end s1=p(1) s2=p(2) s3=p(3) s4=p(4) for k=1:4
for z=x(k):0.001:x(k+1)
q=eval(subs(p(k),'X','z')); plot(z,q,'b') holdon end end gridon legend('èy′??ùì??ú??') title('2??μ') xlabel('x') ylabel('p')
计算结果为:
s1 =
-(705394867539368680*X^3 - 529046150654530286*X^2 + 23087199381998953*X - 40023205577025431)/112589990684262400 s2 =
(257361898089296225*X^3)/136796838681378816 - (160081743506404901*X^2)/60798594969501696 + (404209705972252727*X)/202661983231672320 + 429142243010323951/3166593487994880000 s3 =
- (3495912536773825*X^3)/7599824371187712 + (1437374409830143*X^2)/13510798882111488 + (10427488839800859*X)/11258999068426240 + 12358231431982259943/45035996273704960000 s4 =
(38626753769033575*X^3)/18014398509481984 - (245666153971053021*X^2)/72057594037927936 + (3614707928905781673*X)/1441151880758558720 + 26732501105874704913/720575940379279360000
通过绘制曲线,发现为3次样条曲线,且数据拟合较好。