数值计算day8-数值积分

上节课主要介绍了计算微分的几种数值方法,对一阶微分,最简单的莫过于两点前向差分、后向差分和中心差分这三种方法,其中中心差分的精度最高,这三种差分公式都可以通过推导泰勒展开式得到,而通过泰勒展开式还可以推导出三点前向差分和三点后向差分。对二阶微分,则可以推导出三点中心差分、三点前向差分、三点后向差分公式。这些数值方法均可以推广到数值偏微分的领域,且对于两个精度不高的数值算法,还可以使用Richardson外推加速算法得到一个精度更高的算法。本节课主要介绍计算积分的数值方法。

1. 矩形和中点法

image.png

区间上的积分表示的是曲线在内的面积,若将分解为个小区间:,区间内的面积可以估算为一个黎曼和:,其中表示内的某个点,,此时积分即为黎曼求和的极限: 其中

1.1 左点法与右点法

image.png

左点法是一种矩形法,每个小区间的面积使用矩形公式来计算,其中高取为区间左端点的函数值,宽取为小区间的长度,
右点法与左点法相对,高取为右端点的函数值,

1.2 中点法

image.png

中点法使用区间的中点函数值作为高,区间的面积为:
下面以线性函数为例,讨论如何控制左点法的结果在想要的误差范围内。区间取为,在该区间上积分的真实值为,利用左点法计算,将等分为个小区间,每个小区间的宽度为,左点法计算的积分值为: 计算误差为,如果想要将误差的量级控制在以内,则 至少需要将区间分为个小区间。而在该例中,中点法计算是没有误差的。

2. 梯形法(Trapezoidal Method)

image.png

梯形法是对矩形法的一种改善,采用线性函数来计算积分,对点(a,f(a))(b,f(b))进行牛顿插值,得到:f(x)=f(a)+(x-a)f[a,b]=f(a)+(x-a)\frac{f(b)-f(a)}{b-a} 将区间[a,b]的积分值计算为梯形面积,有:I(f)=f(a)(b-a)+\frac{1}{2}[f(b)-f(a)](b-a)=\frac{f(a)+f(b)}{2}(b-a) 对区间进行划分a=x_0<x_1<x_2<\cdots<x_n=b,则:I_i(f)=\int_{x_{i-1}}^{x_i}f(x)dx=\frac{f(x_{i-1})+f(x_i)}{2}(x_i-x_{i-1}) I(f)=\sum^n_{i=1}I_i(f)=\sum^n_{i=1}\frac{f(x_{i-1})+f(x_i)}{2}(x_i-x_{i-1})=\frac{h}{2}[f(a)+f(b)]+h\sum^{n-1}_{i=2}f(x_i)

image.png

MATLAB实现:

function I = trapezoidal(Fun,a,b,N)
% trapezoidal numerically integrate using the Composite Trapezoidal Method.
% Input Variables:
% Fun Name for the function to be integrated.
% (Fun is assumed to be written with element-by-element calculations.).
% a  Lower limit of integration.
% b  Upper limit of integration.
% N  Number of subintervals.
% Output Variable:
% I  Value of the integral.
 
h = (b-a)/N;
x=a:h:b;
F=Fun(x);
I=h*(F(1)+F(N+1))/2+h*sum(F(2:N));

3. Simpson 方法

3.1 Simpson 1/3 法

梯形法是采用线性函数插值,来估算区间面积,一种改进方法是使用非线性插值,包括二次插值、三次插值。Simpson 1/3 法采用二次插值法。

image.png

取区间内的三点,使用牛顿多项式插值法,过这三个点的曲线方程为: 求解可得 , 使用曲线代替原曲线计算积分,
现在将区间等分成(偶数)个小区间:,则区间内的积分值为:
image.png

对所有小区间积分值求和,有:

3.2 Simpson 3/8 法

在Simpson 3/8 法中,对每四个点进行三次插值:p(x)=c_3x^3+c_2x^2+c_1x+c_0

