MATLAB上机实验——最佳平方逼近
研究生数值分析(20)函数的最佳平方逼近
L1 )
3 2
1 1
t 2
1 tdt 6 15
可知
q1(t)
2 3
L0 (x)
6 15
L1 ( x)
2 3
6 15
t,
1 t 1
把 t =2x-1代人 q1(t) 得 x 在区间[0，1]上的一次最佳平方 逼近多项式
p1(t)
2 3
6 15
(2x
1)
4 15
12 15
m
m
[ * ( xi )
i 1
yi ]2
min ( x )
[ ( xi )
i 1
yi ]2
其中 (x) 为函数类Φ 中任意函数。
因此，用最小二乘法解决实际问题包含 如下2个基本环节：
（1）确定函数类Φ ，即确定 (x) 的形式。 这不是一个单纯的数学问题，还与其
12 15
所求的最佳平方逼近元素为
p(x) 4 12 x, 15 15
0 x 1
二、正交函数系在最佳平方逼近中的应用 对于一般的基底 0 (x),1(x),,n (x)
当 n 较大时，计算法方程中的 (k , j ) 以及求解法方程的计算量都是很大的。 1, x, x2 ,, xn 作基底，当ρ(x)≡1时， 虽然 (k , j ) (xk , x j ) 容易计算，但由此形成 的法方程系数矩阵G在 n稍大时是病态矩阵， 在计算机上求解法方程，其结果不太可靠。
§6 函数的最佳平方逼近 一、最佳平方逼近的概念与解法
用简单函数 p (x)去近似一个给定区间[a， b]上的连续函数 f (x)，是函数逼近要研究的 问题。度量逼近误差标准有许多种，这里 介绍一种称为平方逼近的函数逼近。
MATLAB语言与数学实验(江世宏编著)PPT模板
第二篇 数值分析
第7章 插值法
03 实验题
02
7.8 拉格朗日插值 多项式与三次样条
插值函数的比较
01
7.7 拉格朗日插值 多项式与埃米特插
值多项式的比较
第二篇 数值分析
第8章 函数的数值逼 近
0 1 8.1 伯恩斯坦多项式逼近连续函数 的动画演示
02
8.2 式
函数的最佳平方逼近多项
0 3 8.3 希尔伯特矩阵的病态性
4.1 二维数据曲线图
4.1.1 绘制单根二维曲线 4.1.2 绘制多根二维曲线 4.1.3 设置曲线样式 4.1.4 图形标注 4.1.5 坐标控制 4.1.6 图形窗口的分割
第4章 绘图
4.2 其他二维图形
4.2.1 其他坐标系下的二维数据曲线图 4.2.2 隐函数作图
第4章 绘图
12.5 选
05
方阵列主
元的lu分
解
04
12.4 方阵
03 1 2 . 3 高 斯全主元 消去法
的lu分解
第二篇 数值分析
第12章 线性方程组求 解
0 1
12.7
法
平方根
0 2
12.8 改进的
平方根法
0 3
第3章 编程
3.2 程序控制结构
3.2.1 顺序结构 3.2.2 选择结构 3.2.3 循环结构
第3章 编程
3.3 函数文件
3.3.1 函数文件的基本结构 3.3.2 函数调用 3.3.3 字符串函数
第4章 绘图
4.1 二维数据曲线图 4.2 其他二维图形 4.3 三维图形 习题
第4章 绘图
第2章 常用的数学函数
2.4 matlab符号计算
2.4 正交多项式和最佳平方逼近
dn Ln ( x ) e x n ( x n e x ) dx
第二章 插值与拟合
L2 ( x ) x 2 4 x 2, L3 ( x ) x 3 9 x 2 18 x 6, L4 ( x ) x 4 16 x 3 72x 2 96 x 24 L4 ( x ) x 5 25x 4 200x 3 600x 2 600x 120
（） (sinix, cos jx) 0, i , j 1,2,, n 3
（） （， dx 2 ； 4 11 ）
(1, sinix ) 0, (1, cosix ) 0, i 1,, n。
第二章 插值与拟合
正交多项式的三项递推公式:
n 设 {k ( x)}k 0 为 a，b]具有权函数 ( x) 的正交多项式组，i ( x) [
0，当i j , 且i , j 1 （）(cosix, cos jx ) cosix cos jxdx 1 ； ，当i j 0 0, 当i j , 且i , j 1 （） (sinix, sin jx ) 2 ； , 当i j 0
(2)
1 连续区间上正交多项式
0, i j (第二章 i , j ) 插值与拟合 ai 0, i j
连续区间上的正交多项式的概念与离散 点集上的正交多项式概念相似，只要将内积 的定义作相应的改变 。 定义2.10 函数f (x)和 g (x)在连续意义下的内积定义为
( f , g ) ( x) f ( x) g ( x)dx, f , g C[a, b]
特别地，
第二章 插值与拟合
最佳平方逼近第二次
由(3.8)和(3.9)确定的sn* ( x)满足
lim ||
n
f
(x) sn* (x) ||2
0.
特别地,勒让德多项式pn( x)在[1,1]上带权( x) 1正交, 且
11
pm ( x)
pn( x)dx
0, 2
2n
m n, , m n. 1
于是, 系数 从而
a*j
2
j 1( 2
f
,
...
...
...
bn n
an
cn
定理 最小二乘拟合多项式存在唯一 (n < m)。
证明：记法方程组为 Ba = c .
则有
B ΦTΦ c ΦT y
其中
Φ 1...
x1 ...
x12 ...
1 xm xm2
对任意 u 0 Rn1 ，必有 Φ u 0。
... ...
x1n ...
1 3
a0
a1
2 3
2 5
,
aa01
4 15
,
4 5
.
p1*
(x)
4 15
4 5
x.
平方误差: || (x) ||22
1
0
xdx
(
4 15
2 3
4 5
2 5
)
1 450
0.0022.
3.3.2、用正交函数族求最佳平方逼近
设f ( x) C[a,b], span{0,,n} C[a,b]，0,,n是
S * (x) a0 a1x an xn
直接通过解法方程得到 H n 中的最佳平方逼近多项式是一致 的.
只是当 n 较大时法方程出现病态，计算误差较大，不能 使用，而用勒让德展开不用解线性方程组，不存在病态问题， 因此通常都用这种方法求最佳平方逼近多项式.
平方逼近——精选推荐
第四章 平方逼近教学目的及要求：掌握最小二乘法、初步掌握平方逼近及直交多项式、平方度量意义下函数的逼近问题。
本书第二章是用数量)()(m ax x f x p f p bx a -=-≤≤来度量逼近多项式)(x p 与已知函数)(x f 的近似程度。
若,,0)()(∞→→-n x f x p n 则意味着序列{})(x p n 在区间],[b a 上一致收敛到)(x f 。
一致逼近度量，亦称Tchebyshev 度量是很重要的一类度量标准。
然而由于它的非线性特性，使得最佳一致逼近多项式的构造问题十分困难。
对于许多问题来讲，人们需要求出的是在确定意义下的“整体”近似。
本章讨论一类新的度量----平方度量意义下函数的逼近问题。
先从最小二乘法谈起。
§1. 最小二乘法最小二乘法起源于以测量和观测为基础的天文学。
Gauss 在1794年利用最小二乘法解决了多余观测问题，当时他只有十七岁。
可以用下面的简单例子描述这类问题。
假定通过观测或实验得到如下一组数据（即列表函数）：我们的目的是一简单的式子表出这些数据间的关系。
从分析数据看出，这些点差不多分布在一条直线上，因此我们自然想到用线性式b ax y +=表示它们之间的关系。
这就须定出参数a 和b 的值来。
这实际上是多余观测问题，用插值法不能确定出a 和b 的值。
代定参数的确定归结为矛盾方程组的求解问题。
假定有某方法可以定出a 和b ，则按bx a y +=，给出一个x 便可以算出一个y 。
我们记).8,,1( =+=k bxa y kk y 称为k y 的估计值，显然它们不会是完全相同的，它们之间的差（通常称为残差）)8,,1( =-=k y y kk k ε无疑是衡量被确定的参数a 和b （也就是近似多项式b ax y +=）好坏的重要标志。
可以规定许多原则来确定参数b a ,。
例如（1） 参数的确定，将使残差绝对值中最大的一个达到最小，即kkT εmax =为最小；（2） 参数的确定，将使残差绝对值之和达到最小，即∑kk ε为最小；（3） 参数的确定，将使残差的平方和达到最小，即∑2k ε为最小。
计算方法最佳平方逼近-最小二乘法
i0
i0
为最小。
这种要求误差（偏差）平方和最小的拟合称为 曲线拟合的最小二乘法。
（1）直线拟合
设已知数据点 (xi , yi ), i 1, 2, … , m 分布大致为 一条直线。作拟合直线 y(x) a0 a1x, 该直线
不是通过所有的数据点 (xi , yi ) ,而是使偏差平方和
F(a0 , a1 )
计算方法 (Numerical Analysis)
第5次 最佳平方逼近与曲线拟合的最小二乘法
主要内容
• 最佳平方逼近 • 曲线拟合的最小二乘法
最佳平方逼近
函数逼近的类型
• 最佳一致逼近：使用多项式对连续函数进行一致 逼近。逼近误差使用范数
||
f(x)
-
s(x)
||
max
a x b
|
f(x)
-
1(1 x2 )dx 1(0.934 0.426x) 1 x2dx
0
0
0.0026
最大值误差 :
同学们自己求一下
|| δ(x) || max | 1 x2 (0.934 0.426x) | 0.066
例题 求f(x) x在[1/4,1]上的在Φ span{1, x} 中的关于ρ(x) 1的最佳平方逼近多项式。
10 27
88 x 135
平方误差 :|| δ(x) ||22
1xdx
1
( 10 27
7 12
31 80
) 88
135
4
1.02
p1* (x)
10 27
88 x. 135
1
f(x) x
平方误差 : || δ(x) ||22 0.0001082.
计算方法 第五章第二节最佳平方逼近
n
2
i 0
a
i 0
上述方程组称为正规方程组。也可以写为
( p, j ) ( f , j )，j 0,1,..., n.
由于0 ( x), 1 ( x),..., n ( x) 线性无关，由性质5.2.3，该方程组 的系数矩阵非奇异，因而方程组存在惟一解。
可以证明，最佳平方问题的解存在惟一且就是正规方程组的解。
b
j i,
j i,
则称多项式族 {g n ( x)} 在[a, b] 上带权 ( x) 正交，并称 g n ( x)是[a, b] 上带权 ( x)的 n 次正交多项式。
一般情况下，当权函数 ( x)及区间[a, b] 给定后，人们 可通过Gram－Schmidt正交化过程，由{1, x,..., x n }构造 出相应的正交多项式。
2
的最小值。
由多元函数取极值的必要条件 S 0,
a j
j 0,1,..., n,
得
n aii ( x) f ( x) j ( x)dx 0， j 0,1,..., n. a ( x) i 0
b
于是有
)，j 0,1,..., n. ( , ) a ( f , S (a0 , ai1 ,...,jan ) :i ( x) j aii ( x) f ( x) dx
2
2
2
2
f g f g 2 f g
2
2
2
2
,
f , g Y.
二、 函数的最佳平方逼近
已知函数 f ( x) C[a,b] 及C[a,b]中的一个子集 span{0 , 1 ,..., n }，如果 p( x) span{0 , 1,..., n}，使得
Matlab数值实验定积分的近似计算教程
Matlab数值实验定积分的近似计算教程一、问题背景与实验目的二、相关函数（命令）及简介三、实验内容1．矩形法2．梯形法3．抛物线法4.直接应用Matlab命令计算结果四、自己动手一、问题背景与实验目的利用牛顿—莱布尼兹公式虽然可以精确地计算定积分的值，但它仅适用于被积函数的原函数能用初等函数表达出来的情形．如果这点办不到或者不容易办到，这就有必要考虑近似计算的方法．在定积分的很多应用问题中，被积函数甚至没有解析表达式，可能只是一条实验记录曲线，或者是一组离散的采样值，这时只能应用近似方法去计算相应的定积分．本实验将主要研究定积分的三种近似计算算法：矩形法、梯形法、抛物线法．对于定积分的近似数值计算，Matlab有专门函数可用．二、相关函数（命令）及简介1．sum(a)：求数组a的和．2．format long：长格式，即屏幕显示15位有效数字．（注：由于本实验要比较近似解法和精确求解间的误差，需要更高的精度）．3．double()：若输入的是字符则转化为相应的ASCII码；若输入的是整型数值则转化为相应的实型数值．4．quad()：抛物线法求数值积分．格式： quad(fun,a,b) ，注意此处的fun是函数，并且为数值形式的，所以使用*、/、^等运算时要在其前加上小数点，即 .*、./、.^等．例：Q = quad('1./(x.^3-2*x-5)',0,2);5．trapz()：梯形法求数值积分．格式：trapz(x,y)其中x为带有步长的积分区间；y为数值形式的运算（相当于上面介绍的函数fun）例：计算x=0:pi/100:pi;y=sin(x);trapz(x,y)6．dblquad()：抛物线法求二重数值积分．格式：dblquad(fun,xmin,xmax,ymin,ymax)，fun可以用inline 定义，也可以通过某个函数文件的句柄传递．例1：Q1 = dblquad(inline('y*sin(x)'), pi, 2*pi, 0, pi)顺便计算下面的Q2，通过计算，比较Q1 与Q2结果（或加上手工验算），找出积分变量x、y的上下限的函数代入方法．Q2 = dblquad(inline('y*sin(x)'), 0, pi, pi, 2*pi) 例2：Q3 = dblquad(@integrnd, pi, 2*pi, 0, pi)这时必须存在一个函数文件integrnd.m：function z = integrnd(x, y)z = y*sin(x);7．fprintf（文件地址，格式，写入的变量）：把数据写入指定文件．例：x = 0:.1:1;y = [x; exp(x)];fid = fopen('exp.txt','w'); %打开文件fprintf(fid,'%6.2f %12.8f\n',y); %写入fclose(fid) %关闭文件8．syms 变量1 变量2 …：定义变量为符号．9．sym('表达式')：将表达式定义为符号．解释：Matlab中的符号运算事实上是借用了Maple的软件包，所以当在Matlab中要对符号进行运算时，必须先把要用到的变量定义为符号．10．int(f,v,a,b)：求f关于v积分，积分区间由a到b．11．subs(f，'x'，a)：将 a 的值赋给符号表达式 f 中的 x，并计算出值．若简单地使用subs(f)，则将f的所有符号变量用可能的数值代入，并计算出值．三、实验内容1．矩形法根据定积分的定义，每一个积分和都可以看作是定积分的一个近似值，即在几何意义上，这是用一系列小矩形面积近似小曲边梯形的结果，所以把这个近似计算方法称为矩形法．不过，只有当积分区间被分割得很细时，矩形法才有一定的精确度．针对不同的取法，计算结果会有不同，我们以为例（取），（1）左点法：对等分区间，在区间上取左端点，即取，0.78789399673078，理论值，此时计算的相对误差（2）右点法：同（1）中划分区间，在区间上取右端点，即取，0.78289399673078，理论值，此时计算的相对误差（3）中点法：同（1）中划分区间，在区间上取中点，即取，0.78540024673078，理论值，此时计算的相对误差如果在分割的每个小区间上采用一次或二次多项式来近似代替被积函数，那么可以期望得到比矩形法效果好得多的近似计算公式．下面介绍的梯形法和抛物线法就是这一指导思想的产物．2．梯形法等分区间，相应函数值为（）．曲线上相应的点为（）将曲线的每一段弧用过点，的弦（线性函数）来代替，这使得每个上的曲边梯形成为真正的梯形，其面积为，．于是各个小梯形面积之和就是曲边梯形面积的近似值，，即，称此式为梯形公式．仍用的近似计算为例，取，0.78539399673078，理论值，此时计算的相对误差很显然，这个误差要比简单的矩形左点法和右点法的计算误差小得多．3．抛物线法由梯形法求近似值，当为凹曲线时，它就偏小；当为凸曲线时，它就偏大．若每段改用与它凸性相接近的抛物线来近似时，就可减少上述缺点，这就是抛物线法．将积分区间作等分，分点依次为，，对应函数值为（），曲线上相应点为（）．现把区间上的曲线段用通过三点，，的抛物线来近似代替，然后求函数从到的定积分：由于，代入上式整理后得同样也有……将这个积分相加即得原来所要计算的定积分的近似值：，即这就是抛物线法公式，也称为辛卜生（Simpson）公式．仍用的近似计算为例，取，=0.78539816339745，理论值，此时计算的相对误差4. 直接应用Matlab命令计算结果（1）数值计算方法1：int('1/(1+x^2)','x',0,1) (符号求积分)方法2：quad('1./(1+x.^2)',0,1) （抛物线法求数值积分）方法3：x=0:0.001:1;y=1./(1+x.^2);trapz(x,y) （梯形法求数值积分）（2）数值计算方法1：int(int('x+y^2','y',-1,1),'x',0,2) （符号求积分）方法2：dblquad(inline('x+y^2'),0,2,-1,1) （抛物线法二重数值积分）四、自己动手1．实现实验内容中的例子，即分别采用矩形法、梯形法、抛物线法计算，取，并比较三种方法的精确程度．2．分别用梯形法与抛物线法，计算，取．并尝试直接使用函数trapz()、quad()进行计算求解，比较结果的差异．3．试计算定积分．（注意：可以运用trapz()、quad()或附录程序求解吗？为什么？）4．将的近似计算结果与Matlab中各命令的计算结果相比较，试猜测Matlab中的数值积分命令最可能采用了哪一种近似计算方法？并找出其他例子支持你的观点．5．通过整个实验内容及练习，你能否作出一些理论上的小结，即针对什么类型的函数（具有某种单调特性或凹凸特性），用某种近似计算方法所得结果更接近于实际值？6．学习fulu2sum.m的程序设计方法，尝试用函数 sum 改写附录1和附录3的程序，避免for 循环．。
§3 最佳平方逼近多项式
若f ( x) C[a, b], span{0 ( x),1 ( x),n ( x)},
若函数组 0 ( x),1 ( x), n ( x)满足条件
则
方程组
0， 当i j ( i , j ) ( j , j )， 当i j
( 0 , 0 ) a 0 ( f , 0 ) a ( f , ) ( , ) 1 1 1 1 ( n , n ) an ( f , n )
* a 展开，而系数 k (k 1,2,, n)
按下式计算
ak ( f ( x),k ( x)) /(k ( x),k ( x)) ; (k 0,1,, n)
得级数
a k k ( x) k 0
* 称为f(x)的广义傅立叶(Foureir)级数,系数 ak (k 1,2,, n)
( x)
1 1 x2
，k 0 ( (Tk Tk ) ) ,k 0 2
1 1 x
2
由切比雪夫(Chebyshev)镇多项式作最佳平方逼近
Cn ( x ) ak Tk ( x ) k 0 n
其中
( f , T0 ) 1 a (T0 , T0 )
f P*
2 2
事实上,
f P
2 2
( f P, f P )
*
( f P* P* P, f Pn P* P)
( f P * , f P* ) ( P* P, P * P ) 2( f P * , P * P )
* 因为（f P ，P P） （f P ， （ a j a j） j ( x) ） * * * j 0
函数逼近研究
函数逼近研究一、课程目标知识目标：1. 让学生理解函数逼近的基本概念，掌握常用的函数逼近方法；2. 使学生掌握函数逼近的误差分析，并能够运用到实际问题的解决中；3. 引导学生运用逼近理论分析具体函数的性质，提高对函数的认知。
技能目标：1. 培养学生运用数学软件进行函数逼近实验的能力，提高数据处理和分析技能；2. 培养学生运用逼近方法解决实际问题的能力，提高解决问题的策略和方法；3. 提高学生的逻辑思维能力和团队合作能力，通过讨论、分享和交流，提升解决问题的效率。
情感态度价值观目标：1. 培养学生对数学学科的兴趣，激发他们探索函数逼近领域的热情；2. 培养学生严谨的科学态度，注重实证分析，增强学生的实证意识；3. 引导学生认识到数学知识在实际应用中的价值，提高学生的应用意识。
课程性质分析：本课程为高中数学选修课程，旨在通过函数逼近的研究，提高学生对数学知识的理解和应用能力。
学生特点分析：高中学生已经具备了一定的数学基础和逻辑思维能力，对新鲜事物充满好奇，但需要进一步引导和激发。
教学要求：1. 注重理论与实践相结合，提高学生的实际操作能力；2. 创设情境，引导学生主动探究，培养学生的自主学习能力；3. 注重团队合作，培养学生的沟通与协作能力。
二、教学内容1. 函数逼近的基本概念：- 逼近函数的定义与性质；- 常用的逼近函数类型及其特点。
2. 常用函数逼近方法：- 插值法：拉格朗日插值、牛顿插值；- 曲线拟合：最小二乘法、正交多项式；- 最佳逼近：最佳一致逼近、最佳平方逼近。
3. 函数逼近的误差分析：- 逼近误差的定义及计算方法；- 误差估计与控制策略。
4. 函数逼近的应用：- 实际问题中的函数逼近案例分析；- 数学软件在函数逼近中的应用。
5. 教学内容的安排与进度：- 第一章：函数逼近的基本概念（1课时）；- 第二章：常用函数逼近方法（2课时）；- 第三章：函数逼近的误差分析（1课时）；- 第四章：函数逼近的应用（2课时）。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
一. 任务：
用MATLAB语言编写连续函数最佳平方逼近的算法程序（函数式M文件）。并用
此程序进行数值试验，写出计算实习报告。
二. 程序功能要求：
在后面的附一leastp.m的基础上进行修改，使其更加完善。要求算法程序可以适应不同
的具体函数，具有一定的通用性。所编程序具有以下功能：
1. 用Lengendre多项式做基，并适合于构造任意次数的最佳平方逼近多项式。
可利用递推关系 0112()1,()()(21)()(1)()/2,3,.....nnnPxPxxPxnxPxnPxnn
2. 被逼近函数f(x)不用内联函数构造，而改用M文件建立数学函数。这样，此程序可通过
修改建立数学函数的M文件以适用不同的被逼近函数（要学会用函数句柄）。
3. 要考虑一般的情况]1,1[],[)(baxf。因此，程序中要有变量代换的功能。
4. 计算组合系数时，计算函数的积分采用变步长复化梯形求积法（见附三）。
5. 程序中应包括帮助文本和必要的注释语句。另外，程序中也要有必要的反馈信息。
6. 程序输入：（1）待求的被逼近函数值的数据点0x(可以是一个数值或向量)
（2）区间端点：a,b。
7. 程序输出：（1）拟合系数：012,,,...,ncccc
（2）待求的被逼近函数值
00001102200()()()()()nn
sxcPxcPxcPxcPx
三：数值试验要求：
1. 试验函数：()cos,[0,4]fxxxx；也可自选其它的试验函数。
2. 用所编程序直接进行计算，检测程序的正确性，并理解算法。
3. 分别求二次、三次、。。。最佳平方逼近函数）xs(。
4. 分别作出逼近函数）xs(和被逼近函数)(xf的曲线图进行比较。
（分别用绘图函数plot(0x,s(0x))和fplot(‘xcosx’,[x1 x2,y1,y2])）
四：计算实习报告要求：
1.简述方法的基本原理，程序功能，使用说明。
2.程序中要加注释。
3.对程序中的主要变量给出说明。
4.附源程序及计算结果。
一、 程序代码
1) legendre(N)函数程序
%legendre(N)函数
function P=legendre(N)
syms t x;%定义符号变量t x
for n=1:N
PP(n)=diff((t^2-1)^(n-1),n-1);%diff函数,求函数的n阶导数
Q(n)=2^(n-1)*prod([1:n-1]);%prod函数，计算数组元素的连乘积
end
PP(1)=1;
Q=sym(Q);
P=PP*(inv(diag(Q)));%inv函数，求逆
2) 采用M文件建立被逼近函数
%用M文件建立被逼近函数
function F=creat(x)
n=length(x);
F=x.*cos(x(1:n));%数组乘法运算
3) 区间变换函数程序
%区间变换函数程序
function f=convert(a,b,F)
syms x t; %定义符号变量t x
s=2\((b-a)*t+a+b); %实现区间转换
f=subs(F,x,s); %用s置换表达式中x，然后将置换完的表达式赋给f，符号替换
4) 变步长复化梯形求积公式程序
%变步长复化梯形求积公式
function I=tx(g)
m=1;
h=1-(-1);%积分区间
T=zeros(1,100);%赋予T初值为0的1*100行向量
T(1)=h*(feval(g,-1)+feval(g,1))/2;% feval函数，执行函数句柄
i=1;
while i<100%进行变步长计算
m=2*m;
h=h/2;
s=0;
for k=1:m/2
x=-1+h*(2*k-1);
s=s+feval(g,x);
end
T(i+1)=T(i)/2+h*s;
if abs(T(i+1)-T(i))<0.00001%设置精度值0.00001，终止条件
I=T(i+1);
break;
end
i=i+1;
end
5) 主程序
% 最佳平法逼近函数leastp
function [c s]= leastp(a,b,N)
syms t x;
F=creat(x); %用M文件建立的被逼近函数
P=legendre(N); %legendre(N)函数
f=convert(a,b,F); %区间变换函数程序
f=P*diag(f);
for i=1:N
g=inline(f(i));
I= tx(g); %变步长复化梯形求积公式
u(i)=I;
Q(i)=2\(2*(i-1)+1);
end
Q=sym(Q);
c=double(u*diag(Q));
S=c*P';
s=subs(S,t,(2*x-a-b)/(b-a));
subplot(211),
ezplot(s,[a:0.01:b]);
subplot(212),
ezplot(F,[a,b]);
二、 运行程序
实验函数为： ()cos,[0,4]fxxxx
输入：
>> a=0;b=4;N=2;
>> [c s]= leastp(a,b,N);
>>
>> a=0;b=4;N=4;
>> [c s]= leastp(a,b, N);
>>
>> a=0;b=4;N=7;
>> [c s]= leastp(a,b, N);
三、 实验结果
N=2
N=4
N=7
由上图可以发现：当阶数N取值不断增加时，图像的拟合程度越高。
附：
一、 参考程序
Lengendre多项式作基的函数最佳平方逼近算法程序LEASTP.m
（此程序只适用于对函数xxexf)( 构造最佳平方逼近多项式）
function [c,s]=leastp(x)
%LEASTP.m:least-square fitting with legendre polynomials
p1=1;
p2=inline('x','x');
p3=inline('(3*x^2-1)','x');
pp1=1;
pp2=inline('x.^2','x');
pp3=inline( '(3*x.^2-1)/2.*(3*x.^2-1)/2','x');
fp1=inline('x.*exp(x)','x');
fp2=inline('(x.^2).*exp(x)','x');
fp3=inline('(x.*exp(x)).*(3*x.^2-1)/2','x');
c(1)=quad(fp1,-1,1)/2;
c(2)=quad(fp2,-1,1)/quad(pp2,-1,1);
c(3)=quad(fp3,-1,1)/quad(pp3,-1,1);
s=c(1)+c(2)*p2(x)+c(3)*p3(x);
二、被逼近函数用M文件建立（例如下面）
function f=fun(x)
f=1./(1+x.^2);
三、变步长复化梯形求积公式的算法
1
11
1
1.,(()())22.0,23.(),3.15.*26.,7.,,2bahbaTfafbhHxaHHfxxxhxbTThHTTITIhhTT
4.若，则转
若则，输出,停机。
转2.
1
210[()()]211222((21))221,2,nnnn
njbaTfafbHTTTbabafajnnn