Write a program that will numerically integrate a given function at a given interval using the trapezoidal rule method

[emergency9177]

Plase can someone help me with this computational problem:

Write a program that will numerically integrate a given function at a given interval using the trapezoidal rule method. Try to write the algorithm in such a way that Integral achieves the desired accuracy or that you evaluate the accuracy of the integration (eg. perform two calculations with different steps). Solve the problem in C programming language. Submit the text file with the program containing the solution.

f(x) = exp(-x2/10)•sin(3x) / 5x

Calculate integral of f(x) in interval between 0 and +∞.

I'v try to solve the problem:

1 #include<conio.h>
2 #include<stdio.h>
3 void main()
4 {
5 float a,b,h,x,y,s=0;
6 float x0,y0,xn,yn,r;
7 int i,n;
8 float f(float);
9 printf("\nEnter lower limit value:: ");
10 scanf("%f",&a);
11 printf("\nEnter upper limit value:: ");
12 scanf("%f",&b);
13 printf("\nEnter the number of intervals:: ");
14 scanf("%d",&n);
15 h=(b-a)/n;
16 x0=a;
17 y0=f(a);
18 yn=f(b);
19 x=x0+h;
20 for(i=1;i<=n-1;i++)
21 {
22 y=f(x);
23 s=s+y;
24 x=x+h;
25 }
26 r=(h/2)*(y0+yn+2*s);
27 printf("\nResult is:: %f",r);
28 getch();
29 }
30 float f(float x)
31 {
32 return exp(-x*x/10.0)*sin(3.0*x) / (5.0*x);
33 }

After running the program I receive this:

main.c:32:9: warning: implicit declaration of function ‘exp’ [-Wimplicit-function-declaration]
32 | return exp(-x*x/10.0)*sin(3.0*x) / (5.0*x);
| ^~~
main.c:32:9: warning: incompatible implicit declaration of built-in function ‘exp’
main.c:3:1: note: include ‘<math.h>’ or provide a declaration of ‘exp’
2 | #include<stdio.h>
+++ |+#include <math.h>
3 | void main()
main.c:32:24: warning: implicit declaration of function ‘sin’ [-Wimplicit-function-declaration]
32 | return exp(-x*x/10.0)*sin(3.0*x) / (5.0*x);
| ^~~
main.c:32:24: warning: incompatible implicit declaration of built-in function ‘sin’
main.c:32:24: note: include ‘<math.h>’ or provide a declaration of ‘sin’

Enter lower limit value:: 0

Enter upper limit value:: +∞

Enter the number of intervals::
Result is:: -nan

...Program finished with exit code 0
Press ENTER to exit console.

_____________________________________________________________________________________________________________________________________________________________________

I've already asked my prof. if the expression f(x) = exp(-x2/10)•sin(3x) / 5x is correct, he told me yes, and I also aksed him if I have to use the Gauss-Laguerre quadrature since I have to Calculate integral of f(x) in interval between 0 and +∞.

He answered me, that we defined the value of a function in lectures when x = 0. We looked at the function limit when x goes to 0. And also that we have to program this as follows: when calculating the function we check if x is equal to 0. In this case we return the value of the limit, otherwise we calculate the value.​

I'm stuck, please some help, any hint?

Thank you very much.

 

[lev888]

1. Can you describe your algorithm?
2. Seems like you need to include math.h.

 

[emergency9177]

1. Can you describe your algorithm?
2. Seems like you need to include math.h.

Thank you.

1. The algorithm use defined function to calculate the value of function f(x) =exp(-x2/10)•sin(3x) / 5x. The variable data type used in the program are integer and float types.

2. I include math.h and when I run the program it comes:

Enter lower limit value:: 0

Enter upper limit value:: +∞

Enter the number of intervals::
Result is:: -nan

...Program finished with exit code 0
Press ENTER to exit console.

Still not working.

Please help me.

Thank you.

 

[lev888]

Thank you.

1. The algorithm use defined function to calculate the value of function f(x) =exp(-x2/10)•sin(3x) / 5x. The variable data type used in the program are integer and float types.

2. I include math.h and when I run the program it comes:

Enter lower limit value:: 0

Enter upper limit value:: +∞

Enter the number of intervals::
Result is:: -nan

...Program finished with exit code 0
Press ENTER to exit console.

Still not working.

Please help me.

Thank you.

This is not a description of the algorithm. What are the actual steps? Do you understand how numerical integration algorithms work?

 

[Cubist]

Enter upper limit value:: +∞

Attempting to set variable "b" to infinity won't work. What value would you expect "h" to contain after this calculation (on line 15)...

h=(b-a)/n;

Basically, you'll have to choose some finite value for b when performing a numerical integration. Perhaps a sketch or graph of the function will help you to choose a suitable value.

 

[emergency9177]

Attempting to set variable "b" to infinity won't work. What value would you expect "h" to contain after this calculation (on line 15)...

h=(b-a)/n;

Basically, you'll have to choose some finite value for b when performing a numerical integration. Perhaps a sketch or graph of the function will help you to choose a suitable value.

