a fortran environment

[danbaron]

Here is a shorter version of the previous program.

I don't think Simpson's Rule is the fastest way to numerically integrate, but, you can get pretty good accuracy if you divide the interval into tiny sub-intervals.

Below the interval is divided into ten million equal sub-intervals, and the result is accurate to 1 part in 3,141,592,653,589 - approximately one part in three trillion.

' code ----------------------------------------------------------------------------------------------------------------

module integrate

contains

double precision function xsquaredcosxdiv4(x)
double precision::x
xsquaredcosxdiv4=x*x*cos(x)/4
end function

double precision function simpsonintegrate(x1,x2,ndivisions,f)
double precision::x1,x2,f
integer::ndivisions
double precision::division,sum
integer::i
division=(x2-x1)/ndivisions
sum=(f(x1)+f(x2))/2
do i= 1,ndivisions-1
sum=sum+f(x1+division*i)
enddo
simpsonintegrate=division*sum
end function

end module

!-----------------------------------------------------------------------------------------------------------

program test
use integrate
implicit none
double precision::pi
pi=4*atan(1d0)
print*,"pi               = ",pi
print*,"integral approx. = ",simpsonintegrate(0d0,2*pi,10000000,xsquaredcosxdiv4)
end program

' output --------------------------------------------------------------------------------------------------------------

 pi               =    3.1415926535897931
 integral approx. =    3.1415926535892416

[danbaron]

I'm at the point now in my "top secret" project, where I needed to find out if Fortran would do what I want with respect to recursion.

So, I made a test.

Amazingly, things worked the way I hoped.

' code -------------------------------------------------------------------------------------------------------------------------------------

recursive subroutine test(i)
integer::i
call dog(i)
if(mod(i,100)==0)print*,i
endsubroutine

!--------------------------------------------------------------

subroutine dog(i)
integer::i
if(i<4000) call test(i+1)
endsubroutine

!--------------------------------------------------------------

program x15
implicit none

call test(0)

end program

' output -----------------------------------------------------------------------------------------------------------------------------------

        4000
        3900
        3800
        3700
        3600
        3500
        3400
        3300
        3200
        3100
        3000
        2900
        2800
        2700
        2600
        2500
        2400
        2300
        2200
        2100
        2000
        1900
        1800
        1700
        1600
        1500
        1400
        1300
        1200
        1100
        1000
         900
         800
         700
         600
         500
         400
         300
         200
         100
           0