数值积分

最新推荐文章于 2023-01-01 11:05:01 发布

Brunon

最新推荐文章于 2023-01-01 11:05:01 发布

阅读量852

点赞数 1

分类专栏：数值分析文章标签：数值分析算法

本文链接：https://blog.csdn.net/brunon/article/details/66479797

版权

数值分析专栏收录该内容

5 篇文章 1 订阅

订阅专栏

本文探讨了基于插值方法的数值积分，重点讲解了Romberg算法，通过实验展示了其工作原理。文章提供了详细的代码实现，并阐述了数值积分在解决实际问题中的应用，尤其是在面对无法直接求解原函数或函数复杂时的重要性。

摘要由CSDN通过智能技术生成

一、写在前面
本篇文章介绍了以插值方法为基础的数值积分算法，其中最经典的便是Romberg提出的思想理念，同上篇文章一样，还是用一个实验来展示Romberg算法思想。
实验内容：
这里写图片描述

二、实验过程
在代码片中，我给出了较为详细的注释，在这里不在多说，最好的语言是代码！^V^
这里写图片描述

三、实验结果：
这里写图片描述
【参考代码】

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

double function(double x)    //原函数（可以不用给出，只给出几组已知数据放在数值中）
{
    if(x == 0)
        return 1;     //当x->0,函数极限是1
    return sin(x)/x;
}

double CompositeTrapezoidalRule(double a, double b, int n)     //复化梯形公式在区间[a, b]上n等分
{
    int i;
    double h = (b-a)/n;    //h为步长
    double sum = function(a) + function(b);
    for(i = 1; i<n; i++)    //运用公式
    {
        sum += 2*function(a + h*i);
    }
    return sum*h/2;
}

double CompositeSimpsonRule(double a, double b, int n)    //复化Simpson公式在区间[a, b]上n等分
{
    int i;
    double h = (b-a)/n;
    double sum = function(a) + function(b) + 4*function(a + h/2);
    for(i = 1; i<n; i++)
    {
        sum += 2*function(a + h*i);
        sum += 4*function(a + h*i + h/2);
    }
    return sum*h/6;
}

double CompositeCotesRule(double a, double b, int n)    //复化Cotes公式在区间[a, b]上n等分
{
    int i;
    double h = (b-a)/n;
    double sum = 7*(function(a) + function(b));
    sum += (32*function(a + h/4) + 12*function(a + h/2) + 32*function(a + 0.75*h));
    for(i = 1; i<n; i++)
    {
        sum += (32*function(a + h*i + h/4) + 12*function(a + h*i + h/2) + 32*function(a + h*i + 0.75*h));
    }
    return sum*h/90;
}

int main()
{
    double answer;
    printf("------------------Romberg Experiment--------------------\n\n");
    printf("k   Trapezoidal Simpson     Cotes     Romberg\n");
    printf("k=0: %lf \n", CompositeTrapezoidalRule(0.0, 1.0, 1));    //第一行显示

    printf("k=1: %lf \t", CompositeTrapezoidalRule(0.0, 1.0, 2));    //第二行显示
    printf("%lf\n", CompositeSimpsonRule(0.0, 1.0, 1));

    printf("k=2: %lf \t", CompositeTrapezoidalRule(0.0, 1.0, 4));    //第三行显示
    printf("%lf   ", CompositeSimpsonRule(0.0, 1.0, 2));
    printf("%lf\n", CompositeCotesRule(0.0, 1.0, 1));

    printf("k=3: %lf \t", CompositeTrapezoidalRule(0.0, 1.0, 8));    //第四行显示
    printf("%lf   ", CompositeSimpsonRule(0.0, 1.0, 4));
//    printf("%lf\n", CompositeCotesRule(0.0, 1.0, 2));    当n=2时，Cotes算法误差较大，此处用Romberg
    double S2 = CompositeSimpsonRule(0.0, 1.0, 2);
    double S4 = CompositeSimpsonRule(0.0, 1.0, 4);
    double C1 = CompositeCotesRule(0.0, 1.0, 1);
    double C2 = 16*S4/15 - S2/15;
    double R1 = 64*C2/63 - C1/63;
    printf("%lf   ", C2);
    printf("%lf\n", R1);

    printf("k=4: %lf \t", CompositeTrapezoidalRule(0.0, 1.0, 16));    //第五行显示
    printf("%lf   ", CompositeSimpsonRule(0.0, 1.0, 8));
//  C4和R2都用Romberg
    double S8 = CompositeSimpsonRule(0.0, 1.0, 8);
    double C4 = 16*S8/15 - S4/15;
    double R2 = 64*C4/63 - C2/63;
    printf("%lf   ", C4);
    printf("%lf\n", R2);

    printf("the 4 result is %lf\n", R2);
    return 0;
}

四、写在后面
数值积分的应用：在许多实际问题中，常常需要计算定积分的值。但是在实际使用中，往往遇到如下困难，而不能使用牛顿-莱布尼兹公式。
（1）找不到用初等函数表示的原函数
（2）虽然找到了原函数，但因表达式过于复杂而不便计算
（3） f(x)是由测量或计算得到的表格函数
由于以上种种困难，有必要研究积分的数值计算问题。