码迷,mamicode.com
首页 > 编程语言 > 详细

龙贝格算法

时间:2015-11-13 00:55:32      阅读:243      评论:0      收藏:0      [点我收藏+]

标签:

求积分的龙贝格算法

计算f(x)=1/x在[1,3]上的积分:

 1 #include <iostream>
 2 #include <cstring>
 3 #include <cmath>
 4 using namespace std;
 5 
 6 double f(double x){
 7     return 1.0/x;
 8 }
 9 
10 int main(){
11     double a,b;//积分区间的上界和下届 
12     a=1;
13     b=3;
14     
15     double t[100];
16     memset(t,-1,sizeof(t));
17     t[1]=(b-a)/2.0*( f(a)+f(b) ) ;
18     for(int n=1;n<=8;n=2*n )
19     {
20         double temp=0;
21         for(int i=0;i<=n-1;i++){
22             temp+=f( a+(2*i+1)*(b-a)/(2*n)  );
23         }
24         temp=temp*(b-a)/(2*n);
25         t[2*n]=t[n]/2.0 + temp;
26     }
27     
28     for(int i=1;i<=16;i=i*2){
29         printf("t(%d)=%2.6f\n",i,t[i]);
30     }
31     printf("\n");
32     
33     double s[100];
34     memset(s,-1,sizeof(s));
35     for(int n=1;n<=8;n*=2)
36     {
37         s[n]=4.0/3.0*t[2*n]-t[n]/3.0;
38     }
39     for(int i=1;i<=8;i=i*2){
40         printf("s(%d)=%2.6f\n",i,s[i]);
41     }
42     printf("\n");
43     
44     double c[100];
45     memset(c,-1,sizeof(c));
46     for(int n=1;n<=4;n*=2)
47     {
48         c[n]=16.0/15.0*s[2*n]-s[n]/15.0;  
49     } 
50     for(int i=1;i<=4;i=i*2){
51         printf("c(%d)=%2.6f\n",i,c[i]);
52     }
53     printf("\n");
54     
55     double r[100];
56     memset(r,-1,sizeof(r));
57     for(int n=1;n<=2;n*=2)
58     {
59         r[n]=64.0/63.0*c[2*n]-c[n]/63.0;
60     }
61     for(int i=1;i<=2;i=i*2){
62         printf("r(%d)=%2.8f\n",i,r[i]);
63     }
64     
65     
66 }

输出结果为:

技术分享

龙贝格算法

标签:

原文地址:http://www.cnblogs.com/liugl7/p/4960706.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!