[問題]使用辛普森法做積分
我想用辛普森法做出積分且增加一容差參數。
計算N=2與N=4之結果的差值,若差值比容差大,則增加N值。
重複第二步,直至差值小於容差時,輸出積分結果與此時差值。
寫到檢驗n=4和n=2的差時就卡關了
我不知道要如何表示result(n=4)-result(n=2)這一句
我用的編譯器是Dev-C++
以下是我的程式碼
#include <stdio.h>
int main(void)
{
double a0, a1, a2, a3; // coefficients, g'(x) = a3*x^3 + a2*x^2 + a1*x + a0
double a, b; // integration interval, [a, b]
int N; // number of segments, N
double h; // length of each segments, dx
double x0, x1, x2; // 3 sample points for simpson's method, x
double f0, f1, f2; // value at 3 sample points, f(x)
int counter; // number of segment count
double result = 0; // result of integral
//get coefficient and interval from user
printf("**** Integration by simpson's method **** \n"
"Integral g'(x) dx at interval [a,b]\n\n");
printf("Please input coeficients for g'(x) = a3*x^3 + a2*x^2 + a1*x + a0
\n>");
scanf("%lf%lf%lf%lf", &a3, &a2, &a1, &a0);
printf("Please input parameter interval [a,b] \n>");
scanf("%lf%lf", &a, &b);
//get number of segment
printf("Please input the number of segments \n>");
scanf("%d", &N);
while( N > 0 )
{
//clear result for each loop
result = 0;
//check N
if( N % 2 != 0)
{
printf("N must be a even number!\nProgram if force to end.\n");
return 0;
}
//calculate length of segment
h = (b - a) / N;
//initialize x to first segment
x0 = a;
x1 = a + h;
x2 = a + 2*h;
//calculation by simpson's method for each segment
for( counter = 1 ; counter <= N/2 ; counter++ )
{
//calculate value at 3 sample points
f0 = a3*x0*x0*x0 + a2*x0*x0 + a1*x0 + a0;
f1 = a3*x1*x1*x1 + a2*x1*x1 + a1*x1 + a0;
f2 = a3*x2*x2*x2 + a2*x2*x2 + a1*x2 + a0;
//add to result by simpson's rule
result = result + h/3 * (f0 + 4.0*f1 + f2);
//move the 3 sample point to next segment
x0 = x0 + 2*h;
x1 = x1 + 2*h;
x2 = x2 + 2*h;
}//end for
//print out result
printf("The result of the integration is %lf\n", result);
//check for re-loop
printf("\nDo you want to rechoose the number of segments?(enter -1 to
leave)\n>");
scanf("%d", &N );
}//end while
return 0;
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 140.114.201.34
→
03/23 15:03, , 1F
03/23 15:03, 1F
→
03/24 18:39, , 2F
03/24 18:39, 2F
討論串 (同標題文章)