(五)实验题目:曲线拟合的最小二乘法
一、程序功能
编程实现用最小二乘法做曲线拟合
二、实验算例选择
y=sinx
三、算法
步骤1: 计算一组x,y的数据; 步骤2: 求出法方程
步骤3: 用递推公式拟合曲线方程
四、重要标识符说明
x 是插入点横坐标 y 为所求函数 M是拟合次数
五、程序运行实例
六、源程序
pi=3.1415926;
x=(pi/180)*linspace(0,2*180,2*18); y=sin(x); M=4;
[m,n]=size(x); w=ones(1,n);
15
c=zeros(1,n); a=zeros(1,n); b=zeros(1,n); p=zeros(M,n); p0=1;
c(1)=sum(w.*y)/sum(w); a(1)=sum(w.*x)/sum(w); p(1,:)=p0*(x-a(1));
c(2)=sum((w.*y).*p(1,:))/sum(w.*p(1,:).^2); s=c(1)+c(2)*(x-a(1));
a(2)=sum((w.*x).*p(1,:).^2)/sum(w.*p(1,:).^2); b(1)=sum(w.*p(1,:).^2)/sum(w.*p0^2); p(2,:)=(x-a(2)).*p(1,:)-b(1)*p0;
c(3)=sum((w.*y).*p(2,:))/sum(w.*p(2,:).^2); s=s+c(3)*p(2,:);
for k=3:M
a(k)=sum((w.*x).*p(k-1,:).^2)/sum(w.*p(k-1,:).^2); b(k-1)=sum(w.*p(k-1,:).^2)/sum(w.*p(k-2,:).^2); p(k,:)=(x-a(k)).*p(k-1,:)-b(k-1)*p(k-2,:); c(k+1)=sum((w.*y).*p(k,:))/sum(w.*p(k,:).^2); s=s+c(k+1)*p(k,:); end
plot(x,s,x,y,'*');
xlabel('\\theta (degree)'); ylabel('y');
legend('4次拟合曲线','sin()'); pause; close;
七、个人实验总结
最小二乘法只用递推公式,而且拟合次数越多越可以得到较好的所求方程。
16
(六)实验题目:复合梯形求积
一、程序功能
用复合梯形求积法求积
二、实验算例选择
编程实现区间为[0 1]对4/(1+x^2)的积分 三、算法
步骤1: 计算步长h;
步骤2: 计算各个结点对应的函数值;
步骤3: 将第二步所得结果代入复合梯形公式。
四、重要标识符说明
x 拟合的x的值
linspace(X1,X2,n) 生成n维向量,其中x1、x2、n分别为起始值、终止值、元素个数
n 为所取的x,y的长度 max() 取最大值
五、程序运行实例
:
六、源程序
function[I]=ComTrapz(n,x,y) n=8;
x=linspace(0,1,n) y=4./(1+(x).^2);
h=(max(x)-min(x))/(n-1); f(1:n)=y(1:n);
I=(f(1)+2*sum(f(2:n-1))+f(n))*h/2; %disp(I);
17
(六)实验题目:复合辛普森
一、程序功能
用复合辛普森求积法求积
二、实验算例选择
编程实现区间为[0 1]对4/(1+x^2)的积分 三、算法
步骤1:计算步长h;
步骤2:计算各个结点对应的函数值; 步骤3: 将第二
步所得结果n?1n?1Sn=h/3[f(a)+f(b)+2
?f(a?ih)+i?1?f(a?(2i?1)h)] i?0四、重要标识符说明
Sum为求和函数
五、程序运行实例
六、源程序
function[I]=ComSimpson(m,x,y) m=4; n=2*m+1;
x=linspace(0,1,n); y=4./(1+(x).^2);
h=(max(x)-min(x))/(n-1); f(1:n)=y(1:n);
I=(f(1)+4*sum(f(2: 2: n-1))+2*sum(f(3: 2: n-2))+f(n))*h/3; %disp(I);
(六)实验题目:龙贝格积分
18
代入
一、程序功能
用龙贝格求积法求积
二、实验算例选择
编程实现区间为[0 1]对4/(1+x^2)的积分
三、算法
步骤1: 计算初值T1,令h=(b-a)/21(i=0,1,2…)计算T2n; 步骤2: 求加速值Sn,Cn,Rn;
步骤3: 若满足精度则输出,否则转向步2。
四、重要标识符说明
a 给定下积分 b 给定上积分 t(1,1) 计算T1
五、程序运行实例
六、源程序
clc; clear; a=0; b=1; fa=a; fb=1/b;
19
err=10e-5; N=6;
t=zeros(N+1);
t(1,1)=(b-a)*(fa+fb)/2; j=1;
for i=j+1:N n=2^(i-2);
h=(b-a)*2^(2-i); for k=1:n
x(k)=a+(2*k-1)*h/2; end
y=4./(1+(x).^2); H=h*sum(y);
t(i,j)=(t(i-1,j)+H)/2; end
for n=j:N
for m=n:N-1
t(m+1,n+1)=((4^(n))*t(m+1,n)-t(m,n))/(4^(n)-1); end
if abs(t(n+1,n+1)-t(n,n)) disp(t(n,n)); 七、个人实验总结 三者比较龙贝格的效果是最好的。 20