FORTRAN: Loop over variable file names

Hi guys

I'm a beginner in fortran. So excuse me for my naivety, let me briefly describe what I was trying to do. I have let's say 2 files named reac-1 and reac-2. After opening these files I've to do some calculations, close these files and open the same files again in a loop. So my faulty code looks like this:

do j= 1,29  !temp loop

write(*,*)"This is for temperature: ",T(j)
    
  do i=1,r	! file loop
  write(filename,100)i
  100 format('reac-',I1)											

	write(*,*)"This is for reac-",i
   
   //open file, do something, close file.....//	

  end do	! end file loop

   //do something...//

end do  ! end temp loop

I get a run time error when it enters temperature loop the 2nd time, it stops after printing "This is for reac-1". My guess is it already had reac-1 & reac-2 in 'filename' variable, so it can't write reac-1 again and that's why it crashes.

Could you guys shed some light on it where I'm going wrong here and how can I overcome it?

Thanks a lot!

Saleheen

Where is T() declared? Is it an array or a function? If it is an array, where are the elements of T set?

Where is r defined? Where is r set?

If it stops after printing reac-1 without printing reac-2 , what makes you think it is dying in the "temp loop" instead of dying in the "file loop"?

Hi Don

Thanks a lot for you reply. 'T' is an array, I've defined the elements initially. 'r' is number of input files of that species (here 'reac') I have. It prompts me for how many species of reactant I have. It doesn't stop without printing 'reac-2', what it does is for a single temperature it prints for reac-1, then reac-2. And then when it goes for a new temperature, it prints

This is for temperature: 298

, then

This is for reac- 1

and than it stops. That's why I was thinking it's dying in the "file loop". I'm posting a part of the code in case you want to have a look at it.

program check
implicit none


! Setting up parameters and necessary conversion factors
real,parameter::h=4.135667516e-15;            		    ! planck's constant in eVs
real,parameter::c=3.0e10;                           	            ! velocity of light in cm/s
real,parameter::kb=8.6173324e-5;                            ! boltzmann's constant in eV/K

							
integer, parameter:: dp = selected_real_kind(15, 307)
real(dp),allocatable::zper(:)					! zero point energy array 
real(dp),allocatable::zpcr(:)				! zero point corrected energy
real(dp),allocatable::qvibr(:)
real(dp)::zpcrtot=0.0					! total zero point corrected energy 

integer::i,j
integer::r
character(len=20)::filename	
integer::ierror,alloc_err,nlines,nless,array_t
real(dp)::scf,vibmulti,freqtot
real(dp),allocatable::freq(:),vib(:)
real(dp)::T(29) =(/(array_t,array_t=273,973,25)/)



write(*,*)"How many reactants do you have?"
read(*,*)r


!allocate proper dimension to all zero point arrays
allocate(zper(r))
allocate(zpcr(r))
allocate(qvibr(r))

! Loop through reactant files

do j= 1,29 !temp loop

write(*,*)"This is for temperature: ",T(j)
    
  do i=1,r	! file loop
  write(filename,100)i
  100 format('reac-',I1)											! get filenames

	write(*,*)"This is for reac-",i
     
	freqtot=0														! initializing total frequency value
	vibmulti=1														!
    ! Count no of lines	
	nlines=0
	OPEN(unit=10,file=filename,status='old',action='read',&
	 iostat=ierror)
			do	! count loop
            read(10,*,iostat=ierror)
             if (ierror /= 0) exit
  			nlines=nlines+1
            end do !end count loop
	CLOSE(unit=10)
	
	! Count no of frequency elements
	nless=nlines-1
	
	!allocate a dimension to frquency array
	allocate(freq(nless))
	allocate(vib(nless))
	
	OPEN(unit=10,file=filename,status='old',action='read',&
	 iostat=ierror)

			read(10,*,iostat=ierror) scf
			
			do  	! freq loop
            read(10,*,iostat=ierror) freq(nless)
            if (ierror /= 0) exit
          	if (freq(nless) < 100) then
              	freq(nless) = 100
            end if            
            
			freqtot=freqtot+freq(nless)
            			           
			vib(nless)= 1/(1-exp(-h*c*freq(nless)/kb/T(j)))
            
            vibmulti=vib(nless)*vibmulti
			
	     	end do	!end freq loop
     
     deallocate (freq,stat=alloc_err) 
     deallocate (vib,stat=alloc_err)      		
     CLOSE(10)


