龙贝格积分(c++)

用龙贝格算法计算积分
 

#include <iostream>

#include<cmath>

#include <iomanip>

using namespace std;

   

int power(int a, int b)

   

{

      int result = 1;

      if (b == 0)

            return result;

      while (b != 0)

      {

            result *= a;

            b--;

      }

      return result;

}

   

   

   

int main()

   

{

      int M = 5;            //M表示二分的次数

      double** T = new double*[M + 1];

      for (int i = 0; i < M + 1; i++)

            T[i] = new double[i + 1];

      double a = 0, b = 1;            //a表示积分下限,b表示积分上限

      double fa = pow(a, 1.5), fb = pow(b, 1.5);            //fa表示在a点处被积函数的函数值,fb表示在b点处被积函数的函数值

      double h = b - a;      //h表示迭代的步长

      double sum;      //用于计算积分内循环累加的部分

      T[0][0] = (fa + fb)*h / 2;      //用梯形公式进行初始化

      for (int i = 1; i < M + 1; i++)

   

      {

            h /= 2;

            sum = 0;

            for (int j = 1; j <= power(2, i - 1); j++)

                  sum += pow(a + (2 * j - 1)*h, 1.5);            //计算产生的新的结点部分的和

            T[i][0] = T[i - 1][0] / 2 + h*sum;                  //利用原来计算得到的积分值和新的节点计算得到新的积分值

            for (int k = 1; k <= i; k++)

                  T[i][k] = T[i][k - 1] + (T[i][k - 1] - T[i - 1][k - 1]) / (power(4, k) - 1);      //利用龙贝格积分公式计算后面的积分

      }

   

      

      cout << setiosflags(ios::left)

            << setw(10) << "T0"

            << setw(10) << "T1"

            << setw(10) << "T2"

            << setw(10) << "T3"

            << setw(10) << "T4"

            << setw(10) << "T5" << endl;

      for (int i = 0; i < M + 1; i++)

      {

            for (int j = 0; j <= i; j++)

                  cout << setw(10) << T[i][j];

            cout << endl;

      }

      cout << endl << endl;

      cout << "最后计算得到的积分值为:" << T[M][M] << endl;

      

      for (int i = 0; i < M + 1; i++)

            delete[] T[i];

      delete[] T;

      return 0;

   

}


posted @ 2015-12-04 01:53  硫酸亚铜  阅读(775)  评论(0编辑  收藏  举报