pot_io.F90 9.74 KB
Newer Older
1
!--------------------------------------------------------------------------------
2
! Copyright (c) 2017 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!
!!! This module is a wrapper for the potential I/O
!!!
!!!                             GM'17
!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

MODULE m_pot_io

   USE m_types
   USE m_juDFT
   USE m_loddop
   USE m_wrtdop
21
   USE m_cdnpot_io_hdf
22
   USE m_cdnpot_io_common
23 24 25
#ifdef CPP_HDF
   USE hdf5
#endif
26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57
   IMPLICIT NONE

   PRIVATE
   PUBLIC readPotential, writePotential
   PUBLIC POT_ARCHIVE_TYPE_TOT_const, POT_ARCHIVE_TYPE_COUL_const
   PUBLIC POT_ARCHIVE_TYPE_X_const

   INTEGER,          PARAMETER :: POT_ARCHIVE_TYPE_TOT_const = 1
   INTEGER,          PARAMETER :: POT_ARCHIVE_TYPE_COUL_const = 2
   INTEGER,          PARAMETER :: POT_ARCHIVE_TYPE_X_const = 3

   INTEGER,          PARAMETER :: POT_DIRECT_MODE = 1
   INTEGER,          PARAMETER :: POT_STREAM_MODE = 2
   INTEGER,          PARAMETER :: POT_HDF5_MODE   = 3

   CONTAINS

   SUBROUTINE readPotential(stars,vacuum,atoms,sphhar,input,sym,archiveType,&
                            iter,fr,fpw,fz,fzxy)

      TYPE(t_stars),INTENT(IN)  :: stars
      TYPE(t_vacuum),INTENT(IN) :: vacuum
      TYPE(t_atoms),INTENT(IN)  :: atoms
      TYPE(t_sphhar),INTENT(IN) :: sphhar
      TYPE(t_input),INTENT(IN)  :: input
      TYPE(t_sym),INTENT(IN)    :: sym

      INTEGER, INTENT (OUT)     :: iter
      INTEGER, INTENT (IN)      :: archiveType

      !     ..
      !     .. Array Arguments ..
58
      COMPLEX, INTENT (OUT) :: fpw(stars%ng3,input%jspins), fzxy(vacuum%nmzxyd,stars%ng2-1,2,input%jspins)
59 60 61
      REAL,    INTENT (OUT) :: fr(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins), fz(vacuum%nmzd,2,input%jspins)

      ! local variables
62
      INTEGER           :: mode, iUnit
63 64 65
      LOGICAL           :: l_exist
      CHARACTER(len=30) :: filename

66 67 68 69
#ifdef CPP_HDF
      INTEGER(HID_T)    :: fileID
#endif
      INTEGER           :: currentStarsIndex,currentLatharmsIndex
70
      INTEGER           :: currentStructureIndex,currentStepfunctionIndex
71 72 73
      INTEGER           :: potentialType
      CHARACTER(LEN=30) :: archiveName

74 75 76
      CALL getMode(mode)

      IF(mode.EQ.POT_HDF5_MODE) THEN
77
#ifdef CPP_HDF
78 79
         INQUIRE(FILE='pot.hdf',EXIST=l_exist)
         IF (l_exist) THEN
80 81
            CALL openPOT_HDF(fileID,currentStarsIndex,currentLatharmsIndex,&
                             currentStructureIndex,currentStepfunctionIndex)
82

83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101
            archiveName = 'illegalPotentialArchive'
            IF (archiveType.EQ.POT_ARCHIVE_TYPE_TOT_const) THEN
               archiveName = 'pottot'
            END IF
            IF (archiveType.EQ.POT_ARCHIVE_TYPE_COUL_const) THEN
               archiveName = 'potcoul'
            END IF
            IF (archiveType.EQ.POT_ARCHIVE_TYPE_X_const) THEN
               archiveName = 'potx'
            END IF

            potentialType = POTENTIAL_TYPE_IN_const

            l_exist = isPotentialEntryPresentHDF(fileID,archiveName,potentialType)

            CALL closeCDNPOT_HDF(fileID)
         END IF

         IF(l_exist) THEN
