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!!** There should be two separate routines for reading and writing 10 11 module m_ioxv 12 13 logical, public, save :: xv_file_read = .false. 14 15 public :: ioxv 16 17 CONTAINS 18 19 subroutine ioxv( task, cell, vcell, 20 . na, isa, iza, xa, va, found ) 21c ******************************************************************* 22c Saves positions and velocities. 23c J.M.Soler. July 1997. 24c ********** INPUT ************************************************** 25c character task*(*) : 'read' or 'write' (or 'READ' or 'WRITE') 26c ********** INPUT OR OUTPUT (depending of task) ******************** 27c real*8 cell(3,3) : Unit cell vectors 28c real*8 vcell(3,3) : Unit cell vector velocities (Parrinello-Rahman) 29c integer na : Number of atoms 30c integer isa(na) : Atomic species index 31c integer iza(na) : Atomic numbers 32c real*8 xa(3,na) : Atomic positions 33c real*8 va(3,na) : Atomic velocities 34c ********** OUTPUT ************************************************* 35c logical found : Has input file been found 36c (only for task='read') 37c ********** UNITS ************************************************** 38c Units are arbitrary, but the use with task='write' and task='read' 39c must be consistent 40c ******************************************************************* 41 42C 43C Modules 44C 45 use files, only : slabel, label_length 46 use parallel, only : Node 47 use precision, only : dp 48#ifdef MPI 49 use mpi_siesta 50#endif 51 52 implicit none 53 54 logical found 55 integer na, isa(na), iza(na) 56 real(dp) cell(3,3), va(3,na), vcell(3,3), xa(3,na) 57 character(len=*) :: task 58 59 external io_assign, io_close 60 61c Internal variables and arrays 62 character(len=label_length+3) :: fname 63 64 integer ia, iu, iv, ix 65#ifdef MPI 66 integer MPIerror 67#endif 68 69C Find name of file 70 fname = trim(slabel) // '.XV' 71 72C Only do reading and writing for IOnode 73 if (Node.eq.0) then 74 75C Choose between read or write 76 if (task.eq.'read' .or. task.eq.'READ') then 77 78C Check if input file exists 79 inquire( file=fname, exist=found ) 80 if (found) then 81 82C Open file 83 call io_assign( iu ) 84 open( iu, file=fname, status='old' ) 85 86C Read data 87 write(6,'(/,a)') 88 . 'ioxv: Reading coordinates and velocities from file' 89 do iv = 1,3 90 read(iu,*) (cell(ix,iv),ix=1,3),(vcell(ix,iv),ix=1,3) 91 enddo 92 read(iu,*) na 93 do ia = 1,na 94 read(iu,*) 95 . isa(ia),iza(ia),(xa(ix,ia),ix=1,3),(va(ix,ia),ix=1,3) 96 enddo 97 98C Close file 99 call io_close( iu ) 100 101 else 102C If input file not found, go to exit point 103 goto 999 104 endif 105 106 elseif (task.eq.'write' .or. task.eq.'WRITE') then 107 108C Open file 109 call io_assign( iu ) 110 open( iu, file=fname, form='formatted', status='unknown' ) 111 112C Write data on file 113 write(iu,'(2(3x,3f18.9))') 114 . ((cell(ix,iv),ix=1,3),(vcell(ix,iv),ix=1,3),iv=1,3) 115 write(iu,*) na 116 do ia = 1,na 117 write(iu,'(i3,i6,3f18.9,3x,3f18.9)') 118 . isa(ia),iza(ia),(xa(ix,ia),ix=1,3),(va(ix,ia),ix=1,3) 119 enddo 120 121 call io_close( iu ) 122 123 endif 124 endif 125 126 999 continue 127 128 if (task.eq.'read' .or. task.eq.'READ') then 129 130C If data has been read in then broadcast the values to all Nodes 131#ifdef MPI 132 call MPI_Bcast(found,1,MPI_logical,0,MPI_Comm_World,MPIerror) 133#endif 134 if (found) then 135 xv_file_read = .true. 136#ifdef MPI 137 call MPI_Bcast(na,1,MPI_integer,0,MPI_Comm_World,MPIerror) 138 call MPI_Bcast(cell(1,1),9,MPI_double_precision,0, 139 . MPI_Comm_World,MPIerror) 140 call MPI_Bcast(vcell(1,1),9,MPI_double_precision,0, 141 . MPI_Comm_World,MPIerror) 142 call MPI_Bcast(isa,na,MPI_integer,0,MPI_Comm_World,MPIerror) 143 call MPI_Bcast(iza,na,MPI_integer,0,MPI_Comm_World,MPIerror) 144 call MPI_Bcast(xa(1,1),3*na,MPI_double_precision,0, 145 . MPI_Comm_World,MPIerror) 146 call MPI_Bcast(va(1,1),3*na,MPI_double_precision,0, 147 . MPI_Comm_World,MPIerror) 148#endif 149 endif 150 endif 151 152 end subroutine ioxv 153 154 end module m_ioxv 155 156