1

Context: I have the intention to put a certain Fortran subroutine to some tests to see what is wrong with it. It's a numerical simulation and the results are not matching with theory. I use the write statement to do some simple debugging.

The problem: While a file is open in my main routine, I can't seem to write anything to the screen (so I can't check which stuff are being incorrectly passed to a certain chain of subroutines and etc.). It writes all fine when I do it before opening the file, but not inside it (after opening and after closing it).

Here is the code I'm referring to:

WRITE(*,*) 'BLA !<---------------------------------------
WRITE(*,*) 'BLA BLA' !<---------------------------------------

do p=1,N !open files
    write(posvel, "(a,i0,a)") "1Dposveldatacomelasticaxy1", p, ".dat"
    OPEN(unit=p, file=trim(posvel), status="unknown")
end do

WRITE(*,*) 'bla' !<---------------------------------------

t = tmin
cont = 0
do while ((t + dt) < (tmax))
    t = t+dt
    cont = cont+1
    do i = 1, N
        forcax(i) = 0.0d0
        forcay(i) = flagy(i)*gravidade(m(i))
        do j = 1, N
    call coefficients(m(i), m(j), gama_n, k_n)
            Fx_elastica(j,i) = 0.0d0
            Fy_elastica(j,i) = 0.0d0
    Fx_viscosa(j,i) = 0.0d0
    Fy_viscosa(j,i) = 0.0d0
    WRITE(*,*) 'inside j loop' !<---------------------------------------
            if (i .NE. j) then
                if ( (abs(sqrt(((xold(i)-xold(j))**2)+(yold(i)-yold(j))**2))).LE. (a(i)+a(j)) ) then
        WRITE(*,*) 'inside collision' !<---------------------------------------
                    call forca_elastica(k_n, a(i), a(j), xold(i), xold(j), yold(i), yold(j), Fx_elastica(j,i),&
        Fy_elastica(j,i))

                if (Fx_elastica(j,i) .GT. 0.0d0) then
                        Fx_elastica(i,j) = -Fx_elastica(j,i)
            WRITE(*,*) 'elastic x is being passed' !<---------------------------------------
                    end if
        if (Fy_elastica(j,i) .GT. 0.0d0) then
                            Fy_elastica(i,j) = -Fy_elastica(j,i)
            WRITE(*,*) 'elastic y is being passed' !<---------------------------------------
                    end if


            forcax(i) = forcax(i) + flagex(i)*Fx_elastica(j,i)
            forcay(i) = forcay(i) + flagey(i)*Fy_elastica(j,i) 
        end do
        call integracao_Euler_xy (xold(i),xnew(i),vxold(i),vxnew(i),forcax(i),yold(i),ynew(i),vyold(i),vynew(i),forcay(i),m(i))

        if (mod (cont,5000).eq. 0) then
            WRITE(p, *) int(cont/5000), t, xold(i), yold(i), forcax(i), forcay(i) !<---------------------------------------

        end if
    end do
end do

do p = 1,N !close files
    close(unit=p)
end do

Just look at the WRITE statements. The first two appear on the screen alright. After OPENing the files, though... It doesn't. The WRITE statements that depend on conditions are the ones I want to see, but Fortran is not even writing the ones that don't depend on those conditions. Also, take a look at the last WRITE statement - it writes to the file with no problems.

Any ideas on how to fix/contour this problem?

I'm using Fortran 90.

DrHAL
  • 129
  • 1
  • 8
  • Instead of the fact that you use Fortran 90 (which you probably don't use anyway, or at least the compiler does not, if it is not an old obsolete version) it is more useful to know which compiler you are using. – Vladimir F Героям слава Feb 20 '17 at 15:46

1 Answers1

3

You should not use small numbers for unit numbers. You are looping from 1 with a step of 1. You are almost guaranteed to hit the pre-connected units for standard output and standard input. See also Standard input and output units in Fortran 90?

Loop from some larger number, say from 100, or use newunit= and store the unit numbers in some array.

Also note that p has value N+1 at WRITE(p, *) ....

  • Wow, I totally forgot about this. I know there is something about the default unit number for the read and write statements to be 5 and 6, but didn't know this would cause any problems. – DrHAL Feb 20 '17 at 15:53
  • 1
    If you have a Fortran 2008 compliant compiler, you can also use the syntax `open(newunit=u(p),...)`, where `u` is an integer array. This way, you can automatically request a "safe" new unit number when opening the file. – jabirali Feb 20 '17 at 23:47