102 103
            CALL openPOT_HDF(fileID,currentStarsIndex,currentLatharmsIndex,&
                             currentStructureIndex,currentStepfunctionIndex)
104 105 106 107 108

            CALL readPotentialHDF(fileID, archiveName, potentialType,&
                                  iter,fr,fpw,fz,fzxy)

            CALL closeCDNPOT_HDF(fileID)
109
         ELSE
110
            WRITE(*,*) 'Potential entry or pot.hdf file not found.'
111 112 113
            WRITE(*,*) 'Falling back to stream access file pot.str.'
            mode = POT_STREAM_MODE
         END IF
114 115 116 117 118
#else
         WRITE(*,*) 'Not compiled for pot.hdf file usage.'
         WRITE(*,*) 'Falling back to stream access file pot.str.'
         mode = POT_STREAM_MODE
#endif
119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156
      END IF

      IF(mode.EQ.POT_STREAM_MODE) THEN
         INQUIRE(FILE='pot.str',EXIST=l_exist)
         IF (l_exist) THEN
            !load density from cdn.str and exit subroutine

            RETURN
         ELSE
            WRITE(*,*) 'pot.str file not found.'
            WRITE(*,*) 'Falling back to direct access file.'
            mode = POT_DIRECT_MODE
         END IF
      END IF

      IF (mode.EQ.POT_DIRECT_MODE) THEN
         filename = 'illegalPotentialArchive'
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_TOT_const) THEN
            filename = 'pottot'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_COUL_const) THEN
            filename = 'potcoul'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_X_const) THEN
            filename = 'potx'
         END IF

         INQUIRE(file=TRIM(ADJUSTL(filename)),EXIST=l_exist)
         IF(.NOT.l_exist) THEN
            CALL juDFT_error("potential file "//TRIM(ADJUSTL(filename))//" missing",calledby ="readPotential")
         END IF

         iUnit = 11
         OPEN (iUnit,file=TRIM(ADJUSTL(filename)),form='unformatted',status='unknown')

         CALL loddop(stars,vacuum,atoms,sphhar,input,sym,&
                     iUnit,iter,fr,fpw,fz,fzxy)

157
         CLOSE(iUnit)
158 159 160 161
      END IF

   END SUBROUTINE readPotential

162
   SUBROUTINE writePotential(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,&
163 164 165 166 167
                             iter,fr,fpw,fz,fzxy)

      TYPE(t_stars),INTENT(IN)  :: stars
      TYPE(t_vacuum),INTENT(IN) :: vacuum
      TYPE(t_atoms),INTENT(IN)  :: atoms
168
      TYPE(t_cell), INTENT(IN)  :: cell
169 170 171
      TYPE(t_sphhar),INTENT(IN) :: sphhar
      TYPE(t_input),INTENT(IN)  :: input
      TYPE(t_sym),INTENT(IN)    :: sym
172
      TYPE(t_oneD),INTENT(IN)   :: oneD
173 174 175 176 177

      INTEGER, INTENT (IN)      :: iter
      INTEGER, INTENT (IN)      :: archiveType
      !     ..
      !     .. Array Arguments ..
178
      COMPLEX, INTENT (IN) :: fpw(stars%ng3,input%jspins), fzxy(vacuum%nmzxyd,stars%ng2-1,2,input%jspins)
179 180 181 182
      REAL,    INTENT (IN) :: fr(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins), fz(vacuum%nmzd,2,input%jspins)

      ! local variables
      INTEGER           :: mode, iUnit
183
      LOGICAL           :: l_exist, l_storeIndices
184 185
      CHARACTER(len=30) :: filename

186 187 188 189
#ifdef CPP_HDF
      INTEGER(HID_T)    :: fileID
#endif
      INTEGER           :: currentStarsIndex,currentLatharmsIndex
190
      INTEGER           :: currentStructureIndex,currentStepfunctionIndex
191 192 193
      INTEGER           :: potentialType
      CHARACTER(LEN=30) :: archiveName

194 195 196
      REAL              :: fzTemp(vacuum%nmzd,2,input%jspins)
      COMPLEX           :: fzxyTemp(vacuum%nmzxyd,stars%ng2-1,2,input%jspins)

197 198 199
      CALL getMode(mode)

      IF(mode.EQ.POT_HDF5_MODE) THEN
200
#ifdef CPP_HDF
201 202
         CALL openPOT_HDF(fileID,currentStarsIndex,currentLatharmsIndex,&
                          currentStructureIndex,currentStepfunctionIndex)
203

204 205 206
         CALL checkAndWriteMetadataHDF(fileID, input, atoms, cell, vacuum, oneD, stars, sphhar, sym,&
                                       currentStarsIndex,currentLatharmsIndex,currentStructureIndex,&
                                       currentStepfunctionIndex,l_storeIndices)
207 208 209 210 211 212 213 214 215 216 217 218 219 220

         archiveName = 'illegalPotentialArchive'
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_TOT_const) THEN
            archiveName = 'pottot'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_COUL_const) THEN
            archiveName = 'potcoul'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_X_const) THEN
            archiveName = 'potx'
         END IF

         potentialType = POTENTIAL_TYPE_IN_const