[math]\int_{0}^{+\infty}\dfrac{\mathrm{e}^{-\frac{x^2}{10}}\sin\left(3x\right)}{5x}[/math]
The result of the above integral is: 0.3141592

The graph of the functions is attached (three different zooms).







I've changed te algorithm:

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
double f (double x) { double myresult = exp(-x*x/10.0)*sin(3.0*x) / (5.0*x); return myresult; }
int main()
{
 float lower, upper, integration=0.0, stepSize, k;
 int i, subInterval;
system("clear");
 printf("Enter lower limit of integration: ");
 scanf("%f", &lower);
 printf("Enter upper limit of integration: ");
 scanf("%f", &upper);
 printf("Enter number of sub intervals: ");
 scanf("%d", &subInterval);
 stepSize = (upper - lower)/subInterval;
 integration = f(lower) + f(upper);
 for(i=1; i<= subInterval-1; i++)
 {
  k = lower + i*stepSize;
  integration = integration + 2 * f(k);
 }
 integration = integration * stepSize/2;
 printf("\nIntegral of given function for given values is: %.3f", integration);
 getchar();
 return 0;
}

I am receiving again the same answers:

Enter lower limit of integration: 0
Enter upper limit of integration: +∞
Enter number of sub intervals:
Integral of given function for given values is: -nan

...Program finished with exit code 0
Press ENTER to exit console.

__________________________________________________________________________________________________________________________________________________

The program seems do not like the fact that I have the upper limit of integration: +∞.

Any help how to code correctly the program?

Thanks!

 

[blamocur]

[math]\int_{0}^{+\infty}\dfrac{\mathrm{e}^{-\frac{x^2}{10}}\sin\left(3x\right)}{5x}[/math]
The result of the above integral is: 0.3141592

The graph of the functions is attached (three different zooms).

View attachment 30497

View attachment 30498

View attachment 30499

I've changed te algorithm:

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
double f (double x) { double myresult = exp(-x*x/10.0)*sin(3.0*x) / (5.0*x); return myresult; }
int main()
{
 float lower, upper, integration=0.0, stepSize, k;
 int i, subInterval;
system("clear");
 printf("Enter lower limit of integration: ");
 scanf("%f", &lower);
 printf("Enter upper limit of integration: ");
 scanf("%f", &upper);
 printf("Enter number of sub intervals: ");
 scanf("%d", &subInterval);
 stepSize = (upper - lower)/subInterval;
 integration = f(lower) + f(upper);
 for(i=1; i<= subInterval-1; i++)
 {
  k = lower + i*stepSize;
  integration = integration + 2 * f(k);
 }
 integration = integration * stepSize/2;
 printf("\nIntegral of given function for given values is: %.3f", integration);
 getchar();
 return 0;
}

I am receiving again the same answers:

Enter lower limit of integration: 0
Enter upper limit of integration: +∞
Enter number of sub intervals:
Integral of given function for given values is: -nan

...Program finished with exit code 0
Press ENTER to exit console.

__________________________________________________________________________________________________________________________________________________

The program seems do not like the fact that I have the upper limit of integration: +∞.

Any help how to code correctly the program?

Thanks!

Reread the post by @Cubist which says that you cannot use infinity; instead use a large enough upper limit that will give good approximation. Try varying that limit and see if the result changes significantly.
Also, your definition of function 'f' won't work for x = 0; you'll need to make a special case there.

 

[lev888]

Try varying that limit

I think you meant try varying the number of intervals.
Edit: I guess it makes sense to play with both parameters.

 

[blamocur]

I think you meant try varying the number of intervals.

Actually both. Varying the upper limit will show whether the upper limit is sufficiently large. I.e., if increasing the upper limit doesn't change the result then the limit is a good enough approximation of infinity. To keep it clean I would adjust the number of intervals to keep the step the same when varying the upper limit.

 

[lev888]

Actually both. Varying the upper limit will show whether the upper limit is sufficiently large. I.e., if increasing the upper limit doesn't change the result then the limit is a good enough approximation of infinity. To keep it clean I would adjust the number of intervals to keep the step the same when varying the upper limit.

Already edited my reply. I finally noticed that the upper limit is infinity.

 

[emergency9177]

Actually both. Varying the upper limit will show whether the upper limit is sufficiently large. I.e., if increasing the upper limit doesn't change the result then the limit is a good enough approximation of infinity. To keep it clean I would adjust the number of intervals to keep the step the same when varying the upper limit.

I've tried different combination varying the lower and upper limits as the number of steps. The program gives results for all type of boundaries except for the lower limit 0 and the upper limit ∞. When gradually varying the combinations, the program gives minimal differences.

What now? I am asked to calculate integral of f(x) in interval between 0 and +∞ ? Not from ex. 1E-9 to 1E9 !

Thanks

 

[lev888]

I've tried different combination varying the lower and upper limits as the number of steps. The program gives results for all type of boundaries except for the lower limit 0 and the upper limit ∞. When gradually varying the combinations, the program gives minimal differences.

What now? I am asked to calculate integral of f(x) in interval between 0 and +∞ ? Not from ex. 1E-9 to 1E9 !

Thanks

How do you propose to handle ∞ as a numeric input???

 

