FloatingFeather
Programmer
Hi. I have this doubt. I have written a code that integrates by means of the composite Simpsons rule:
As can be seen in the wikipedia article, the weights are distributed according to the point where we are evaluating the function. So, for the first and last terms we have a certain weight, let's say w1,
then for the even points we have a weight w2, and for the odd points of evaluation we have a third weight w3. Suppose we want to integrate f(x)=x^2, between xmin and xmax.
One easy way of coding this in fortran would be like this:
The problem with this approach is that, if we are making heavy use of this, the recurrent use of if's inside the loops can become computationally expensive. A more efficient way of coding this would be:
My question is if I should actually write the code in this last form (which is more laborious) in order to avoid the recursive if statements, or if the compiler will automatically realize of this, and automatically unwrap the first do loop with the if's into something like the last code.
Thanks in advance.
As can be seen in the wikipedia article, the weights are distributed according to the point where we are evaluating the function. So, for the first and last terms we have a certain weight, let's say w1,
then for the even points we have a weight w2, and for the odd points of evaluation we have a third weight w3. Suppose we want to integrate f(x)=x^2, between xmin and xmax.
One easy way of coding this in fortran would be like this:
Code:
[indent]
dx=(xmax-xmin)/N
rint=0.d0
do i=1,N+1
x=xmin+(i-1)*dx
a=4.d0/3.d0 !this one would be w2
if(mod(i,2).ne.0) a=2.d0/3.d0 !this one is w3
if(i.eq.1 .or. i.eq.N+1) a=1.d0/3.d0 !this one is w1
rint=x**2*a+rint
enddo
rint=dx*rint
[/indent]The problem with this approach is that, if we are making heavy use of this, the recurrent use of if's inside the loops can become computationally expensive. A more efficient way of coding this would be:
Code:
[indent]
dx=(xmax-xmin)/N
rint=0.d0
x=xmin
a=1.d0/3.d0 !this one is w1
rint=x**2*a
x=xmax
rint=x**2*a+rint
do i=2,N,2
x=xmin+(i-1)*dx
a=4.d0/3.d0 !this one would be w2
rint=x**2*a+rint
enddo
do i=3,N-1,2
x=xmin+(i-1)*dx
a=2.d0/3.d0 !this one would be w3
rint=x**2*a+rint
enddo
rint=dx*rint
[/indent]My question is if I should actually write the code in this last form (which is more laborious) in order to avoid the recursive if statements, or if the compiler will automatically realize of this, and automatically unwrap the first do loop with the if's into something like the last code.
Thanks in advance.
