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!***************************************************************************** 10! module iosockets 11!***************************************************************************** 12! Handles communications between siesta-as-a-subroutine and a driver 13! program, transfering coordinates and forces trough POSIX sockets. 14! The routines that handle the other side of the communication are 15! in file fsiesta_sockets.f90 and/or in the i-PI code. 16! J.M.Soler, A.Garcia, and M.Ceriotti, 2014 17!***************************************************************************** 18! 19! Modules used: 20! use precision, only: dp 21! use parallel, only: IOnode 22! use fdf 23! use f90sockets, only: open_socket, writebuffer, readbuffer 24! use sys, only: die, bye 25! use m_mpi_utils, only: broadcast 26! use cellSubs, only: volcel 27! use mpi_siesta 28! 29! Public derived types: 30! none 31! 32! Public variables and arrays: 33! none 34! 35! Public procedures: 36! coordsFromSocket, &! receive coords and cell vectors from master program 37! forcesToSocket ! send energy, forces, and stress to master program 38! 39!***************************************************************************** 40! subroutine coordsFromSocket( na, xa, cell ) 41! integer, intent(in) :: na ! Number of atoms 42! real(dp), intent(out):: xa(3,na) ! Atomic coordinates (bohr) 43! real(dp), intent(out):: cell(3,3) ! Lattice vectors (bohr) 44!***************************************************************************** 45! subroutine forcesToSocket( na, energy, forces, stress ) 46! integer, intent(in) :: na ! Number of atoms 47! real(dp), intent(in) :: energy ! Total energy (Ry) 48! real(dp), intent(in) :: forces(3,na) ! Atomic forces (Ry/bohr) 49! real(dp), intent(in) :: stress(3,3) ! Stress tensor (Ry/bohr^3) 50!***************************************************************************** 51 52module iosockets 53 54! Used modules 55 use precision, only: dp 56 use parallel, only: IOnode 57 use fdf 58 use f90sockets, only: open_socket, writebuffer, readbuffer, close_socket 59 use sys, only: die, bye 60 use m_mpi_utils, only: broadcast 61 use cellSubs, only: volcel 62#ifdef MPI 63 use mpi_siesta 64#endif 65 66 implicit none 67 68! Public procedures 69PUBLIC :: & 70 coordsFromSocket, &! receive coords and cell vectors from master program 71 forcesToSocket ! send energy, forces, and stress to master program 72 73PRIVATE ! Nothing is declared public beyond this point 74 75! Private module parameters 76! WARNING: IPI_MSGLEN and MSGLEN must be equal in i-PI and fsiesta, respectively 77 integer, parameter :: IPI_MSGLEN = 12 78 integer, parameter :: MSGLEN = 80 79 integer, parameter :: unit_len = 32 80 character(len=unit_len), parameter :: siesta_xunit = 'bohr' 81 character(len=unit_len), parameter :: siesta_eunit = 'ry' 82 character(len=unit_len), parameter :: ipi_xunit = 'bohr' 83 character(len=unit_len), parameter :: ipi_eunit = 'hartree' 84 85! Private module variables 86 integer, save :: socket 87 real(dp),save :: cellv 88 character(len=32), save:: master='unknown', stype='unknown' 89 character(len=unit_len),save:: master_eunit='unknown', & 90 master_xunit='unknown' 91#ifdef MPI 92 integer :: MPIerror 93#endif 94 95CONTAINS 96 97!***************************************************************************** 98 99subroutine coordsFromSocket( na, xa, cell ) 100 ! Reads coordinates from socket 101 implicit none 102 integer, intent(in) :: na ! Number of atoms 103 real(dp), intent(out) :: xa(3,na) ! Atomic coordinates (bohr) 104 real(dp), intent(out) :: cell(3,3) ! Lattice vectors (bohr) 105 106! Local parameters 107 character(len=*),parameter:: myName='coordsFromSocket ' 108 109! Local variables and arrays 110 logical, save :: firstTime = .true. 111 CHARACTER(len=MSGLEN) :: header, message 112 CHARACTER(len=1024) :: host 113 INTEGER :: inet, port, n 114 REAL(dp) :: aux(9), c(9) 115 REAL(dp),ALLOCATABLE :: x(:) 116 117! Open the socket, shared by the receive and send sides of communication 118 if (firstTime) then 119 master = fdf_get( "Master.code", "fsiesta") 120 host = fdf_get( "Master.address", "localhost") 121 port = fdf_get( "Master.port", 10001) 122 stype = fdf_get( "Master.socketType", "inet") 123 124 if (leqi(stype,'unix')) then 125 inet = 0 126 elseif (leqi(stype,'inet')) then 127 inet = 1 128 else 129 call die(myName//'ERROR: unknown socket type '//trim(stype)) 130 endif 131 132 if (IOnode) then 133 write(*,'(/,a,i4,i8,2x,a)') & 134 myName//'opening socket: inet,port,host=',inet,port,trim(host) 135 call open_socket(socket, inet, port, host) 136 endif 137 138 firstTime = .false. 139 end if ! (firstTime .and. IOnode) 140 141! Read header from socket 142 header = '' 143 if (IOnode) then 144 do 145 if (leqi(master,'i-pi')) then 146 call readbuffer(socket, header, IPI_MSGLEN) 147 148 ! Immediately stop if requested 149 if (trim(header)=='EXIT') then ! we are done! 150 call close_socket(socket) 151 call bye(myName//'STOP requested by driver') 152 end if 153 154 if (trim(header)/='STATUS') exit ! do loop 155 message = 'READY' 156 call writebuffer(socket, message, IPI_MSGLEN) 157 elseif (leqi(master,'fsiesta')) then 158 call readbuffer(socket, header) 159 160 ! Immediately stop if requested 161 if (trim(header)=='quit') then 162 call writebuffer(socket,'quitting') 163 call close_socket(socket) 164 call bye(myName//'STOP requested by driver') 165 elseif (trim(header)=='wait') then 166 call writebuffer(socket,'ready') 167 cycle ! do loop 168 else 169 exit ! do loop 170 endif 171 else 172 call die(myName//'ERROR: unknown master') 173 endif ! leqi(master) 174 enddo 175 endif ! IOnode 176#ifdef MPI 177 call MPI_Bcast(header, MSGLEN, MPI_Character, 0, MPI_Comm_World, MPIerror) 178#endif 179 180! Read cell vectors from socket, as a single buffer vector 181 if (IOnode) then 182 if (leqi(master,"i-pi") .and. trim(header)=="POSDATA") then 183 call readbuffer(socket, c, 9) 184 ! Please see below for the cell array 185 ! If this is to be used, it should *ALSO* be transposed! 186 call readbuffer(socket, aux, 9) ! not used in siesta 187 master_xunit = ipi_xunit 188 master_eunit = ipi_eunit 189 elseif (leqi(master,"fsiesta") .and. trim(header)=="begin_coords") then 190 call readbuffer(socket, master_xunit) 191 call readbuffer(socket, master_eunit) 192 call readbuffer(socket, c, 9) 193 else 194 call die(myName//'ERROR: unexpected header: '//trim(header)) 195 end if 196 endif 197 198! Broadcast and copy cell from buffer 199#ifdef MPI 200 call MPI_Bcast(master_xunit, unit_len, MPI_Character,0, & 201 MPI_Comm_World, MPIerror) 202 call MPI_Bcast(master_eunit, unit_len, MPI_Character,0, & 203 MPI_Comm_World, MPIerror) 204 call MPI_Bcast(c,9, MPI_Double_Precision,0, MPI_Comm_World, MPIerror) 205#endif 206 ! I-Pi assumes row-major order of *only* the cell parameters 207 ! So we need to transpose the cell 208 ! This is a very bad practice since the same is happening in 209 ! ASE I-pi implementation. 210 ! The cell is transposed when transfering, but not the coordinates. 211 ! Essentially one could leave out *any* transposes in all implementations 212 ! and it would work! However, this is now a legacy implementation 213 ! that should not be changed, neither in I-pi, nor ASE. 214 cell = TRANSPOSE( RESHAPE( c, (/3,3/) ) ) 215 216! Read and check number of atoms 217 if (ionode) then 218 call readbuffer(socket, n) 219 if (n/=na) call die(myName//'ERROR: unexpected number of atoms') 220 endif 221 222! Read and broadcast atomic coordinates 223 allocate(x(3*na)) 224 if (IOnode) call readbuffer(socket, x, 3*na) 225#ifdef MPI 226 call MPI_Bcast(x,3*na, MPI_Double_Precision,0, MPI_Comm_World, MPIerror) 227#endif 228 xa = RESHAPE( x, (/3,na/) ) 229 deallocate(x) 230 231! Read trailing message 232 if (IOnode .and. leqi(master,'fsiesta')) then 233 call readbuffer(socket, message) 234 if (message/='end_coords') then 235 call die(myName//'ERROR: unexpected trailer:'//trim(message)) 236 end if 237 endif 238 239 if ( IONode ) then 240 ! Print coordinates and cell vectors received 241 write(*,'(/,4a,/,(3f12.6))') myName,'cell (',trim(master_xunit),') =', cell 242 write(*,'(4a,/,(3f12.6))') myName,'coords (',trim(master_xunit),') =', xa 243 end if 244 245! Convert physical units 246 cell = cell * fdf_convfac( master_xunit, siesta_xunit ) 247 xa = xa * fdf_convfac( master_xunit, siesta_xunit ) 248 249! Compute and save the cell volume, to be used by forcesToSocket 250 cellv = volcel(cell) 251 252end subroutine coordsFromSocket 253 254!***************************************************************************** 255 256subroutine forcesToSocket( na, energy, forces, stress ) 257! Writes stress and forces to socket 258 implicit none 259 integer, intent(in) :: na ! Number of atoms 260 real(dp), intent(in) :: energy ! Total energy (Ry) 261 real(dp), intent(in) :: forces(3,na) ! Atomic forces (Ry/bohr) 262 real(dp), intent(in) :: stress(3,3) ! Stress tensor (Ry/bohr^3) 263 264! Local parameters 265 character(len=*),parameter:: myName='forcesToSocket ' 266 267! Local variables and arrays 268 character(len=MSGLEN):: header, message 269 real(dp) :: e, s(9), vir(9) 270 real(dp),allocatable :: f(:) 271 272! Copy input to local variables 273 allocate(f(3*na)) 274 e = energy 275 f(:) = reshape( forces, (/3*na/) ) 276 s(:) = reshape( stress, (/9/) ) 277 278! Convert physical units 279 e = e * fdf_convfac( siesta_eunit, master_eunit ) 280 f(:) = f(:) * fdf_convfac( siesta_eunit, master_eunit ) & 281 / fdf_convfac( siesta_xunit, master_xunit ) 282 s(:) = s(:) * fdf_convfac( siesta_eunit, master_eunit ) & 283 / fdf_convfac( siesta_xunit, master_xunit )**3 284 285! Find virial tensor for i-pi, in master's units 286 vir(:) = -s * cellv * fdf_convfac( siesta_xunit, master_xunit )**3 287 288! Write forces to socket 289 if (IOnode) then 290 if (leqi(master,'i-pi')) then 291 do 292 call readbuffer(socket, header, IPI_MSGLEN) 293 if (trim(header)=='STATUS') then ! inform i-pi of my status 294 message = 'HAVEDATA' 295 call writebuffer(socket,message,IPI_MSGLEN) 296 cycle ! do loop 297 elseif (trim(header)=='GETFORCE') then ! proceed to send forces 298 exit ! do loop 299 else 300 call die(myName//'ERROR: unexpected header from i-pi') 301 endif 302 enddo 303 message = 'FORCEREADY' 304 call writebuffer(socket,message,IPI_MSGLEN) 305 call writebuffer(socket,e) 306 call writebuffer(socket,na) 307 call writebuffer(socket,f,3*na) 308 call writebuffer(socket,vir,9) 309 ! i-pi can also receive an arbitrary string, that will be printed 310 ! out to the "extra" trajectory file. This is useful if you want 311 ! to return additional information, e.g. atomic charges, wannier 312 ! centres, etc. one must return the number of characters, then 313 ! the string. here we just send back zero characters. 314 call writebuffer(socket,0) 315 elseif (leqi(master,'fsiesta')) then 316! do 317! call readbuffer(socket, header) 318! if (trim(header)/='wait') exit 319! enddo 320 call writebuffer(socket,'begin_forces') 321 call writebuffer(socket,e) 322 call writebuffer(socket,s,9) 323 call writebuffer(socket,na) 324 call writebuffer(socket,f,3*na) 325 call writebuffer(socket,'end_forces') 326 else 327 call die(myName//'ERROR: unknown master') 328 endif ! leqi(master) 329 end if ! IOnode 330 331 if ( IONode ) then 332 ! Print energy, forces, and stress tensor sent to master 333 write(*,'(/,a,f12.6)') myName// 'energy ('//trim(master_eunit)//') =', e 334 write(*,'(a,/,(3f12.6))') myName//'stress ('//trim(master_eunit)//'/'//trim(master_xunit)//'^3) =', s 335 write(*,'(a,/,(3f12.6))') myName//'forces ('//trim(master_eunit)//'/'//trim(master_xunit)//') =', f 336 end if 337 deallocate(f) 338 339end subroutine forcesToSocket 340 341end module iosockets 342