221 222 223 224 225 226 227 228 229 230
         fzTemp(:,:,:) = fz(:,:,:)
         fzxyTemp(:,:,:,:) = fzxy(:,:,:,:)
         IF(vacuum%nvac.EQ.1) THEN
            fzTemp(:,2,:)=fzTemp(:,1,:)
            IF (sym%invs) THEN
               fzxyTemp(:,:,2,:) = CONJG(fzxyTemp(:,:,1,:))
            ELSE
               fzxyTemp(:,:,2,:) = fzxyTemp(:,:,1,:)
            END IF
         END IF
231 232
         CALL writePotentialHDF(input, fileID, archiveName, potentialType,&
                                currentStarsIndex, currentLatharmsIndex, currentStructureIndex,&
233
                                currentStepfunctionIndex,iter,fr,fpw,fzTemp,fzxyTemp)
234 235 236

         IF(l_storeIndices) THEN
            CALL writePOTHeaderData(fileID,currentStarsIndex,currentLatharmsIndex,&
237
                                    currentStructureIndex,currentStepfunctionIndex)
238 239 240 241
         END IF

         CALL closeCDNPOT_HDF(fileID)
#endif
242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270
      ELSE IF(mode.EQ.POT_STREAM_MODE) THEN
         ! Write potential to pot.str file
         STOP 'POT_STREAM_MODE not yet implemented!'
      ELSE
         ! Direct mode
         filename = 'illegalPotentialArchive'
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_TOT_const) THEN
            filename = 'pottot'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_COUL_const) THEN
            filename = 'potcoul'
         END IF
         IF (archiveType.EQ.POT_ARCHIVE_TYPE_X_const) THEN
            filename = 'potx'
         END IF

         iUnit = 11
         OPEN (iUnit,file=TRIM(ADJUSTL(filename)),form='unformatted',status='unknown')
         CALL wrtdop(stars,vacuum,atoms,sphhar,input,sym,&
                     iUnit,iter,fr,fpw,fz,fzxy)
         CLOSE(iUnit)
      END IF

   END SUBROUTINE writePotential

   SUBROUTINE getMode(mode)
      INTEGER, INTENT(OUT) :: mode

      mode = POT_DIRECT_MODE
271
      IF (juDFT_was_argument("-stream_cdn")) mode=POT_STREAM_MODE
272
      IF (.NOT.juDFT_was_argument("-no_cdn_hdf")) THEN !juDFT_was_argument("-hdf_cdn")) THEN
273 274 275
#ifdef CPP_HDF
         mode=POT_HDF5_MODE
#else
276 277
!         WRITE(*,*) 'Code not compiled with HDF5 support.'
!         WRITE(*,*) 'Falling back to direct access.'
278 279
#endif
      END IF
280 281 282
   END SUBROUTINE getMode

END MODULE m_pot_io