NaN output outside but not inside DO loop

[Vahid9]

Dear Fortran Users,

In the following code, I can print all the values for conc_h inside the first DO loop over i before ENDDO and the total conc_h right after the ENDDO.
When I try to print conc_e inside the second DO loop over i, I get REAL values but when I put the print statement right after the ENDDO,
I get NaN. Th two DO loops over i seem identical. I have attached the input file just in case. I may be missing something obvious but I cannot find it. Any hints would be appreciated.

Thank you,
Vahid

PROGRAM conc_fermi
IMPLICIT NONE

INTEGER, PARAMETER :: max_rows = 200000
INTEGER :: row, rowc_min, rowc_max, rowv_min, rowv_max
INTEGER :: i,j,io,total_row
REAL, DIMENSION(max_rows) :: E, DOS, DOSTOT
REAL :: EVmin, EVmax, ECmin, ECmax, E_F, conc, conc_e, conc_h
REAL :: DOS_ECmin, DOS_ECmax, DOS_EVmin, DOS_EVmax, V, dE, kT
CHARACTER (20) :: input_file
REAL, PARAMETER :: temp=300 ! K
REAL, PARAMETER :: kb = 1.3806488e-23 ! J/K
REAL, PARAMETER :: qq = 1.60217657e-19 ! eV
REAL, PARAMETER :: bohr_2_cm = 5.2917721067e-9

kT = kb*temp/qq ! eV
EVmax=6.28
EVmin=-5.83
ECmax=21.47
ECmin=6.8
V=262.8707

PRINT*, "Enter DOS file name"
READ*, input_file

! Convert V from bohr^3 to cm^3
V=V*(bohr_2_cm)**3

OPEN(UNIT=100, FILE=input_file, STATUS='old', ACTION='read')

! Read the DOS file
E)=0
DOS)=0
DOSTOT)=0
DO row=1,max_rows
READ(100,*,IOSTAT=io) E(row), DOS(row), DOSTOT(row)

IF (io .LT. 0) THEN
total_row = row - 1
EXIT
ENDIF

! Define row number and DOS at EVmax
IF ( E(row) .EQ. EVmax ) THEN
rowv_max = row
DOS_EVmax = DOS(rowv_max)
ENDIF

! Define row number and DOS at EVmin
IF ( E(row) .EQ. EVmin ) THEN
rowv_min = row
DOS_EVmin = DOS(rowv_min)
ENDIF

! Define row number and DOS at ECmax
IF ( E(row) .EQ. ECmax ) THEN
rowc_max = row
DOS_ECmax = DOS(rowc_max)
ENDIF

! Define row number and DOS at ECmin
IF ( E(row) .EQ. ECmin ) THEN
rowc_min = row
DOS_ECmin = DOS(rowc_min)
ENDIF

ENDDO
CLOSE(100)

! energy interval in DOS
dE = ABS(E(1)-E(2))

OPEN(1,file='concentration-Fermi.dat')
WRITE(1,*) "Fermi Level (eV) Electron concentration (in 1/cm^3)"

! Vary E_F from 0.3eV below to 0.3eV above ECmin

DO j = -3, 3

E_F=ECmin+(REAL(j)/10)

! Trapezoid Rule
conc_h = 0.0
DO i = rowv_min, rowv_max
conc_h = conc_h + DOS(i) * (1-(1/(1+exp((E(i)-E_F)/kT))))
END DO
PRINT*, conc_h
conc_h=dE*(conc_h-0.5*(DOS_EVmin*(1-(1/(1+exp((EVmin-E_F)/kT)))) &
+ DOS_EVmax*(1-(1/(1+exp((EVmax-E_F)/kT))))))

conc_e = 0.0
DO i = rowc_min, rowc_max
conc_e = conc_e + DOS(i) * (1/(1+exp((E(i)-E_F)/kT)))
! PRINT*, conc_e ! prints all the conc_e
END DO
PRINT*, conc_e ! prints NaN
conc_e = dE*(conc_e-0.5*(DOS_ECmin*(1/(1+exp((ECmin-E_F)/kT))) &
+ DOS_ECmax*(1/(1+exp((ECmax-E_F)/kT)))))
conc=conc_h-conc_e

WRITE(1,*) E_F, " ", conc/V
ENDDO !j

CLOSE(1)
END

 

[mikrom]

I modified the print statements and tried your program

...
     DO i = rowv_min, rowv_max
       conc_h = conc_h + DOS(i) * (1-(1/(1+exp((E(i)-E_F)/kT)))) 
     END DO
     [highlight]PRINT*, 'conc_h =', conc_h[/highlight]
...
     DO i = rowc_min, rowc_max
       conc_e = conc_e + DOS(i) * (1/(1+exp((E(i)-E_F)/kT)))
       ! PRINT*, conc_e ! prints all the conc_e
     END DO
     [highlight]PRINT*, 'conc_e =', conc_e ! prints NaN[/highlight]
...

and got normal output without an error:

$ gfortran vahid9.f95 -o vahid9

$ vahid9
 Enter DOS file name
si.dos
 conc_h =   1.76246031E-05
 conc_e =   6.62863044E-07
 conc_h =   3.68310367E-07
 conc_e =   3.17181075E-05
 conc_h =   7.69643638E-09
 conc_e =   1.50968332E-03
 conc_h =   1.60829211E-10
 conc_e =   5.95530570E-02
 conc_h =   3.36078019E-12
 conc_e =  0.701390862
 conc_h =   7.02287620E-14
 conc_e =   2.32466984
 conc_h =   1.46748299E-15
 conc_e =   4.66214800

 

[Vahid9]

Thanks Mikrom. The NaN seem to be related to compiling with ifort which is the compiler I used. When I compile with gfortran, I get the same result as yours.

Best wishes,
Vahid