Numerical integration: Difference between revisions
Content added Content deleted
(C) |
|||
Line 256:
simpson = h / 6 * (f(a) + f(b) + 4 * sum1 + 2 * sum2)
END FUNCTION
=={{header|C}}==
{{trans|Java}}
<c>#include <math.h>
double int_leftrect(double from, double to, double n, double (*func)())
{
double h = (to-from)/n;
double sum = 0, x;
for(x=from; x <= (to-h); x += h)
sum += func(x);
return h*sum;
}
double int_rightrect(double from, double to, double n, double (*func)())
{
double h = (to-from)/n;
double sum = 0, x;
for(x=from+h; x <= (to-h); x += h)
sum += func(x);
return h*sum;
}
double int_midrect(double from, double to, double n, double (*func)())
{
double h = (to-from)/n;
double sum = 0.0, x;
for(x=from; x <= (to-h); x += h)
sum += (func(x) + func(x + h));
return h*sum/2.0;
}
double int_trapezium(double from, double to, double n, double (*func)())
{
double h = (to - from) / n;
double sum = func(from) + func(to);
int i;
for(i = 1;i < n;i++)
sum += func(from + i * h);
return h * sum;
}
double int_simpson(double from, double to, double n, double (*func)())
{
double h = (to - from) / n;
double sum1 = 0.0;
double sum2 = 0.0;
int i;
for(i = 0;i < n;i++)
sum1 += func(from + h * i + h / 2.0);
for(i = 1;i < n;i++)
sum2 += func(from + h * i);
return h / 6.0 * (func(from) + func(to) + 4.0 * sum1 + 2.0 * sum2);
}</c>
Usage example:
<c>#include <stdio.h>
#include <stdlib.h>
double f1(double x)
{
return 2.0*x*log(x*x+1.0);
}
double f1a(double x)
{
double y = x*x+1;
return y*(log(y)-1.0);
}
double f2(double x)
{
return cos(x);
}
double f2a(double x)
{
return sin(x);
}
typedef double (*pfunc)(double, double, double, double (*)());
typedef double (*rfunc)(double);
#define INTG(F,A,B) (F((B))-F((A)))
int main()
{
int i,j;
double ic;
pfunc f[5] = { &int_leftrect, &int_rightrect,
&int_midrect, &int_trapezium,
&int_simpson };
const char *names[5] = {"leftrect", "rightrect", "midrect",
"trapezium", "simpson" };
rfunc rf[2] = { &f1, &f2 };
rfunc If[2] = { &f1a, &f2a };
for(j=0; j<2; j++)
{
for(i=0; i < 5 ; i++)
{
ic = (*f[i])(0.0, 1.0, 100.0, rf[j]);
printf("%10s [ 0,1] num: %+lf, an: %lf\n",
names[i], ic, INTG((*If[j]), 0.0, 1.0));
}
printf("\n");
}
}</c>
Output of the test (with a little bit of enhancing by hand):
<pre>
leftrect [ 0,1] num: +0.365865, an: 0.386294
rightrect [ 0,1] num: +0.365865
midrect [ 0,1] num: +0.372628
trapezium [ 0,1] num: +0.393254
simpson [ 0,1] num: +0.386294
leftrect [ 0,1] num: +0.838276, an: 0.841471
rightrect [ 0,1] num: +0.828276
midrect [ 0,1] num: +0.836019
trapezium [ 0,1] num: +0.849165
simpson [ 0,1] num: +0.841471
</pre>
=={{header|C++}}==
| |||