TABLE OF CONTENTS


ABINIT/parsefile [ Functions ]

[ Top ] [ Functions ]

NAME

 parsefile

FUNCTION

  Glue function, to read the given file, put it into a string,
  change everything to uppercase, remove carriage returns and
  non significant blank characters. May also read a XYZ input
  file if specified. Finally read ndtset input variable.

INPUTS

  filnamin= the file to read
  comm=MPI communicator

OUTPUT

  lenstr= the length of the resulting string.
  ndtset= the number of declared datasets.
  string= contains on output the content of the file, ready for parsing.

PARENTS

      abinit,m_ab7_invars_f90,ujdet

CHILDREN

      importxyz,instrng,intagm,inupper,xmpi_bcast

SOURCE

 28 #if defined HAVE_CONFIG_H
 29 #include "config.h"
 30 #endif
 31 
 32 #include "abi_common.h"
 33 
 34 
 35 subroutine parsefile(filnamin,lenstr,ndtset,string,comm)
 36 
 37  use defs_basis
 38  use m_profiling_abi
 39  use m_errors
 40  use m_xmpi
 41 
 42 !This section has been created automatically by the script Abilint (TD).
 43 !Do not modify the following lines by hand.
 44 #undef ABI_FUNC
 45 #define ABI_FUNC 'parsefile'
 46  use interfaces_32_util
 47  use interfaces_42_parser
 48  use interfaces_57_iovars, except_this_one => parsefile
 49 !End of the abilint section
 50 
 51  implicit none
 52 
 53 !Arguments ------------------------------------
 54  character(len=*),intent(in) :: filnamin
 55  integer,intent(in) :: comm
 56  integer,intent(out) :: ndtset,lenstr
 57  character(len=strlen),intent(out) :: string
 58 
 59 !Local variables-------------------------------
 60 !scalars
 61  integer,parameter :: master=0
 62  integer :: option,marr,tread,lenstr_noxyz,ierr
 63  character(len=strlen) :: string_raw
 64  character(len=500) :: message
 65 !arrays
 66  integer :: intarr(1)
 67  real(dp) :: dprarr(1)
 68 
 69 ! *************************************************************************
 70 
 71  ! Read the input file, and store the information in a long string of characters
 72  ! Note: this is done only by me=0, and then string and other output vars are BCASTED
 73 
 74  if (xmpi_comm_rank(comm) == master) then
 75    !strlen from defs_basis module
 76    option=1
 77    call instrng (filnamin,lenstr,option,strlen,string)
 78 
 79    ! Copy original file, without change of case
 80    string_raw=string
 81 
 82    ! To make case-insensitive, map characters of string to upper case:
 83    call inupper(string(1:lenstr))
 84 
 85    ! Might import data from xyz file(s) into string
 86    ! Need string_raw to deal properly with xyz filenames
 87    lenstr_noxyz = lenstr
 88    call importxyz(lenstr,string_raw,string,strlen)
 89    
 90    !6) Take ndtset from the input string
 91    ndtset=0; marr=1
 92    call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),"ndtset",tread,'INT')
 93    if (tread==1) ndtset=intarr(1)
 94    ! Check that ndtset is not negative
 95    if (ndtset<0 .or. ndtset>9999) then
 96      write(message, '(a,i0,a,a,a,a)' )&
 97 &     'Input ndtset must be non-negative and < 10000, but was ',ndtset,ch10,&
 98 &     'This is not allowed.  ',ch10,&
 99 &     'Action : modify ndtset in the input file.'
100      MSG_ERROR(message)
101    end if
102  end if ! master 
103 
104  if (xmpi_comm_size(comm) > 1) then
105    ! Broadcast data.
106    call xmpi_bcast(lenstr,master,comm,ierr)
107    call xmpi_bcast(ndtset,master,comm,ierr)
108    call xmpi_bcast(string,master,comm,ierr)
109  end if 
110 
111 end subroutine parsefile