《计算方法实验报告(附代码)》由会员分享，可在线阅读，更多相关《计算方法实验报告(附代码)（37页珍藏版）》请在金锄头文库上搜索。
1、1实验一 牛顿下山法实验说明：求非线性方程组的解是科学计算常遇到的问题，有很多实际背景各种算法层出不穷，其中迭代是主流算法。只有建立有效的迭代格式，迭代数列才可以收敛于所求的根。因此设计算法之前，对于一般迭代进行收敛性的判断是至关重要的。牛顿法也叫切线法，是迭代算法中典型方法，只要初值选取适当，在单根附近，牛顿法收敛速度很快，初值对于牛顿迭代至关重要。当初值选取不当可以采用牛顿下山算法进行纠正。牛顿下山公式： )(1kkxfx下山因子 L， 32下山条件 |)(|)(|1kkxff实验代码：#include#include#includeusing namespace std;2double
2、newton_downhill(double x0,double x1); /牛顿下山法函数，返回下山成功后的修正初值 double Y; /定义下山因子 Ydouble k; /k 为下山因子 Y允许的最小值double dfun(double x)return 3*x*x-1; /dfun()计算 f(x)的导数值double fun1(double x)return x*x*x-x-1; /fun1()计算 f(x)的函数值double fun2(double x) return x-fun1(x)/dfun(x); /fun2()计算迭代值int N; /N 记录迭代次数double
3、e; /e 表示要求的精度int main()double x0,x1;coutx0;coute;N=1;if(dfun(x0)=0)cout=fabs(fun1(x0) /初值不满足要求时,转入牛顿下山法x1=newton_downhill(x0,x1); /牛顿下山法结束后,转入牛顿迭代法进行计算while(fabs(x1-x0)=e) /当精度不满足要求时N=N+1;x0=x1;if(dfun(x0)=0)coutk;while(fabs(fun1(x1)=fabs(fun1(x0)if(Yk)Y=Y/2;else cout#include #include using namespac
4、e std;class gsdlprivate: int n;5double *a, *b, *x, eps;public:gsdl (int nn)int i;n = nn;a = new double*n; /动态分配内存空间for (i=0; istr1;ifstream fin (str1);if (!fin) cout aij;for (i=0; ibi; /读入常数向量 Bfin eps; /读入 efin.close ();void gsdl:a_gsdl () int i,j;double p,t,s,q;for (i=0; i=fabs(aii)cout =eps) p=0.
5、0;for (i=0; ip) p=q;void gsdl:output () /输出结果到文件并显示int i;char str220;cout str2;ofstream fout (str2);if (!fout) cout #include #include using namespace std;class gaus private: int n;double *a, *b;public:gaus (int nn)int i;n = nn;a = new double*n; /动态分配内存空间for (i=0; istr1;ifstream fin (str1);if (!fin)
6、cout aij;for (i=0; ibi; /读入常数向量 Bfin.close ();void gaus:gauss () int *js,l,k,i,j,is;9double d,t;js = new intn;l=1;for (k=0; kd) d=t; jsk=j; is=i; if (d+1.0=1.0) l=0;else if (jsk!=k)for (i=0;i=0;i-) t=0.0;for (j=i+1;j=0;k-)if (jsk!=k) t=bk; bk=bjsk; bjsk=t;delete js;void gaus:output () /输出结果到文件并显示int
7、 i;char str220;cout str2;ofstream fout (str2);if (!fout) cout #include #include #include using namespace std;class splin3private: int n, m;double *x, *y, *dy, *ddy, *t, *z, *dz, *ddz, integ;double *s;public:splin3 (int nn, int mm)n = nn; m = mm;x = new doublen; /动态分配内存y = new doublen;dy = new double
8、n;ddy = new doublen;s = new doublen;t = new doublem;z = new doublem;dz = new doublem;ddz = new doublem;void input (); /由文件读入 n 个数据点(x, y)以及 m 个指定插值点 tvoid interp (); /执行三次样条函数插值、微商与积分12void output (); /输出 n 个数据点及其一阶与二阶导数值 m 个插值点 t 处的函数值及其一阶与二阶导数值与积分值到文件并显示splin3 () delete x, y, dy, ddy, t, z, dz, dd
9、z, s; ;void splin3:input () /由文件读入 n 个数据点(x, y)以及 m 个指定插值点 tint k;ifstream fin (ain.txt);if (!fin) cout xk; fin yk; for (k=0; ktk; /读入 m 个插值点fin.close ();void splin3:interp () /执行三次样条函数插值、微商与积分 int i,j;double h0,y0,h1,y1,alpha,beta,u;h0=xn-1-xn-2;y0=yn-1-yn-2;dy0=0.0; ddy0=0.0; ddyn-1=0.0;s0=1.0; sn
10、-1=1.0;for (j=1;j=1;j-)13 sj=dyj*sj+1+sj;ddyj=dyj*ddyj+1+ddyj;dyn-2=(beta-alpha*ddy1-(1.0-alpha)*ddyn-2)/(alpha*s1+(1.0-alpha)*sn-2+2.0);for (j=2;j=xn-1) h0=h0-(xn-1-x0);while (h0xi+1) i=i+1;u=h0;h1=(xi+1-u)/si;h0=h1*h1;zj=(3.0*h0-2.0*h0*h1)*yi;zj=zj+si*(h0-h0*h1)*dyi;dzj=6.0*(h0-h1)*yi/si;dzj=dzj+(
11、3.0*h0-2.0*h1)*dyi;ddzj=(6.0-12.0*h1)*yi/(si*si);ddzj=ddzj+(2.0-6.0*h1)*dyi/si;h1=(u-xi)/si;h0=h1*h1;zj=zj+(3.0*h0-2.0*h0*h1)*yi+1;zj=zj-si*(h0-h0*h1)*dyi+1;dzj=dzj-6.0*(h0-h1)*yi+1/si;dzj=dzj+(3.0*h0-2.0*h1)*dyi+1;ddzj=ddzj+(6.0-12.0*h1)*yi+1/(si*si);ddzj=ddzj-(2.0-6.0*h1)*dyi+1/si;14void splin3:ou
12、tput () /输出 n 个数据点及其一阶与二阶导数值 m 个插值点 t处的函数值及其一阶与二阶导数值与积分值到文件并显示int i;ofstream fout (aout.txt);if (!fout) cout #include #include using namespace std;class pir1private: int n;double *x, *y,a,b;public:pir1 (int nn)0021nni iii iyaxx16n = nn;x = new doublen; /动态分配内存y = new doublen;void input (); /由文件读入 n
13、个数据点(x, y)void fit (); /执行最小二乘直线拟合void output (); /输出拟合直线pir1 () delete x, y; ;void pir1:input ()int k;ifstream fin (bin.txt);if (!fin) cout xk;fin yk; fin.close ();void pir1:fit () int k;double s1,s2,s3,s4;for(s1=0,s2=0,s3=0,s4=0,k=0;k0)cout #include #include using namespace std;class romb private: double a, b, eps, integ;public:romb (double aa, double bb, double es)/顺序提供 a,b,eps 值的构造函数 a = aa; b = bb; eps = es;
