TABLE OF CONTENTS


ABINIT/mrgddb [ Programs ]

[ Top ] [ Programs ]

NAME

 mrgddb

FUNCTION

 This code merges the derivative databases.

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, SP)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt

INPUTS

  (main routine)

OUTPUT

  (main routine)

NOTES

 The heading of the constituted database is read,
 then the heading of the temporary database to be added is read,
 the code check their compatibility, and create a new
 database that mixes the old and the temporary ones.
 This process can be iterated.
 The whole database will be stored in
 central memory. One could introduce a third mode in which
 only the temporary DDB is in central memory, while the
 input DDB is read twice : first to make a table of blocks,
 counting the final number of blocks, and second to merge
 the two DDBs. This would save memory.

PARENTS

CHILDREN

      abi_io_redirect,abimem_init,abinit_doctor,ddb_hdr_free
      ddb_hdr_open_read,get_command_argument,herald,mblktyp1,mblktyp5,timein
      wrtout,xmpi_init

SOURCE

 43 #if defined HAVE_CONFIG_H
 44 #include "config.h"
 45 #endif
 46 
 47 #include "abi_common.h"
 48 
 49 program mrgddb
 50 
 51  use defs_basis
 52  use m_build_info
 53  use m_abicore
 54  use m_errors
 55  use m_xmpi
 56  use m_ddb_hdr
 57 
 58  use m_specialmsg,   only : specialmsg_getcount, herald
 59  use m_time ,        only : asctime, timein
 60  use m_io_tools,     only : file_exists
 61  use m_fstrings,     only : sjoin
 62  use m_ddb,          only : DDB_VERSION, mblktyp1, mblktyp5
 63 
 64 !This section has been created automatically by the script Abilint (TD).
 65 !Do not modify the following lines by hand.
 66 #undef ABI_FUNC
 67 #define ABI_FUNC 'mrgddb'
 68 !End of the abilint section
 69 
 70  implicit none
 71 
 72 !Local variables-------------------------------
 73 !scalars
 74  integer,parameter :: mddb=5000,ddbun=2 ! mddb=maximum number of databases (cannot be made dynamic)
 75  integer :: chkopt,ios
 76  integer :: iddb,ii,mblktyp,mblktyptmp,nddb,nfiles_cli,nargs,msym,comm,my_rank
 77  real(dp) :: tcpu,tcpui,twall,twalli
 78  logical :: cannot_overwrite=.True.
 79  character(len=24) :: codename
 80  character(len=fnlen) :: dscrpt
 81  type(ddb_hdr_type) :: ddb_hdr
 82 !arrays
 83  real(dp) :: tsec(2)
 84  character(len=fnlen) :: filnam(mddb+1)
 85  character(len=500) :: msg,arg
 86 
 87 !******************************************************************
 88 
 89  ! Change communicator for I/O (mandatory!)
 90  call abi_io_redirect(new_io_comm=xmpi_world)
 91 
 92  ! Initialize MPI
 93  call xmpi_init()
 94  comm = xmpi_world; my_rank = xmpi_comm_rank(comm)
 95 
 96  ! Initialize memory profiling if it is activated
 97  ! if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
 98  ! note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
 99 #ifdef HAVE_MEM_PROFILING
