A strange fortran code


by phynewb
Tags: code, fortran, strange
phynewb
phynewb is offline
#1
Oct20-10, 02:38 PM
P: 13
Hi guys,
I just wrote a simple code for delta function and print it out. Here it is.

      program main
        implicit none
        integer i
        real*8 del
        real*8,parameter:: ep=1D-2
        real*8,parameter:: pi=3.1415926
        do i=-10,10,1
        write(*,*) i,del(real(i))  ! use function
        enddo
        do i=-10,10,1
        write(*,*) i,(1.0/pi*ep/((real(i))**2+ep**2)) ! no use function
        enddo
        end

      function del(x)
      implicit none
      real*8 del
      real*8 x
      real*8,parameter::ep=1D-2
      real*8,parameter:: pi=3.1415926

      del=(1.0/pi*ep/(x**2+ep**2))

      return
      end
It turns out to be quite wrong by using function del. I cannot figure out what's wrong with the function because the expression is the same. Can you tell me how to debug the function? Thanks!
Phys.Org News Partner Science news on Phys.org
Better thermal-imaging lens from waste sulfur
Hackathon team's GoogolPlex gives Siri extra powers
Bright points in Sun's atmosphere mark patterns deep in its interior
Mark44
Mark44 is offline
#2
Oct21-10, 10:11 AM
Mentor
P: 20,962
Quote Quote by phynewb View Post
Hi guys,
I just wrote a simple code for delta function and print it out. Here it is.

      program main
        implicit none
        integer i
        real*8 del
        real*8,parameter:: ep=1D-2
        real*8,parameter:: pi=3.1415926
        do i=-10,10,1
        write(*,*) i,del(real(i))  ! use function
        enddo
        do i=-10,10,1
        write(*,*) i,(1.0/pi*ep/((real(i))**2+ep**2)) ! no use function
        enddo
        end

      function del(x)
      implicit none
      real*8 del
      real*8 x
      real*8,parameter::ep=1D-2
      real*8,parameter:: pi=3.1415926

      del=(1.0/pi*ep/(x**2+ep**2))

      return
      end
It turns out to be quite wrong by using function del. I cannot figure out what's wrong with the function because the expression is the same. Can you tell me how to debug the function? Thanks!
What output do your two loops produce?
phynewb
phynewb is offline
#3
Oct21-10, 12:17 PM
P: 13
Quote Quote by Mark44 View Post
What output do your two loops produce?
It shows
-10 31.830990148286446
-9 4.42379935910539249E-054
-8 3.96072696157977595E-306
-7 31.830990148286446
-6 4.42379936653263965E-054
-5 0.0000000000000000
-4 31.830990148286446
-3 4.42379937841623373E-054
-2 0.0000000000000000
-1 31.830990148286446
0 4.42380394766179441E-054
1 31.830990148286446
2 4.42380242655975552E-054
3 31.830990148286446
4 4.42380241467614869E-054
5 31.830990148286446
6 4.42380240873434528E-054
7 31.830990148286446
8 4.42380240279254245E-054
9 31.830990148286446
10 4.42380239982164045E-054

-10 3.18309583173281331E-005
-9 3.92974701861929112E-005
-8 4.97358443944407051E-005
-7 6.49610718106421125E-005
-6 8.84191714697638248E-005
-5 1.27323451299340604E-004
-4 1.98942445036508832E-004
-3 3.53673738606087140E-004
-2 7.95754859835665312E-004
-1 3.18278073675496917E-003
0 31.830990148286446
1 3.18278073675496917E-003
2 7.95754859835665312E-004
3 3.53673738606087140E-004
4 1.98942445036508832E-004
5 1.27323451299340604E-004
6 8.84191714697638248E-005
7 6.49610718106421125E-005
8 4.97358443944407051E-005
9 3.92974701861929112E-005
10 3.18309583173281331E-005


The upper part is from the function and the lower part is calculated without function. It shall look like sharp gaussian distribution around 0 instead of random distribution.

Mark44
Mark44 is offline
#4
Oct21-10, 03:56 PM
Mentor
P: 20,962

A strange fortran code


I'm pretty sure this line is causing the problem:
write(*,*) i,del(real(i))  ! use function
In the call to your del function, you are converting an integer value, i, to a four-byte real (a single-precision floating point value), but the function is expecting an eight-byte real (a double-precision floating point value).

I think this will fix the problem:
write(*,*) i,del(dble(i))  ! use function
Here's a link to a page with more information - http://www.liv.ac.uk/HPC/HTMLF90Cour...esnode114.html.
phynewb
phynewb is offline
#5
Oct21-10, 08:55 PM
P: 13
Quote Quote by Mark44 View Post
I'm pretty sure this line is causing the problem:
write(*,*) i,del(real(i))  ! use function
In the call to your del function, you are converting an integer value, i, to a four-byte real (a single-precision floating point value), but the function is expecting an eight-byte real (a double-precision floating point value).

I think this will fix the problem:
write(*,*) i,del(dble(i))  ! use function
Here's a link to a page with more information - http://www.liv.ac.uk/HPC/HTMLF90Cour...esnode114.html.
Thanks Mark44. It really works! I am always troubled by this accuracy stuff...
Thank you a lot :p
Mark44
Mark44 is offline
#6
Oct22-10, 12:45 AM
Mentor
P: 20,962
It's more about the size, in bytes, of the variables you're using. For floating point variables, the choices are real*4 (single precision - 4 bytes) and real*8 (double precision - 8 bytes). real() converts whatever to a 4-byte floating point value and dble() converts whatever to an 8-byte floating point value.

Make sure that what you pass to a function is what it is expecting - otherwise you'll have problems.


Register to reply

Related Discussions
Help me write Fortran code Programming & Computer Science 2
Looking for fortran code or subroutine ,please help!!! Programming & Computer Science 0
Problem with Fortran g77 code Programming & Computer Science 4
Code::Blocks and Fortran Programming & Computer Science 1
Convert latex code to fortran code? Math & Science Software 1