[Cubist]

Why not try to get an idea of the error introduced by NOT going all the way to +infinity

[math]\frac{\mathrm{d}}{\mathrm{d}x}\left(-\mathrm{e}^{\frac{-x^2}{10}}\right)[/math][math]= \frac{x \mathrm{e}^{\frac{-x^2}{10}}}{ 5}[/math]
The above will always be bigger than your f(x) for x>1. What is the integral of the above from, say, 100 to infinity?

 

[emergency9177]

Why not try to get an idea of the error introduced by NOT going all the way to +infinity

[math]\frac{\mathrm{d}}{\mathrm{d}x}\left(-\mathrm{e}^{\frac{-x^2}{10}}\right)[/math][math]= \frac{x \mathrm{e}^{\frac{-x^2}{10}}}{ 5}[/math]
The above will always be bigger than your f(x) for x>1. What is the integral of the above from, say, 100 to infinity?

Thanks for hint!

The numerical integration of

[math]\dfrac{x\mathrm{e}^{-\frac{x^2}{10}}}{5}[/math]
from 0 to infinity gives: 1

from 100 to infinity gives: 0

I slightly modified the code in C program:

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
double f (double x) { double myresult = (x*exp(-x*x/10))/5; return myresult; }
int main()
{
 float lower, upper, integration=0.0, stepSize, k;
 int i, subInterval;
system("clear");
 printf("Enter lower limit of integration: ");
 scanf("%f", &lower);
 printf("Enter upper limit of integration: ");
 scanf("%f", &upper);
 printf("Enter number of sub intervals: ");
 scanf("%d", &subInterval);
 stepSize = (upper - lower)/subInterval;
 integration = f(lower) + f(upper);
 for(i=1; i<= subInterval-1; i++)
 {
  k = lower + i*stepSize;
  integration = integration + 2 * f(k);
 }
 integration = integration * stepSize/2;
 printf("\nIntegral of given function for given values is: %.3f", integration);
 getchar();
 return 0;
}

But when I run the program I obtain (as always):

Enter lower limit of integration: 0
Enter upper limit of integration: ∞
Enter number of sub intervals:
Integral of given function for given values is: -nan

...Program finished with exit code 0
Press ENTER to exit console.

Or:

Enter lower limit of integration: 100
Enter upper limit of integration: ∞
Enter number of sub intervals:
Integral of given function for given values is: -nan

...Program finished with exit code 0
Press ENTER to exit console.

________________________________________________________________________________________________________________________________________________

At this point I think there is something wrong with the codes, and not with the conceptual approach to the problem .

 

[Cubist]

Thanks for hint!

The numerical integration of
[math]\dfrac{x\mathrm{e}^{-\frac{x^2}{10}}}{5}[/math]from 100 to infinity gives: 0

Actually, instead of "100" let's use "15". The integration (of the expression from post#13) between 15 and +inf gives exactly [imath]e^{-22.5}[/imath] which approximates to 0.000000000169. This is the maximum amount of error that would be introduced by performing a numerical integration from 0 to 15 instead of 0 to +inf. Is this an acceptable amount of error to introduce?

So, instead of typing in, "Enter upper limit of integration: ∞"
Type this, "Enter upper limit of integration: 15"
...and accept that the output might be different by up to ≈ 0.000000000169

The errors introduced by the numerical integration itself will probably be much more significant than the above.

 

[emergency9177]

Actually, instead of "100" let's use "15". The integration (of the expression from post#13) between 15 and +inf gives exactly [imath]e^{-22.5}[/imath] which approximates to 0.000000000169. This is the maximum amount of error that would be introduced by performing a numerical integration from 0 to 15 instead of 0 to +inf. Is this an acceptable amount of error to introduce?

So, instead of typing in, "Enter upper limit of integration: ∞"
Type this, "Enter upper limit of integration: 15"
...and accept that the output might be different by up to ≈ 0.000000000169

The errors introduced by the numerical integration itself will probably be much more significant than the above.

Accordingly to the problem "write the algorithm in such a way that Integral achieves the desired accuracy or that you evaluate the accuracy of the integration (eg. perform two calculations with different steps)." it should be and acceptable amount of error.

I did not understand if I have also to insert one more code step into algorithm " accept that the output might be different by up to ≈ 0.000000000169 " , or it was just a comment from you to alert me, that this is the maximum amount of error by performing a numerical integration from 0 to 15?

Thanks!

 

[blamocur]

that this is the maximum amount of error by performing a numerical integration from 0 to 15?

I think you got it right. I.e., if you integrate from 0 to 15 the resulting error will be less then the one @Cubist computed for a limiting function. And if the error of [imath]1.69\times 10^{-10}[/imath] is not considered acceptable you can recompute the integration limit by using @Cubist's technique to get a different error limit.

 

[emergency9177]

I think you got it right. I.e., if you integrate from 0 to 15 the resulting error will be less then the one @Cubist computed for a limiting function. And if the error of [imath]1.69\times 10^{-10}[/imath] is not considered acceptable you can recompute the integration limit by using @Cubist's technique to get a different error limit.

I got it.

Thank you!