100  call abimem_init(0)
101 #endif
102 
103  call timein(tcpui,twalli)
104 
105  codename='MRGDDB'//repeat(' ',18)
106  call herald(codename,abinit_version,std_out)
107 
108  ABI_CHECK(xmpi_comm_size(comm)==1, "mrgddb not programmed for parallel execution")
109 
110  nargs = command_argument_count()
111 
112  chkopt = 1; nfiles_cli = 0
113  do ii=1,nargs
114    call get_command_argument(ii, arg)
115    if (arg == "-v" .or. arg == "--version") then
116      write(std_out,"(a)") trim(abinit_version); goto 100
117 
118    else if (arg == "--nostrict") then
119      ! Disable consistency checks
120      chkopt = 0
121 
122    else if (arg == "-f") then
123      cannot_overwrite = .False.
124 
125    else if (arg == "-h" .or. arg == "--help") then
126      ! Document the options.
127      write(std_out,*)"Usage:"
128      write(std_out,*)"    mrgddb                           Interactive prompt."
129      write(std_out,*)"    mrgddb < run.files               Read arguments from run.files."
130      write(std_out,*)"    mrgddb out_DDB in1_DDB in2_DDB   Merge list of input DDB files, produce new out_DDB file."
131      write(std_out,*)"    mrgddb out_DDB in*_DDB           Same as above but use shell wildcards instead of file list."
132      write(std_out,*)" "
133      write(std_out,*)"Available options:"
134      write(std_out,*)"    -v, --version      Show version number and exit."
135      write(std_out,*)"    -f                 Overwrite output DDB if file already exists."
136      write(std_out,*)"    --nostrict         Disable consistency checks"
137      write(std_out,*)"    -h, --help         Show this help and exit."
138      goto 100
139 
140    else
141      ! Save filenames passed via command-line.
142      nfiles_cli = nfiles_cli + 1
143      if (nfiles_cli > mddb+1) then
144        write(msg, '(a,i0,2a)')&
145        'Number of files should be lower than mddb+1= ',mddb+1,ch10,&
146        'Action: change mddb in mrgddb.f90 and recompile.'
147        MSG_ERROR(msg)
148      end if
149      filnam(nfiles_cli) = trim(arg)
150    end if
151  end do
152 
153  if (nfiles_cli == 0) then
154    ! Read names of files, operating mode (also check its value),
155    ! and short description of new database.
156 
157    ! Read the name of the output ddb
158    write(std_out,*)' Give name for output derivative database : '
159    read(std_in, '(a)' ) filnam(1)
160    write(std_out,'(a,a)' )' ',trim(filnam(1))
161 
162    ! Read the description of the derivative database
163    write(std_out,*)' Give short description of the derivative database :'
164    read(std_in, '(a)' )dscrpt
165    write(std_out,'(a,a)' )' ',trim(dscrpt)
166 
167    ! Read the number of input ddbs, and check its value
168    ! MG NOTE: In the documentation of mrgddb_init I found:
169    !
170    ! nddb = (=1 => will initialize the ddb, using an input GS file)
171    !        (>1 => will merge the whole set of ddbs listed)
172    !    if nddb==1,
173    !     (2) Formatted input file for the Corning ground-state code
174    !
175    ! but the case nddb=1 with input file is not supported anymore!
176 
177    write(std_out,*)' Give number of input ddbs, or 1 if input GS file'
178    read(std_in,*)nddb
179    write(std_out,*)nddb
180    if (nddb<=0 .or. nddb>mddb) then
181      write(msg, '(a,a,i0,a,i0,a,a,a)' )&
182 &     'nddb should be positive, >1 , and lower',&
183 &     'than mddb= ',mddb,' while the input nddb is ',nddb,'.',ch10,&
184 &     'Action: change mddb in mrgddb.F90 and recompile.'
185      MSG_ERROR(msg)
186    end if
187 
188    ! Read the file names
189    if (nddb==1) then
190      write(std_out,*)' Give name for ABINIT input file : '
191      read(std_in, '(a)' ) filnam(2)
192      write(std_out,'(a,a)' )' ',trim(filnam(2))
193    else
194      do iddb=1,nddb
195        !Added to catch error message if the number of input ddbs is greater than the
196        !actually number of ddb files entered by the user.
197        read(std_in, '(a)',IOSTAT =ios ) filnam(iddb+1)
198        if (ios < 0) then
199          write(msg, '(a,i0,a,a,a,a)' )&
200 &         'The number of input ddb files: ',nddb,' exceeds the number ',&
201 &         'of ddb file names.', ch10, &
202 &         'Action: change the number of ddb files in the mrgddb input file.'
203          MSG_ERROR(msg)
204        else
205          write(std_out,*)' Give name for derivative database number',iddb,' : '
206          write(std_out,'(a,a)' )' ',trim(filnam(iddb+1))
207        end if
208      end do
209    end if
210 
211  else
212    ! Command-line interface.
213    if (nfiles_cli == 1) then
214      MSG_ERROR("Need more than one argument")
215    end if
216    if (cannot_overwrite .and. file_exists(filnam(1))) then
217      MSG_ERROR(sjoin("Cannot overwrite existing file:", filnam(1)))
218    end if
219    nddb = nfiles_cli - 1
220    dscrpt = sjoin("Generated by mrgddb on:", asctime())
221  end if
222 
223  ! Evaluate the mblktyp of the databases
224  ! msym = maximum number of symmetry elements in space group
225  mblktyptmp=1
226  do iddb=1,nddb
227    call ddb_hdr_open_read(ddb_hdr,filnam(iddb+1),ddbun,DDB_VERSION,dimonly=1)
228 
229    if (ddb_hdr%mblktyp > mblktyptmp) mblktyptmp = ddb_hdr%mblktyp
230    call ddb_hdr_free(ddb_hdr)
231  end do
232 
233  mblktyp = mblktyptmp
234  ! write(std_out,*),'mblktyp',mblktyp
235 
236  if (mblktyp==5) then
237    ! Memory optimized routine
238    call mblktyp5(chkopt,ddbun,dscrpt,filnam,mddb,msym,nddb,DDB_VERSION)
239  else
240    ! Speed optimized routine
241    call mblktyp1(chkopt,ddbun,dscrpt,filnam,mddb,msym,nddb,DDB_VERSION)
242  end if
243 
244 !**********************************************************************
245 
246  call timein(tcpu,twall)
247 
248  tsec(1)=tcpu-tcpui
249  tsec(2)=twall-twalli
250 
251  write(std_out, '(a,a,a,f13.1,a,f13.1)' ) &
252 & '-',ch10,'- Proc.   0 individual time (sec): cpu=',tsec(1),'  wall=',tsec(2)
253  call wrtout(std_out,'+mrgddb : the run completed successfully ','COLL', do_flush=.True.)
254 
255  call abinit_doctor("__mrgddb")
256 
257  100 call xmpi_end()
258 
259  end program mrgddb