c语言实现定积分求解

1. #include <stdio.h>  
2. #include <stdlib.h>  
3. #include <math.h>  
4. #include <time.h>  
5.   
6. #define N       3  
7.   
8. double fun(double x)  
9. {  
10.     double y;  
11.       
12.     y = sqrt(4-(x)*(x));  
13.     //y = sin(x);  
14.       
15.     return y;  
16. }  
17.   
18. /*随机点法求定积分*/  
19. double Darts(int n)  
20. {  
21.     double x, y;  
22.     time_t t;  
23.     int i = 0;  
24.     int count = 0;  
25.       
26.     srand((unsigned)time(&t));  
27.       
28.     for (i=0; i<n; i++)  
29.     {  
30.         x = rand()%100/100.0;  
31.         y = rand()%100/100.0;  
32.           
33.         if (y <= 1-pow(x,2))  
34.         {  
35.             count++;      
36.         }     
37.     }  
38.       
39.     return (double)count/(double)n;  
40. }  
41.   
42. /*左矩形法求定积分*/  
43. double LeftRect(double down, double up, int n)  
44. {  
45.     double h, s;  
46.     int i;  
47.       
48.     /*计算步长*/  
49.     h = (up-down)/n;  
50.       
51.     s = fun(down)*h;  
52.       
53.     for (i=1; i<n; i++)  
54.     {  
55.         s = s + fun(down+i*h)*h;      
56.     }   
57.       
58.     return s;  
59. }  
60.   
61. /*梯形公式求定积分*/  
62. double Trape(double down, double up, int n)  
63. {  
64.     double h, s;  
65.     int i = 0;  
66.       
67.     /*计算步长*/  
68.     h = (up-down)/n;  
69.       
70.     s = 0.5*(fun(down)+fun(down+h))*h;  
71.       
72.     for (i=1; i<n; i++)  
73.     {  
74.         s = s + 0.5 * (fun(down+i*h) + fun(down+(i+1)*h))*h;  
75.     }  
76.       
77.     return s;     
78. }   
79.   
80. /*复合梯形公式*/  
81. double T(double x, double y, int z)  
82. {  
83.     double h, Tn;  
84.     int i = 0;  
85.       
86.     h = (y-x)/z;  
87.     Tn = (fun(x)+fun(y))/2;  
88.       
89.     for (i=0; i<z; i++)  
90.     {  
91.         Tn = Tn+fun(x+i*h);   
92.     }     
93.       
94.     Tn = Tn*h;  
95.       
96.     return Tn;  
97. }   
98.   
99. /*辛普生公式求定积分，公式为：S[n]=(4*T[2*n]-T[n])/3，其中T[2n],T[n]为梯形公式计算结果*/  
100. double Simposn(double down, double up, int n)  
101. {  
102.     double s;  
103.       
104.     /*辛普生公式*/  
105.     s = (4*T(down, up, 2*n) - T(down, up, n))/3;  
106.       
107.     return s;     
108. }  
109.   
110. /*高斯公式求定积分*/  
111. double Gass(double (*func)(double x), double a, double b, int n)  
112. {  
113.     int i = 0;  
114.       
115.     //高斯点及其求积系数列表  
116.     float x1[1]={0.0};                                                    
117.     float A1[1]={2};  
118.     float x2[2]={-0.5573503,0.5573503};                                  
119.     float A2[2]={1,1};  
120.     float x3[3]={-0.7745967,0.0,0.7745967};                               
121.     float A3[3]={0.555556,0.888889,0.555556};  
122.     float x4[4]={0.3399810,-0.3399810,0.8611363,-0.8611363};              
123.     float A4[4]={0.6521452,0.6521452,0.3478548,0.3478548};  
124.     float x5[5]={0.0,0.5384693,-0.5384693,0.9061799,-0.9061799};         
125.     float A5[5]={0.5688889,0.4786287,0.4786287,0.2369269,0.2369269};     
126.   
127.   
128.     float *p, *t;  
129.     switch (n)  
130.     {  
131.         case 1:            
132.             p = x1;  
133.             t = A1;  
134.             break;  
135.         case 2:            
136.             p = x2;  
137.             t = A2;  
138.             break;  
139.         case 3:            
140.             p = x3;  
141.             t = A3;  
142.             break;  
143.         case 4:            
144.             p = x4;  
145.             t = A4;  
146.             break;  
147.         case 5:            
148.             p = x5;  
149.             t = A5;  
150.             break;  
151.         default :     
152.             printf("intput wrong!");  
153.     }  
154.       
155.     float g = 0;  
156.     for (i=0; i<n; i++)  
157.     {  
158.         g += (*func)((b-a)*p[i]/2+(a+b)/2)*t[i];  
159.     }   
160.     g *= (b-a)/2;  
161.       
162.     return g;  
163. }  
164.   
165. int main(int argc, char *argv[])  
166. {  
167.     printf("随机点法积分值%f\n", Darts(10000));  
168.       
169.     double down, up;  
170.     int n;  
171.     double sum = 0;  
172.       
173.     printf("积分下限:\n");  
174.     scanf("%lf", &down);  
175.     printf("积分上限:\n");  
176.     scanf("%lf", &up);  
177.     printf("分隔数目:\n");  
178.     scanf("%d", &n);  
179.       
180.     sum = LeftRect(down, up, n);  
181.     printf("左矩形法积分值为:%f\n", sum);  
182.       
183.     sum = Trape(down, up, n);  
184.     printf("梯形公式积分值为:%f\n", sum);  
185.       
186.     sum = Simposn(down, up, n);  
187.     printf("辛普生公式积分值为:%f\n", sum);  
188.       
189.     sum = Gass(fun, down, up, N);  
190.     printf("高斯公式积分值为:%f\n", sum);  
191.       
192.     return 0;  
193. }