C数值计算算法编程数值积分类Integral的设计.docx

1、C数值计算算法编程数值积分类Integral的设计C#数值计算算法编程数值积分类Integral的设计namespace CSharpAlgorithm.Algorithm /* * 计算数值积分的类Integral * * author 周长发 * version 1.0 */ public abstract class Integral /* * 抽象函数：计算积分函数值，必须在派生类中覆盖该函数 * * param x - 函数变量 * return double型，对应的函数值 */ public abstract double Func(double x); /* * 基本构造函数

2、*/ public Integral() /* * 变步长梯形求积法 * * 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x) * * param a - 积分下限 * param b - 积分上限，要求ba * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValueTrapezia(double a, double b, double eps) int n,k; double fa,fb,h,t1,p,s,x,t=0; / 积分区间端点的函数值 fa=Func(a); fb=Fun

3、c(b); / 迭代初值 n=1; h=b-a; t1=h*(fa+fb)/2.0; p=eps+1.0; / 迭代计算 while (p=eps) s=0.0; for (k=0;ka * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValueSimpson(double a, double b, double eps) int n,k; double h,t1,t2,s1,s2=0,ep,p,x; / 计算初值 n=1; h=b-a; t1=h*(Func(a)+Func(b)/2.0; s1=t1; ep=ep

4、s+1.0; / 迭代计算 while (ep=eps) p=0.0; for (k=0;ka * param d - 对积分区间进行分割的最小步长，当子区间的宽度 * 小于d时，即使没有满足精度要求，也不再往下进行分割 * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValueATrapezia(double a, double b, double d, double eps) double h,f0,f1,t0,z; double t = new double2; / 迭代初值 h=b-a; t0=0.0; f0

5、=Func(a); f1=Func(b); t0=h*(f0+f1)/2.0; / 递归计算 ppp(a,b,h,f0,f1,t0,eps,d,t); z=t0; return(z); /* * 内部函数 */ private void ppp(double x0, double x1, double h, double f0, double f1, double t0, double eps, double d, double t) double x,f,t1,t2,p,g,eps1; x=x0+h/2.0; f=Func(x); t1=h*(f0+f)/4.0; t2=h*(f+f1)/4

6、.0; p=Math.Abs(t0-(t1+t2); if (peps)|(h/2.0a * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValueRomberg(double a, double b, double eps) int m,n,i,k; double h,ep,p,x,s,q=0; double y = new double10; / 迭代初值 h=b-a; y0=h*(Func(a)+Func(b)/2.0; m=1; n=1; ep=eps+1.0; / 迭代计算 while (ep=eps)&(

7、m=9) p=0.0; for (i=0;i=n-1;i+) x=a+(i+0.5)*h; p=p+Func(x); p=(y0+h*p)/2.0; s=1.0; for (k=1;ka * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValuePq(double a, double b, double eps) int m,n,k,l,j; double hh,t1,s1,ep,s,x,t2,g=0; double h = new double10; double bb = new double10; / 迭代初值

8、 m=1; n=1; hh=b-a; h0=hh; t1=hh*(Func(a)+Func(b)/2.0; s1=t1; bb0=s1; ep=1.0+eps; / 迭代计算 while (ep=eps)&(m=9) s=0.0; for (k=0;k=n-1;k+) x=a+(k+0.5)*hh; s=s+Func(x); t2=(t1+hh*s)/2.0; m=m+1; hm-1=hm-2/2.0; g=t2; l=0; j=2; while (l=0)&(j=2;j-) g=bbj-2-hj-2/g; ep=Math.Abs(g-s1); s1=g; t1=t2; hh=hh/2.0;

9、n=n+n; return(g); /* * 高振荡函数求积法 * * 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x) * * param a - 积分下限 * param b - 积分上限，要求ba * param m - 被积函数中振荡函数的角频率 * param n - 给定积分区间两端点上的导数最高阶数 * param fa - 一维数组，长度为n，存放f(x)在积分区间端点x=a处的各阶导数值 * param fb - 一维数组，长度为n，存放f(x)在积分区间端点x=b处的各阶导数值 * param s - 一维数组，长度为，其中s(1)返回f(

10、x)cos(mx)在积分区间的积分值， * s(2) 返回f(x)sin(mx)在积分区间的积分值 * return double 型，积分值 */ public double GetValuePart(double a, double b, int m, int n, double fa, double fb, double s) int mm,k,j; double sma,smb,cma,cmb; double sa = new double4; double sb = new double4; double ca = new double4; double cb = new doubl

11、e4; / 三角函数值 sma=Math.Sin(m*a); smb=Math.Sin(m*b); cma=Math.Cos(m*a); cmb=Math.Cos(m*b); / 迭代初值 sa0=sma; sa1=cma; sa2=-sma; sa3=-cma; sb0=smb; sb1=cmb; sb2=-smb; sb3=-cmb; ca0=cma; ca1=-sma; ca2=-cma; ca3=sma; cb0=cmb; cb1=-smb; cb2=-cmb; cb3=smb; s0=0.0; s1=0.0; mm=1; / 循环迭代 for (k=0;k=4) j=j-4; mm=

12、mm*m; s0=s0+(fbk*sbj-fak*saj)/(1.0*mm); s1=s1+(fbk*cbj-fak*caj)/(1.0*mm); s1=-s1; return s0; /* * 勒让德高斯求积法 * * 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x) * * param a - 积分下限 * param b - 积分上限，要求ba * param eps - 积分精度要求 * return double 型，积分值 */ public double GetValueLegdGauss(double a, double b, double e

13、ps) int m,i,j; double s,p,ep,h,aa,bb,w,x,g=0; / 勒让德高斯求积系数 double t=-0.9061798459,-0.5384693101,0.0, 0.5384693101,0.9061798459; double c=0.2369268851,0.4786286705,0.5688888889, 0.4786286705,0.2369268851; / 迭代初值 m=1; h=b-a; s=Math.Abs(0.001*h); p=1.0e+35; ep=eps+1.0; / 迭代计算 while (ep=eps)&(Math.Abs(h)

14、s) g=0.0; for (i=1;i=m;i+) aa=a+(i-1.0)*h; bb=a+i*h; w=0.0; for (j=0;j=4;j+) x=(bb-aa)*tj+(bb+aa)/2.0; w=w+Func(x)*cj; g=g+w; g=g*h/2.0; ep=Math.Abs(g-p)/(1.0+Math.Abs(g); p=g; m=m+1; h=(b-a)/m; return(g); /* * 拉盖尔高斯求积法 * * 调用时，须覆盖计算函数f(x)值的虚函数double Func(double x) * * return double 型，积分值 */ public

15、double GetValueLgreGauss() int i; double x,g; / 拉盖尔高斯求积系数 double t=0.26355990, 1.41340290, 3.59642600, 7.08580990, 12.64080000; double c=0.6790941054, 1.638487956, 2.769426772, 4.315944000, 7.104896230; / 循环计算 g=0.0; for (i=0; i=4; i+) x=ti; g=g+ci*Func(x); return(g); /* * 埃尔米特高斯求积法 * * 调用时，须覆盖计算函数f

16、(x)值的虚函数double Func(double x) * * return double 型，积分值 */ public double GetValueHermiteGauss() int i; double x,g; / 埃尔米特高斯求积系数 double t=-2.02018200, -0.95857190, 0.0,0.95857190, 2.02018200; double c=1.181469599, 0.9865791417, 0.9453089237, 0.9865791417, 1.181469599; / 循环计算 g=0.0; for (i=0; i=4; i+) x=ti; g=g+ci*Func(x); return(g);