image.png
小区间的积分值为: 因此将区间等分为(3的倍数)个小区间,对每个小区间进行三次插值计算积分后求和,可得:

4. 高斯求积

前面几种算法一般要求小区间的宽度相同,即将区间等距划分,高斯法则不要求等距,其一般形式为:\int^b_{a}f(x)dx \approx \sum^n_{i=1}C_if(x_i) 其中C_i为点x_i对应的权重,点x_i是待定的高斯点。现在考虑[-1,1]区间:\int^1_{-1}f(x)dx \approx \sum^n_{i=1}C_if(x_i)n=2的情况,有:\int^1_{-1}1dx=C_1+C_2=2 \int^1_{-1}xdx=C_1x_1+C_2x_2=0 \int^1_{-1}x^2dx=C_1x^2_1+C_2x^2_2=\frac{2}{3} \int^1_{-1}x^3dx=C_1x^3_1+C_2x^3_2=0 这是一组非线性方程,存在无数个解,因此,额外提出对称性要求,即x_1=-x_2,则可求解得:C_1=1,C_2=1,x_1=-\frac{1}{\sqrt{3}}=-0.57735027,x_2=\frac{1}{\sqrt{3}}=0.57735027\int^1_{-1}f(x)dx \approx = f(-\frac{1}{\sqrt{3}})+ f(\frac{1}{\sqrt{3}}) 同样地,可以推导出n=3的情况,当n越大时,计算的精度越高。下表给出了n=2,3,4,5,6时,对应的高斯点和权重:

image.png

上述推导过程是以区间为例,对积分区间,可以通过如下线性变换,将积分区间变换为:

5. 总结

本节课主要介绍函数积分的数值算法,最简单的有矩形法(左右点法)和中点法,在中点法的基础上,线性插值法提供了一种梯形求解算法,而采用非线性插值,则可以推导出Simpson 1/3 法和 3/8 法。这几种算法都要求对区间等距划分,最后介绍的高斯求积法则是确定高斯点和相应的权重来直接计算积分值。与数值微分类似,在数值积分算法中,也可以用Richardson外推加速算法将两个精度不高的方法组合构建出一个精度更高的数值积分算法,这里不做详细讨论,包括多重积分、多元积分等相关知识,这里暂不做深入研究,如果后面要使用,再进一步讨论。至此,暑期学校数值方法与计算课程的内容已经整理完毕。原本打算两周内整理好,后来因为各种事情耽搁了,以后尽量避免长线作战,趁热打铁方能加深印象。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 205,724评论 6 479
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 88,104评论 2 382
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 152,142评论 0 341
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 55,086评论 1 278
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 64,076评论 5 370
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,914评论 1 283
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,220评论 3 399
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,871评论 0 259
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 43,318评论 1 300
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,834评论 2 323
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,951评论 1 333
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,574评论 4 322
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,162评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,162评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,383评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 45,349评论 2 352
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,652评论 2 343

推荐阅读更多精彩内容

  • 何谓数值分析? 众所周知现实生活中科学技术上面的问题大多数都能够通过建立对应的数学模型把实际问题转化成一个数学问题...
    罗泽坤阅读 8,774评论 0 14
  • 上一节课主要介绍了曲线拟合与插值,曲线拟合主要包括线性拟合(单特征线性回归和非线性拟合(非线性方程特征变换、高阶多...
    xkzhai阅读 2,615评论 0 0
  • 一、前言 本文主要讨论定积分的数值解,即数值积分。为什么需要数值积分:假设在区间上可积,为其原函数;现实中的多数情...
    胜负55开阅读 6,436评论 0 9
  • 研究生时,极其不喜欢数值分析,觉得太抽象。 结果研究生也没研究出什么,没有能够用好数学这个工具,这也成了心里一个结...
    建_5f75阅读 960评论 0 2
  • 1 数值积分概述 1.1 引言   对于许多实际问题的求解往往需要计算积分。在高等数学中计算积分采用的是著名的牛顿...
    Bocchi阅读 4,379评论 0 5