一重积分的C++实现

发布时间:2020-06-21 13:40:27 作者:Dussssss
来源:网络 阅读:1249
#include<iostream>  
#include<iomanip>  
#include<cmath>  

using namespace std;  

//求积:梯形公式  
double Func_Integral_Trapezoid  
(double lo, double hi, double(*Func)(double), int n = 1000)  
{//lo:下限;hi:上限;Func:函数;n:等分数  
   if (n <= 0) n = 100;  
      double x;  
    double step = (hi - lo) / n;  

    double result = 0.0;  
    x = lo;  
    for (int i = 1; i < n; i++) {  
        x += step;  
        result += Func(x);  
    }  
    result += (Func(lo) + Func(hi)) / 2;  
    result *= step;  

    return result;  
}  

//求积:Simpson公式  
double Func_Integral_Simpson  
(double lo, double hi, double(*Func)(double), int n = 1000)  
{//lo:下限;hi:上限;Func:函数;n:等分数  

    if (n <= 0) n = 100;  

    double x;  
    double step = (hi - lo) / n;  

    double result1 = 0.0;  
    x = lo;  
    for (int i = 1; i < n; i++) {  
        x += step;  
        result1 += Func(x);  
    }  
    result1 *= 2;  

    double result2 = 0.0;  
    x = lo + step / 2;  
    for (int i = 0; i < n; i++) {  
        result2 += Func(x);  
        x += step;  
    }  
    result2 *= 4;  

    double result = result1 + result2 + Func(lo) + Func(hi);  
    result *= step / 6;  

    return result;  
}  

//求积:Cotes公式  
double Func_Integral_Cotes  
(double lo, double hi, double(*Func)(double), int n = 1000)  
{//lo:下限;hi:上限;Func:函数;n:等分数  

    if (n <= 0) n = 100;  

    double x;  
    double step = (hi - lo) / n;  

    double result1 = 0.0;  
    x = lo;  
    for (int i = 1; i < n; i++) {  
        x += step;  
        result1 += Func(x);  
    }  
    result1 *= 14;  

    double result2 = 0.0;  
    x = lo + step / 2;  

    double result3 = 0.0;  
    double x1 = lo + step / 4;  
    double x2 = lo + step / 4 * 3;  

    for (int i = 0; i < n; i++) {  
        result2 += Func(x);  
        result3 += Func(x1) + Func(x2);  
        x += step;  
        x1 += step;  
        x2 += step;  
    }  
    result2 *= 12;  
    result3 *= 32;  

    double result4 = (Func(lo) + Func(hi)) * 7;  

    double result = result1 + result2 + result3 + result4;  
    result *= step / 90;  

    return result;  
}  

//求积:Romberg公式  
double Func_Integral_Romberg  
(double lo, double hi, double(*Func)(double), int k = 4)  
{//lo:下限;hi:上限;Func:函数;k:等分指数(区间等分成2^k份)  
    int size = k + 1;  
    double *matrix = new double[size*size];  
    for (int i = 0; i < size*size; i++) matrix[i] = 0.0;  
    double step = hi - lo;  
    matrix[0] = Func_Integral_Trapezoid(lo, hi, Func, 1);  
    for (int i = 1; i < size; i++) {  
        int n = 1 << (i - 1);  
        for (int k = 0; k < n; k++) {  
            matrix[i*size + 0] += Func(lo + (k + 0.5)*step);  
        }  
        matrix[i*size + 0] *= step;  
        matrix[i*size + 0] += matrix[(i - 1)*size];  
        matrix[i*size + 0] /= 2.0;  
        step /= 2.0;  
    }  

    double temp = 1.0;  
    double factor1, factor2;  
    for (int j = 1; j < size; j++) {  
        temp *= 4.0;  
        factor1 = temp / (temp - 1);  
        factor2 = 1 / (temp - 1);  
        for (int i = j; i < size; i++) {  
            matrix[i*size + j] = factor1*matrix[i*size + j - 1]  
                - factor2*matrix[(i - 1)*size + j - 1];  
        }  
    }  

    double result = matrix[k*size + k];  
    delete[] matrix;  

    return result;  
}  

//测试用的被积函数,0到1积分为PI  
double Func_test1(double x)  
{  
    return 4 / (1 + x*x);  
}  

int main() {  
    cout << "Trapezoid Numerical Integration" << endl;  
    cout << setprecision(15) << Func_Integral_Trapezoid(0, 1, Func_test1, 1024) << endl << endl;  
    cout << "Simpson Numerical Integration" << endl;  
    cout << setprecision(15) << Func_Integral_Simpson(0, 1, Func_test1, 1024) << endl << endl;  
    cout << "Cotes Numerical Integration" << endl;  
    cout << setprecision(15) << Func_Integral_Cotes(0, 1, Func_test1, 1024) << endl << endl;  
    cout << "Romberg Numerical Integraion" << endl;  
    cout << setprecision(15) << Func_Integral_Romberg(0, 1, Func_test1, 10) << endl << endl;  
    getchar();  
    return 0;  
}  
推荐阅读:
  1. Python如何实现北京积分落户数据
  2. 怎么用python求一重积分和二重积分

免责声明:本站发布的内容(图片、视频和文字)以原创、转载和分享为主,文章观点不代表本网站立场,如果涉及侵权请联系站长邮箱:is@yisu.com进行举报,并提供相关证据,一经查实,将立刻删除涉嫌侵权内容。

c++ 实现 c+

上一篇:自定义django模板的 tags和filters

下一篇:java如何获取键盘输入值

相关阅读

您好,登录后才能下订单哦!

密码登录
登录注册
其他方式登录
点击 登录注册 即表示同意《亿速云用户服务条款》