zper(i)= 0.5*h*c*freqtot												! zero point energy
zpcr(i)= scf + zper(i)													! zero point correction
zpcrtot= zpcrtot+zpcr(i)												! total sum of zero point corrected energies of reactants												

write(*,*)"ZPER", zper(i)
write(*,*)"ZPCR", zpcr(i)
write(*,*)"ZPCRTOT", zpcrtot
write(*,*)"Vib.P.F for reac",i,vibmulti  
qvibr(i)=vibmulti
write(*,*)"This is test qvib: ", qvibr(i) 

	end do	! end file loop
    
write(*,*) "Try product: ",PRODUCT(qvibr)

deallocate (qvibr,stat=alloc_err)     
deallocate (zper,stat=alloc_err) 
deallocate (zpcr,stat=alloc_err)

end do  ! end temp loop

end program check

Again thanks a lot!

I haven't written much FORTRAN since ~1975, the language is a lot more free-format than the language I used at that time, and I find your indentation hard to follow; but if I am reading your code correctly, you are allocating space for the zper() , zpcr() , and qvibr() arrays before the start of the temp loop and deallocating space for those arrays inside the temp loop. So, my guess would be that your program is dying on the line:

zper(i)= 0.5*h*c*freqtot								! zero point energy

the 2nd time through your temp loop because you are referencing an unallocated array.

If I were writing this code, I think I would also read the reac-X files once shortly after reading r and save the file sizes in an array so I wouldn't need to recalculate these values each time through your temp loop.

Awesome! Thank you so much Don!

Hi.

To make one's code more readable, one can use a tidy or pretty-printing code. Here's an example:

#!/usr/bin/env bash

# @(#) s1	Demonstrate Fortran-90 preprocessor, pretty-print, f90pp4.
# For pre-processor and pretty-printer Fortran-90, see:
# http://www.ifremer.fr/ditigo/molagnon/fortran90/

pe() { for _i;do printf "%s" "$_i";done; printf "\n"; }
pl() { pe;pe "-----" ;pe "$*"; }
db() { ( printf " db, ";for _i;do printf "%s" "$_i";done;printf "\n" ) >&2 ; }
db() { : ; }
C=$HOME/bin/context && [ -f $C ] && $C f90ppr

FILE=${1-check.f90}

pl " Input code, $FILE:"
head $FILE

pl " Results from f90ppr:"
rm -f new.txt
f90ppr < $FILE > new.txt
head new.txt

exit 0

producing:

$ ./s1

Environment: LC_ALL = , LANG = en_US.UTF-8
(Versions displayed with local utility "version")
OS, ker|rel, machine: Linux, 2.6.26-2-amd64, x86_64
Distribution        : Debian 5.0.8 (lenny, workstation) 
bash GNU bash 3.2.39
f90ppr - ( local: ~/executable/f90ppr, 2011-01-30 )

-----
 Input code, check.f90:
$define FPPR_STP_INDENT 2
program check
implicit none


! Setting up parameters and necessary conversion factors
real,parameter::h=4.135667516e-15;            		    ! planck's constant in eVs
real,parameter::c=3.0e10;                           	            ! velocity of light in cm/s
real,parameter::kb=8.6173324e-5;                            ! boltzmann's constant in eV/K


-----
 Results from f90ppr:
This is f90ppr: @(#) fppridnt.f90	V-1.3 00/05/09 Michel Olagnon
( usage: f90ppr < file.F90  > file.f90 )
Program check
    Implicit None
!
!
! Setting up parameters and necessary conversion factors
    Real, Parameter :: h = 4.135667516e-15
! planck's constant in eVs
    Real, Parameter :: c = 3.0e10
! velocity of light in cm/s
    Real, Parameter :: kb = 8.6173324e-5

There are a number of settings ( like $define FPPR_STP_INDENT 2 ) for modifying the output. There is also a unix-like code in C to allow settings from the command line.

See attached file for a completed tidying.

Best wishes ... cheers, drl

Awesome man! You guys are the best! Thanks a lot!