当前位置：文档之家› 三次样条插值课后题集

三次样条插值课后题集

例1 设)(x f 为定义在[0,3]上的函数，有下列函数值表：且2.0)('0=x f ，1)('3-=x f ，试求区间[0,3]上满足上述条件的三次样条插值函数)(x s本算法求解出的三次样条插值函数将写成三弯矩方程的形式：)()6()()6()(6)(6)(211123131j j jj j j j j j j j jj j jj x x h h M y x x h h M y x x h M x x h M x s --+--+-+-=+++++其中，方程中的系数jj h M 6，jj h M 61+，jj j j h h M y )6(2-，jjj j h h M y )6(211++-将由Matlab代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。

以下为Matlab 代码：%============================= % 本段代码解决作业题的例1%============================= clear all clc% 自变量x 与因变量y ，两个边界条件的取值 IndVar = [0, 1, 2, 3]; DepVar = [0, 0.5, 2, 1.5];LeftBoun = 0.2;RightBoun = -1;% 区间长度向量，其各元素为自变量各段的长度h = zeros(1, length(IndVar) - 1);for i = 1 : length(IndVar) - 1h(i) = IndVar(i + 1) - IndVar(i);end% 为向量μ赋值mu = zeros(1, length(h));for i = 1 : length(mu) - 1mu(i) = h(i) / (h(i) + h(i + 1));endmu(i + 1) = 1;% 为向量λ赋值lambda = zeros(1, length(h));lambda(1) = 1;for i = 2 : length(lambda)lambda(i) = h(i) / (h(i - 1) + h(i));% 为向量d赋值d = zeros(1, length(h) + 1);d(1) = 6 * ( (DepVar(2) - DepVar(1) ) / ( IndVar(2) - IndVar(1) ) - LeftBoun) / h(1); for i = 2 : length(h)a = ( DepVar(i) - DepVar(i - 1) ) / ( IndVar(i) - IndVar(i - 1) );b = ( DepVar(i + 1) - DepVar(i) ) / ( IndVar(i + 1) - IndVar(i) );c = (b - a) / ( IndVar(i + 1) - IndVar(i - 1) );d(i) = 6 * c;endd(i + 1) = 6 *( RightBoun - ( DepVar(i + 1) - DepVar(i) ) / ( IndVar(i + 1) - IndVar(i) ) ) / h(i);% 为矩阵A赋值% 将主对角线上的元素全部置为2A = zeros( length(d), length(d) );for i = 1 : length(d)A(i, i) = 2;end% 将向量λ的各元素赋给主对角线右侧第一条对角线for i = 1 : length(d) - 1A(i, i + 1) = lambda(i);end% 将向量d的各元素赋给主对角线左侧第一条对角线for i = 1 : length(d) - 1A(i + 1, i) = mu(i);end% 求解向量MM =A \ d';% 求解每一段曲线的函数表达式for i = 1 : length(h)Coefs_1 = M(i) / (6 * h(i));Part_1 = conv( Coefs_1, ...conv( [-1, IndVar(i + 1)], ...conv( [-1, IndVar(i + 1)], [-1, IndVar(i + 1)] ) ) );S_1 = polyval (Part_1, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_2 = M(i + 1)/(6 * h(i));Part_2 = conv( Coefs_2, ...conv( [1, -IndVar(i)], ...conv( [1, -IndVar(i)], [1, -IndVar(i)] ) ) );S_2 = polyval (Part_2, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_3 = (DepVar(i) - M(i) * h(i)^2 / 6) / h(i);Part_3 = conv(Coefs_3, [-1, IndVar(i + 1)]);S_3 = polyval (Part_3, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_4 = (DepVar(i + 1) - M(i + 1) * h(i)^2 / 6) / h(i);Part_4 = conv(Coefs_4, [1, -IndVar(i)]);S_4 = polyval (Part_4, [IndVar(i) : 0.01 : IndVar(i + 1)]);S = S_1 + S_2 + S_3 + S_4;plot ([IndVar(i) : 0.01 : IndVar(i + 1)], S, 'LineWidth', 1.25)% 在样条插值曲线的相应位置标注该段曲线的函数表达式text(i - 1, polyval(Part_1, 3), ...['\itS', num2str(i), '(x)=', num2str(Coefs_1), '(', num2str( IndVar(i + 1) ),'-x)^{3}+', ...num2str(Coefs_2), '(x-', num2str( IndVar(i) ), ')^{3}+', num2str(Coefs_3), ...'(', num2str( IndVar(i + 1) ), '-x)+', num2str(Coefs_4), '(x-', num2str( IndVar(i) ), ')'], ...'FontName', 'Times New Roman', 'FontSize', 14) hold on end% 过x=1和x=2两个横轴点作垂线 % line([1, 1], [2.5, -0.5], 'LineStyle', '--'); line([2, 2], [2.5, -0.5], 'LineStyle', '--');% 为x 轴和y 轴添加标注xlabel( '\itx', 'FontName', 'Times New Roman', ... 'FontSize', 14, 'FontWeight', 'bold');ylabel( '\its(x)', 'FontName', 'Times New Roman', ...'Rotation', 0, 'FontSize', 14, 'FontWeight', 'bold');最终，三次样条插值函数s(x)表达式为：[][][]⎪⎩⎪⎨⎧∈-+-+-+--∈-+-+---∈+-++--=.3,2,)2(44.1)3(62.2)2(06.0)3(62.0,2,1,)1(62.2)2(08.0)1(62.0)2(42.0,1,0,08.0)1(06.042.0)1(06.0)(333333x x x x x x x x x x x x x x x x s曲线的图像如图所示：例2 已知函数值表：试求在区间[1,5]上满足上述函数表所给出的插值条件的三次自然样条插值函数)(x s本算法求解出的三次样条插值函数将写成三弯矩方程的形式：)()6()()6()(6)(6)(211123131j j jj j j j j j j j jj j jj x x h h M y x x h h M y x x h M x x h M x s --+--+-+-=+++++其中，方程中的系数jj h M 6，jj h M 61+，jj j j h h M y )6(2-，jjj j h h M y )6(211++-将由Matlab代码中的变量Coefs_1、Coefs_2、Coefs_3以及Coefs_4的值求出。

以下为Matlab 代码：%============================= % 本段代码解决作业题的例2%============================= clear all clc% 自变量x 与因变量y 的取值 IndVar = [1, 2, 4, 5]; DepVar = [1, 3, 4, 2];% 区间长度向量，其各元素为自变量各段的长度h = zeros(1, length(IndVar) - 1);for i = 1 : length(IndVar) - 1h(i) = IndVar(i + 1) - IndVar(i);end% 为向量μ赋值mu = zeros(1, length(h));for i = 1 : length(mu) - 1mu(i) = h(i) / (h(i) + h(i + 1));endmu(i + 1) = 0;% 为向量λ赋值lambda = zeros(1, length(h));lambda(1) = 0;for i = 2 : length(lambda)lambda(i) = h(i) / (h(i - 1) + h(i));end% 为向量d赋值d = zeros(1, length(h) + 1);d(1) = 0;for i = 2 : length(h)a = ( DepVar(i) - DepVar(i - 1) ) / ( IndVar(i) - IndVar(i - 1) );b = ( DepVar(i + 1) - DepVar(i) ) / ( IndVar(i + 1) - IndVar(i) );c = (b - a) / ( IndVar(i + 1) - IndVar(i - 1) );d(i) = 6 * c;endd(i + 1) = 0;% 为矩阵A赋值% 将主对角线上的元素全部置为2A = zeros( length(d), length(d) );for i = 1 : length(d)A(i, i) = 2;end% 将向量λ的各元素赋给主对角线右侧第一条对角线for i = 1 : length(d) - 1A(i, i + 1) = lambda(i);end% 将向量d的各元素赋给主对角线左侧第一条对角线for i = 1 : length(d) - 1A(i + 1, i) = mu(i);end% 求解向量MM =A \ d';% 求解每一段曲线的函数表达式for i = 1 : length(h)Coefs_1 = M(i) / (6 * h(i));Part_1 = conv( Coefs_1, ...conv( [-1, IndVar(i + 1)], ...conv( [-1, IndVar(i + 1)], [-1, IndVar(i + 1)] ) ) );S_1 = polyval (Part_1, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_2 = M(i + 1)/(6 * h(i));Part_2 = conv( Coefs_2, ...conv( [1, -IndVar(i)], ...conv( [1, -IndVar(i)], [1, -IndVar(i)] ) ) );S_2 = polyval (Part_2, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_3 = (DepVar(i) - M(i) * h(i)^2 / 6) / h(i);Part_3 = conv(Coefs_3, [-1, IndVar(i + 1)]);S_3 = polyval (Part_3, [IndVar(i) : 0.01 : IndVar(i + 1)]);Coefs_4 = (DepVar(i + 1) - M(i + 1) * h(i)^2 / 6) / h(i);Part_4 = conv(Coefs_4, [1, -IndVar(i)]);S_4 = polyval (Part_4, [IndVar(i) : 0.01 : IndVar(i + 1)]);S = S_1 + S_2 + S_3 + S_4;plot ([IndVar(i) : 0.01 : IndVar(i + 1)], S, 'LineWidth', 1.25)% 在样条插值曲线的相应位置标注该段曲线的函数表达式text(i, polyval(Part_1, 5), ...['\itS', num2str(i), '(x)=', num2str(Coefs_1), '(', num2str( IndVar(i + 1) ),'-x)^{3}+', ...num2str(Coefs_2), '(x-', num2str( IndVar(i) ), ')^{3}+', num2str(Coefs_3), ...'(', num2str( IndVar(i + 1) ), '-x)+', num2str(Coefs_4), '(x-', num2str( IndVar(i) ), ')'], ...'FontName', 'Times New Roman', 'FontSize', 14)hold onend% 过x=2和x=4两个横轴点作垂线 % line([2, 2], [4.5, 0.5], 'LineStyle', '--'); line([4, 4], [4.5, 0.5], 'LineStyle', '--');% 为x 轴和y 轴添加标注xlabel( '\itx', 'FontName', 'Times New Roman', ... 'FontSize', 14, 'FontWeight', 'bold');ylabel( '\its(x)', 'FontName', 'Times New Roman', ...'Rotation', 0, 'FontSize', 14, 'FontWeight', 'bold');最终，三次自然样条插值函数s(x)表达式为：[][][]⎪⎩⎪⎨⎧∈-+-+--∈-+-+----∈-+-+--=.5,4,)4(2)5(375.4)5(375.0,4,2,)2(75.2)4(75.1)2(1875.0)4(0625.0,2,1,)1(125.3)2()1(125.0)(3333x x x x x x x x x x x x x x s曲线的图像如图所示：[][][][]⎪⎪⎩⎪⎪⎨⎧∈-+-+--∈-+-+----∈-+-+----∈-+-+--=.53.0,45.0,)45.0(1.9)53.0(3987.8)53.0(1442.2,45.0,39.0,)39.0(1903.11)45.0(417.10)39.0(859.2)45.0(399.2,39.0,30.0,)3.0(9518.6)39.0(1137.6)3.0(5993.1)39.0(4806.3,30.0,25.0,)25.0(9697.10)3.0(10)25.0(2652.6)(333333x x x x x x x x x x x x x x x x x x x s试求在区间[0.25,0.53]上满足上述函数表所给出的插值条件的三次自然样条插值函数s(x)求解出的三次样条插值函数将写成三弯矩方程的形式：)()6()()6()(6)(6)(211123131j jjj j j jj j j j jj j jj x x h h M y x x h h M y x x h M x x h M x s --+--+-+-=+++++本题采用和例2基本相同的Matlab 代码，只改变初始条件。

相关主题

相关文档推荐：