TABLE OF CONTENTS


ABINIT/ioprof [ Programs ]

[ Top ] [ Programs ]

NAME

 ioprof

FUNCTION

 Tool for profiling and and testing the IO routines used in abinit

COPYRIGHT

 Copyright (C) 2004-2022 ABINIT group (MG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  (main program)

NOTES

SOURCE

 21 #if defined HAVE_CONFIG_H
 22 #include "config.h"
 23 #endif
 24 
 25 #include "abi_common.h"
 26 
 27 program ioprof
 28 
 29  use defs_basis
 30  use m_build_info
 31  use m_errors
 32  use m_xmpi
 33  use m_wfk
 34  use m_abicore
 35  use m_hdr
 36  use netcdf
 37 
 38  use m_specialmsg,     only : specialmsg_getcount, herald
 39  use m_fstrings,       only : lower, sjoin, itoa
 40  use m_io_tools,       only : delete_file, file_exists, iomode_from_fname, get_unit
 41  use m_argparse,       only : get_arg !, get_arg_list, get_start_step_num
 42 
 43  implicit none
 44 
 45 !Local variables-------------------------------
 46 !scalars
 47  integer,parameter :: master=0,MAX_NFILES=50
 48  integer :: comm,my_rank,nprocs,iomode,formeig,ierr,fform
 49  integer :: ii,io,check_iomode,method,feg,ount,nband2read,nargs, abimem_level
 50  logical :: verbose=.FALSE.
 51  real(dp) :: abimem_limit_mb
 52  character(len=24) :: codename
 53  character(len=500) :: msg,command,arg
 54  character(len=fnlen) :: new_fname,wfk_source,wfk_dest,wfk_path
 55  type(hdr_type) :: hdr
 56  type(wfk_t) :: wfk
 57 !arrays
 58  integer,parameter :: formeigs(2) = (/0,1/)
 59  integer,parameter :: io_modes(1) = (/IO_MODE_FORTRAN/)
 60  !integer,parameter :: io_modes(1) = (/IO_MODE_MPI/)
 61  !integer,parameter :: io_modes(1) = (/IO_MODE_ETSF/)
 62  !integer,parameter :: io_modes(2) = (/IO_MODE_FORTRAN, IO_MODE_MPI/)
 63  !integer,parameter :: io_modes(3) = (/IO_MODE_FORTRAN, IO_MODE_MPI, IO_MODE_ETSF/)
 64  real(dp) :: cwtimes(2)
 65  type(kvars_t),allocatable :: Kvars(:)
 66  character(len=fnlen) :: hdr_fnames(MAX_NFILES) = ABI_NOFILE
 67 
 68 ! *************************************************************************
 69 
 70  ! Change communicator for I/O (mandatory!)
 71  call abi_io_redirect(new_io_comm=xmpi_world)
 72 
 73  ! Initialize MPI
 74  call xmpi_init()
 75  comm = xmpi_world; my_rank = xmpi_comm_rank(comm); nprocs  = xmpi_comm_size(comm)
 76 
 77  ! Initialize memory profiling if it is activated
 78  ! if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
 79  ! note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
 80  ABI_CHECK(get_arg("abimem-level", abimem_level, msg, default=0) == 0, msg)
 81  ABI_CHECK(get_arg("abimem-limit-mb", abimem_limit_mb, msg, default=20.0_dp) == 0, msg)
 82 #ifdef HAVE_MEM_PROFILING
 83  call abimem_init(abimem_level, limit_mb=abimem_limit_mb)
 84 #endif
 85 
 86  codename='IOPROF'//REPEAT(' ',17)
 87  call herald(codename,abinit_version,std_out)
 88 
 89  !verbose = .TRUE.
 90  ount = dev_null; if (verbose) ount = std_out
 91 
 92  if (nprocs>1) then
 93    ABI_ERROR("not programmed for parallel execution, Run it in sequential")
 94  end if
 95 
 96  nargs = command_argument_count()
 97 
 98  ! Command line options.
 99  do ii=1,command_argument_count()
100    call get_command_argument(ii, arg)
101    if (arg == "-v" .or. arg == "--version") then
102      write(std_out,"(a)") trim(abinit_version); goto 100
103 
104    else if (arg == "-h" .or. arg == "--help") then
105      ! Document the options.
106      write(std_out,*)"-v, --version                 Show version number and exit."
107      write(std_out,*)"info wfk_file                 Print info on WFK file"
108      write(std_out,*)"nc2fort ncfile fortran_file   Convert WFK ncfile to Fortran WFK file"
109      write(std_out,*)"-h, --help                    Show this help and exit."
110      goto 100
111    end if
112  end do
113 
114  call get_command_argument(1, command)
115 
116  select case (command)
117 
118  case ("info")
119    ! Print info
120    ABI_CHECK(nargs == 2, "Usage: info wfk_file")
121    call get_command_argument(2, wfk_path)
122    if (my_rank == master) then
123      formeig = 0
124      call wfk_open_read(wfk, wfk_path, formeig, iomode_from_fname(wfk_path), get_unit(), xmpi_comm_self)
125      call wfk%print()
126      call wfk%close()
127    end if
128 
129  case ("nc2fort")
130    ! Converter: netcdf --> Fortran
131    ABI_CHECK(nargs == 3, "Usage: nc2form wfk_source wfk_dest")
132    call get_command_argument(2, wfk_source)
133    call get_command_argument(3, wfk_dest)
134 
135    if (file_exists(wfk_dest)) then
136      ABI_ERROR(sjoin("Cannot overwrite:", wfk_dest, ". Remove file and rerun"))
137    end if
138 
139    if (my_rank == master) call wfk_nc2fort(wfk_source, wfk_dest)
140 
141  case ("prof")
142    ! Profile
143    ABI_CHECK(nargs > 1, "Usage: bench wfk_path")
144    call get_command_argument(2, wfk_path)
145    formeig = 0
146    nband2read = 1500
147    call wfk_prof(wfk_path, formeig, nband2read, comm)
148    ! TO CREATE new file
149    !call wfk_create_wfkfile(new_fname,hdr,iomode,formeig,Kvars,cwtimes,xmpi_comm_self)
150 
151  case ("unittests")
152    ! Unit tests.
153    do ii=1,count(hdr_fnames/=ABI_NOFILE)
154 
155      ! Read the header from an external netcdf file
156      ! This trick is needed because the initialization of
157      ! the header is *IMPOSSIBLE* if we don't have a full ABINIT input file!
158      !
159      call hdr_read_from_fname(hdr,hdr_fnames(ii),fform,comm)
160      ABI_CHECK(fform/=0,"fform==0")
161 
162      call hdr%echo(fform,4,unit=std_out)
163 
164      do feg=1,size(formeigs)
165        formeig = formeigs(feg)
166        do io=1,size(io_modes)
167          iomode = io_modes(io)
168 
169          new_fname = "NEW_WFK"
170          if (iomode==IO_MODE_ETSF) new_fname = "NEW_WFK.nc"
171 
172          if (file_exists(new_fname)) then
173            call delete_file(new_fname,ierr)
174          end if
175 
176          ! TODO
177          if (formeig==1 .and. iomode==IO_MODE_ETSF) then
178            ABI_WARNING("iomode==1 with ETSF_IO not coded yet")
179            CYCLE
180          end if
181 
182          ABI_MALLOC(Kvars, (hdr%nkpt))
183 
184          if (my_rank == master) then
185            call wrtout(ount,"Calling wfk_create_wfkfile","COLL")
186            call wfk_create_wfkfile(new_fname,hdr,iomode,formeig,Kvars,cwtimes,xmpi_comm_self)
187            !call wfk_create_wfkfile(new_fname,hdr,IO_MODE_FORTRAN,formeig,Kvars,cwtimes,xmpi_comm_self)
188            call wrtout(ount,"Done wfk_create_wfkfile","COLL")
189          end if
190 
191          if (nprocs > 1) then
192            ABI_ERROR("Not coded")
193            !call kvars_mpicast(Kvars,master,comm)
194          end if
195 
196          method = 0
197 
198          call wrtout(std_out,sjoin("Checking file:",new_fname,", with iomode=",itoa(iomode)))
199          call wfk_check_wfkfile(new_fname,hdr,iomode,method,formeig,Kvars,cwtimes,comm,ierr)
200 
201          write(msg,'(2(a,i2),2(a,f8.2))')&
202          " Read with iomode: ",iomode,", nproc: ",nprocs,", cpu: ",cwtimes(1),", wall:",cwtimes(2)
203          call wrtout(ount,msg,"COLL")
204 
205          if (ierr/=0) then
206            ABI_ERROR(sjoin("wfk_check_wfkfile returned ierr ",itoa(ierr)))
207          end if
208 
209          ! If not netcdf file, try to read the file with the other mode that is compatible with it.
210          check_iomode = -100
211          if (iomode == IO_MODE_FORTRAN) check_iomode = IO_MODE_MPI
212          if (iomode == IO_MODE_MPI)     check_iomode = IO_MODE_FORTRAN
213 
214          if (check_iomode /= -100) then
215            call wrtout(std_out,sjoin("Trying to read file:",new_fname,", with check_iomode=",itoa(check_iomode)))
216 
217            call wfk_check_wfkfile(new_fname,hdr,check_iomode,method,formeig,Kvars,cwtimes,comm,ierr)
218 
219            write(msg,'(2(a,i2),2(a,f8.2))')&
220            " Read with iomode: ",iomode,", nproc: ",nprocs,", cpu: ",cwtimes(1),", wall:",cwtimes(2)
221            call wrtout(ount,msg,"COLL")
222 
223            if (ierr/=0) then
224              ABI_ERROR(sjoin("wfk_check_wfkfile returned ierr ",itoa(ierr)))
225            end if
226          end if
227 
228          ABI_FREE(Kvars)
229        end do ! iomode
230      end do ! formeig
231 
232      call hdr%free()
233    end do
234 
235  case default
236    ABI_ERROR(sjoin("Unknown command:", command))
237  end select
238 
239  !call wrtout(std_out,ch10//" Analysis completed.","COLL")
240 
241 !Writes information on file about the memory before ending mpi module, if memory profiling is enabled
242  !ABI_MALLOC(nband, (2))
243  call abinit_doctor("__ioprof")
244 
245  100 call xmpi_end()
246 
247  end program ioprof