数值分析/计算方法 实验(C 或 Matlab) 拉格朗日插值/埃尔米特插值/最小二乘法/复化求积公式
数值分析/计算方法
Lagrange插值多项式
实验要求和提示
实验代码(C·无画图)
#define N 13 #includeusing namespace std; int main() { double x[N]={0,10,20,30,40,50,60,70,80,90,100,110,120}; double y[N]={5,1,7.5,3,4.5,8.8,15.5,6.5,-5,-10,-2,4.5,7}; double l; double X=65; double Y=0; int i,k; for(k=0;k ) { l=1; for(i=0;i ) if(i!=k) l=l*((X-x[i])/(x[k]-x[i])); Y=Y+y[k]*l; } cout<<Y; }
埃尔米特插值
实验要求和提示
实验代码 (C·无画图)
#define N 10 #includeusing namespace std; int main() { double x[N]={0.10,0.20,0.30,0.40,0.50,0.60,0.70,0.80,0.90,1.00}; double y[N]={0.904837,0.818731,0.740818,0.670320,0.606531,0.548812,0.496585,0.449329,0.406570,0.367879}; double m[N]={-0.904837,-0.818731,-0.740818,-0.670320,-0.606531,-0.548812,-0.496585,-0.449329,-0.406570,-0.367879}; double l; double ldao[N]={0}; //ldao(x)为l(x)的导数 double X=0.55; double Y=0; int j,k,i; for(j=0;j ) { l=1; for(k=0;k ) if(k!=j) { ldao[j]=ldao[j]+(1/(x[j]-x[k])); l=l*((X-x[k])/(x[j]-x[k])); } Y=Y+y[j]*((1-2*(X-x[j])*ldao[j])*l*l)+m[j]*((X-x[j])*l*l); } cout<<Y; return 0; }
最小二乘法
实验要求和提示
实验代码(Matlab·有画图)
x = [0 10 20 30 40 50 60 70 80 90]; y = [68 67.1 66.4 65.6 64.6 61.8 61.0 60.8 60.4 60]; xifsum = 0; yisum = 0; xisum = 0; xycsum = 0; for i = 1:10 xifsum = x(i) * x(i) + xifsum; yisum = y(i) + yisum; xisum = x(i) + xisum; xycsum = x(i) * y(i) + xycsum; end m = 9; b = (xifsum*yisum - xisum*xycsum) / ((m+1)*xifsum - xisum*xisum); a = ((m+1)*xycsum - xisum*yisum) / ((m+1)*xifsum - xisum*xisum); X = 55; S = a * X + b; disp('X的值为55时,表达式的值为'); disp(S); X = 0 : 1 : 90; yy = a * X + b; figure plot(X, yy, 'k') for j = 1:10 text(x(j), y(j), '*', 'color', 'r'); end text(55, S, '*', 'color', 'b');
实验演示
复化求积公式
实验代码(C)
#include#include #include using namespace std; double function(double x); double T(double a,double b,double h,int N); double Fh(double a,double b,double h,int N); void getIntegralValue(double a,double b,int N);//复化梯形求积 int main() { double a,b ;//产生的节点数 int N=1; cout<<"请输入左右区间范围:"; cin>>a; cin>>b; //getIntegralValue(a,b,N); double t,TN,Tn,Hn=0,h=b-a; t=T(a,b,h,N); do { Hn=Fh(a,b,h,N); Tn=t; TN=0.5*(T(a,b,h,N)+Hn); N=2*N; h=h/2; t=TN; }while(fabs(TN-Tn)>=0.000001); cout<<"Tn="< 8)< endl; cout<<"等分数为:"< 2<<endl; return 0; } double function(double x) { return exp(x)*cos(x); } double T(double a,double b,double h,int N) { double f=0; for(int k=1;k<=N-1;k++) f+=function(a+k*h); return (h/2)*(function(a)+2*f+function(b)); } double Fh(double a,double b,double h,int N) { double Hn=0; for(int i=1;i<=N;i++) { h=(b-a)/N; Hn+=function(a+(h*(2*i-1)/2)); } return Hn=h*Hn; }