1! 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt. 6! See Docs/Contributors.txt for a list of contributors. 7!----------------------------------------------------------------------------- 8! 9! module siesta_master 10! 11! Keeps a repository of coordinates, forces, and other properties, and handles 12! their transmission to/from a master program that uses siesta as a subroutine. 13! Since siesta was written as a main program, not as a subroutine, the flow of 14! data to/from the master program is not straighforward as subroutine arguments 15! and it was simpler to adress it with this module. Furthermore, there are two 16! communication modes available with the master program, through unix pipes and 17! through MPI, what justifies the following flow organization: 18! Using unix pipes (master+fsiesta_pipes compiled together): 19! master <--> fsiesta_pipes <--pipe--> iopipes <--> siesta_master <--> siesta 20! Using MPI (master+siesta compiled together): 21! master <--> fsiesta_mpi <--> siesta_master <--> siesta 22! 23!----------------------------------------------------------------------------- 24! 25! Public procedures provided by this module: 26! coordsFromMaster : Get atomic coordinates from master program 27! forcesToMaster : Send atomic forces to master program 28! getForcesForMaster : Get atomic forces to be sent to master program 29! getPropertyForMaster : Get other properties to be sent to master program 30! setPropertyForMaster : Set other properties to be sent to master program 31! setCoordsFromMaster : Set atomic positions sent by master program 32! setMasterUnits : Set physical units used by master program 33! 34! Public variables: 35! siesta_server : Is siesta a server? 36! siesta_subroutine : Is siesta a subroutine? 37! input_file : Name of fdf data file 38! 39!----------------------------------------------------------------------------- 40! 41! Interfaces of public procedures: 42! 43! subroutine coordsFromMaster( na, xa, cell ) 44! integer, intent(in) :: na ! Number of atoms 45! real(dp),intent(out):: xa(3,na) ! Atomic coordinates 46! real(dp),intent(out):: cell(3,3) ! Unit cell vectors 47! end subroutine coordsFromMaster 48! 49! subroutine forcesToMaster( na, Etot, fa, stress ) 50! integer, intent(in):: na ! Number of atoms 51! real(dp),intent(in):: Etot ! Total energy 52! real(dp),intent(in):: fa(3,na) ! Atomic forces 53! real(dp),intent(in):: stress(3,3) ! Stress tensor 54! end subroutine coordsFromMaster 55! 56! subroutine getForcesForMaster( na, fa, stress ) 57! integer, intent(in) :: na ! Number of atoms 58! real(dp),intent(out):: Etot ! Total energy 59! real(dp),intent(out):: fa(3,na) ! Atomic forces 60! real(dp),intent(out):: stress(3,3) ! Stress tensor 61! end subroutine getForcesForMaster 62! 63! subroutine setCoordsFromMaster( na, xa, cell ) 64! integer, intent(in):: na ! Number of atoms 65! real(dp),intent(in):: xa(3,na) ! Atomic coordinates 66! real(dp),intent(in):: cell(3,3) ! Unit cell vectors 67! end subroutine setCoordsFromMaster 68! 69! subroutine setMasterUnits( xunit, eunit ) 70! character(len=*),intent(in):: xunit ! Physical unit of length 71! character(len=*),intent(in):: eunit ! Physical unit of energy 72! end subroutine setMasterUnits 73! 74! subroutine getPropertyForMaster( property, valueSize, value, units, error ) 75! character(len=*),intent(in) :: property ! Property name 76! integer, intent(in) :: valueSize ! Size of value array 77! real(dp), intent(out):: value(:) ! Property value(s) 78! character(len=*),intent(out):: units ! Physical units 79! character(len=*),intent(out):: error ! Error name 80! end subroutine getPropertyForMaster 81! 82! subroutine setPropertyForMaster( property, valueSize, value, units ) 83! character(len=*),intent(in) :: property ! Property name 84! integer, intent(in) :: valueSize ! Size of value array 85! real(dp), intent(in) :: value(*) ! Property value(s) 86! character(len=*),intent(in) :: units ! Physical units 87! end subroutine setPropertyForMaster 88! 89!----------------------------------------------------------------------------- 90 91MODULE siesta_master 92 93! Used module parameters and procedures 94 use precision, only: dp ! Double precision real kind 95 use sys, only: die ! Termination routine 96 use fdf, only: fdf_get ! Reading fdf-options 97 use fdf, only: fdf_convfac ! Conversion of physical units 98 99 use iosockets, only: coordsFromSocket, forcesToSocket 100 use iopipes, only: coordsFromPipe ! Read coordinates from pipe 101 use iopipes, only: forcesToPipe ! Write forces to pipe 102 103 implicit none 104 105! Public procedures 106PUBLIC :: & 107 coordsFromMaster, &! Get atomic coordinates from master program 108 forcesToMaster, &! Send atomic forces to master program 109 getForcesForMaster, &! Get atomic forces to be sent to master program 110 getPropertyForMaster, &! Get other properties to be sent to master program 111 setPropertyForMaster, &! Set other properties to be sent to master program 112 setCoordsFromMaster, &! Set atomic positions sent by master program 113 setMasterUnits ! Set physical units used by master program 114 115! Public variables 116 logical,public,save:: siesta_server = .false. ! Is siesta a server? 117 logical,public,save:: siesta_subroutine = .false. ! Is siesta a subroutine? 118 character(len=132),public,save:: input_file = ' ' ! fdf data file 119 120PRIVATE ! Nothing is declared public beyond this point 121 122! Derived type to hold a physical property 123 type propType 124 private 125 character(len=32):: name = ' ' 126 character(len=32):: units = ' ' 127 integer :: size = 0 128 real(dp),pointer :: value(:) 129 end type propType 130 131! Physical units used by siesta 132 character(len=*), parameter :: siesta_xunit = 'Bohr' 133 character(len=*), parameter :: siesta_eunit = 'Ry' 134 character(len=*), parameter :: siesta_funit = 'Ry/Bohr' 135 character(len=*), parameter :: siesta_sunit = 'Ry/Bohr**3' 136 137! Global module variables and default values 138 integer, parameter:: maxProps = 100 ! Max number of properties 139 integer, save:: nProps = 0 ! Number of stored properties 140 integer, save:: nAtoms = 0 ! Number of atoms 141 real(dp), save:: ucell(3,3) = 0 ! Unit cell vectors 142 real(dp), save:: stressT(3,3) = 0 ! Stress tensor 143 real(dp), save:: energy = 0 ! Total energy 144 real(dp),pointer, save:: coords(:,:) ! Atomic coordinates 145 real(dp),pointer, save:: forces(:,:) ! Atomic forces 146 type(propType), save:: prop(maxProps) ! Stored physical properties 147 character(len=32),save:: master_xunit = 'Ang' ! Master's unit of length 148 character(len=32),save:: master_eunit = 'eV' ! Master's unit of energy 149 character(len=32),save:: master_funit = 'eV/Ang' ! Master's unit of force 150 character(len=32),save:: master_sunit = 'eV/Ang**3' ! Master's unit of stress 151 152CONTAINS 153 154!----------------------------------------------------------------------------- 155 156subroutine coordsFromMaster( na, xa, cell ) 157 158 implicit none 159 integer, intent(in) :: na ! Number of atoms 160 real(dp),intent(out):: xa(3,na) ! Atomic coordinates 161 real(dp),intent(out):: cell(3,3) ! Unit cell vectors 162 character*32 :: iface ! Interface mode 163 164! In the pipes version, this is where we first know that we are a server 165 siesta_server = .true. 166 167 if (siesta_subroutine) then 168 if (na==nAtoms) then 169 xa = coords 170 cell = ucell 171 else 172 call die('coordsFromMaster: ERROR: number-of-atoms mismatch') 173 endif 174 else 175 iface = fdf_get( "Master.interface", "pipe") 176 if ( iface == "pipe") call coordsFromPipe( na, xa, cell ) 177 if ( iface == "socket") call coordsFromSocket (na, xa, cell ) 178 end if ! (siesta_subroutine) 179 180end subroutine coordsFromMaster 181 182!----------------------------------------------------------------------------- 183 184subroutine forcesToMaster( na, Etot, fa, stress ) 185 186 implicit none 187 integer, intent(in):: na ! Number of atoms 188 real(dp),intent(in):: Etot ! Total energy 189 real(dp),intent(in):: fa(3,na) ! Atomic forces 190 real(dp),intent(in):: stress(3,3) ! Stress tensor 191 character*32 :: iface ! Interface mode 192 193 if (siesta_subroutine) then 194 if (na==nAtoms) then 195 energy = Etot 196 forces = fa 197 stressT = stress 198 else ! (na/=nAtoms) 199 call die('coordsFromMaster: ERROR: number-of-atoms mismatch') 200 endif ! (na==nAtoms) 201 else 202 iface = fdf_get( "Master.interface", "pipe") 203 if ( iface == "pipe") call forcesToPipe( na, Etot, fa, stress ) 204 if ( iface == "socket") call forcesToSocket( na, Etot, fa, stress ) 205 end if ! (siesta_subroutine) 206 207end subroutine forcesToMaster 208 209!----------------------------------------------------------------------------- 210 211subroutine getForcesForMaster( na, Etot, fa, stress ) 212 213 implicit none 214 integer, intent(in) :: na ! Number of atoms 215 real(dp),intent(out):: Etot ! Total energy 216 real(dp),intent(out):: fa(3,na) ! Atomic forces 217 real(dp),intent(out):: stress(3,3) ! Stress tensor 218 219! Check that number of atoms is right, then copy data from repository 220! and convert them to master's physical units 221 if (na==nAtoms) then 222 Etot = energy * fdf_convfac( siesta_eunit, master_eunit ) 223 fa = forces * fdf_convfac( siesta_funit, master_funit ) 224 stress = stressT * fdf_convfac( siesta_sunit, master_sunit ) 225 else 226 call die('getForcesForMaster: ERROR: number-of-atoms mismatch') 227 end if 228 229end subroutine getForcesForMaster 230 231!----------------------------------------------------------------------------- 232 233subroutine setCoordsFromMaster( na, xa, cell ) 234 235 implicit none 236 integer, intent(in):: na ! Number of atoms 237 real(dp),intent(in):: xa(3,na) ! Atomic coordinates 238 real(dp),intent(in):: cell(3,3) ! Unit cell vectors 239 240! Allocate arrays for repository the first time 241 if (nAtoms==0) then 242 nAtoms = na 243 allocate( coords(3,na), forces(3,na) ) 244 end if 245 246! Check that number of atoms is right, then copy data to repository 247! after converting them to siesta physical units 248 if (na==nAtoms) then 249 coords = xa * fdf_convfac( master_xunit, siesta_xunit ) 250 ucell = cell * fdf_convfac( master_xunit, siesta_xunit ) 251 else 252 call die('setCoordsFromMaster: ERROR: number-of-atoms mismatch') 253 end if 254 255end subroutine setCoordsFromMaster 256 257!----------------------------------------------------------------------------- 258 259subroutine setMasterUnits( xunit, eunit ) 260 261 implicit none 262 character(len=*),intent(in):: xunit ! Physical unit of length 263 character(len=*),intent(in):: eunit ! Physical unit of energy 264 265 master_xunit = xunit 266 master_eunit = eunit 267 master_funit = trim(eunit)//'/'//trim(xunit) ! Unit of force 268 master_sunit = trim(eunit)//'/'//trim(xunit)//'**3' ! Unit of stress 269 270end subroutine setMasterUnits 271 272!----------------------------------------------------------------------------- 273 274subroutine getPropertyForMaster( property, valueSize, value, units, error ) 275 276 implicit none 277 character(len=*),intent(in) :: property ! Property name 278 integer, intent(in) :: valueSize ! Size of value array 279 real(dp), intent(out):: value(:) ! Property value(s) 280 character(len=*),intent(out):: units ! Physical units 281 character(len=*),intent(out):: error ! Error name 282 283 logical:: found 284 integer:: iProp 285 286 found = .false. 287 do iProp = 1,nProps 288 if (prop(iProp)%name == property) then 289 if (prop(iProp)%size == valueSize) then 290 value = prop(iProp)%value 291 units = prop(iProp)%units 292 error = ' ' 293 else 294 error = 'wrong_size' 295 end if 296 found = .true. 297 exit ! iProp loop 298 end if 299 end do ! iProp 300 if (.not.found) error = 'unknown_property' 301 302end subroutine getPropertyForMaster 303 304!----------------------------------------------------------------------------- 305 306subroutine setPropertyForMaster( property, valueSize, value, units ) 307 308 implicit none 309 character(len=*),intent(in) :: property ! Property name 310 integer, intent(in) :: valueSize ! Size of value array 311 real(dp), intent(in) :: value(*) ! Property value(s) 312 character(len=*),intent(in) :: units ! Physical units 313 314 logical:: found 315 integer:: iProp 316 317! Check if property is already stored 318 found = .false. 319 do iProp = 1,nProps 320 if (prop(iProp)%name == property) then 321 found = .true. 322 exit ! iProp loop 323 end if 324 end do ! iProp 325 326! Increase the number of stored properties, if necessary 327 if (.not.found) then 328 nProps = nProps + 1 329 iProp = nProps 330 end if 331 332! (Re)allocate the array to store value(s) 333 if (prop(iProp)%size /= valueSize) then 334 if (prop(iProp)%size /= 0) deallocate( prop(iProp)%value ) 335 allocate( prop(iProp)%value(valueSize) ) 336 end if 337 338! Store data 339 prop(iProp)%name = property 340 prop(iProp)%size = valueSize 341 prop(iProp)%value = value(1:valueSize) 342 prop(iProp)%units = units 343 344end subroutine setPropertyForMaster 345 346END MODULE siesta_master 347 348