TABLE OF CONTENTS
ABINIT/mrgddb [ Programs ]
NAME
mrgddb
FUNCTION
This code merges the derivative databases.
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA, XG, GMR, SP, GA) 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.
SOURCE
36 #if defined HAVE_CONFIG_H 37 #include "config.h" 38 #endif 39 40 #include "abi_common.h" 41 42 program mrgddb 43 44 use defs_basis 45 use m_abicore 46 use m_errors 47 use m_xmpi 48 use m_ddb_hdr 49 use m_ddb 50 51 use m_build_info, only : abinit_version 52 use m_specialmsg, only : herald 53 use m_time , only : asctime, timein 54 use m_io_tools, only : file_exists 55 use m_fstrings, only : sjoin 56 57 implicit none 58 59 !Local variables------------------------------- 60 !scalars 61 integer,parameter :: ddbun=2 62 integer :: chkopt,ios, mddb 63 integer :: iddb,ii,nddb,nfiles_cli,nargs,comm,my_rank 64 real(dp) :: tcpu,tcpui,twall,twalli 65 logical :: cannot_overwrite=.True. 66 character(len=24) :: codename 67 character(len=fnlen) :: dscrpt, outname 68 !arrays 69 real(dp) :: tsec(2) 70 character(len=fnlen),allocatable :: filnam(:),copy_filnam(:) 71 character(len=500) :: msg,arg 72 73 !****************************************************************** 74 75 ! Change communicator for I/O (mandatory!) 76 call abi_io_redirect(new_io_comm=xmpi_world) 77 78 ! Initialize MPI 79 call xmpi_init() 80 comm = xmpi_world; my_rank = xmpi_comm_rank(comm) 81 82 ! Initialize memory profiling if it is activated 83 ! if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0" 84 ! note that abimem.mocc files can easily be multiple GB in size so don't use this option normally 85 #ifdef HAVE_MEM_PROFILING 86 call abimem_init(0) 87 #endif 88 89 call timein(tcpui,twalli) 90 91 codename='MRGDDB'//repeat(' ',18) 92 call herald(codename,abinit_version,std_out) 93 94 ABI_CHECK(xmpi_comm_size(comm)==1, "mrgddb not programmed for parallel execution") 95 96 nargs = command_argument_count() 97 98 mddb = 5000 ! maximum number of databases (initial guess) 99 100 chkopt = 1; nfiles_cli = 0 101 do ii=1,nargs 102 call get_command_argument(ii, arg) 103 if (arg == "-v" .or. arg == "--version") then 104 write(std_out,"(a)") trim(abinit_version); goto 100 105 106 else if (arg == "--nostrict") then 107 ! Disable consistency checks 108 chkopt = 0 109 110 else if (arg == "-f") then 111 cannot_overwrite = .False. 112 113 else if (arg == "-h" .or. arg == "--help") then 114 ! Document the options. 115 write(std_out,*)"Usage:" 116 write(std_out,*)" mrgddb Interactive prompt." 117 write(std_out,*)" mrgddb < run.files Read arguments from run.files." 118 write(std_out,*)" mrgddb out_DDB in1_DDB in2_DDB Merge list of input DDB files, produce new out_DDB file." 119 write(std_out,*)" mrgddb out_DDB in*_DDB Same as above but use shell wildcards instead of file list." 120 write(std_out,*)" " 121 write(std_out,*)"Available options:" 122 write(std_out,*)" -v, --version Show version number and exit." 123 write(std_out,*)" -f Overwrite output DDB if file already exists." 124 write(std_out,*)" --nostrict Disable consistency checks" 125 write(std_out,*)" -h, --help Show this help and exit." 126 goto 100 127 128 else 129 ! Save filenames passed via command-line. 130 if (.not. allocated(filnam)) then 131 ABI_MALLOC(filnam, (mddb+1)) 132 end if 133 nfiles_cli = nfiles_cli + 1 134 if (nfiles_cli > mddb+1) then 135 ! Extend filnam 136 ABI_MALLOC(copy_filnam, (mddb+1)) 137 copy_filnam = filnam 138 iddb = mddb + 1 139 mddb = 2 * mddb 140 ABI_FREE(filnam) 141 ABI_MALLOC(filnam, (mddb+1)) 142 filnam(:iddb) = copy_filnam(:iddb) 143 ABI_FREE(copy_filnam) 144 end if 145 filnam(nfiles_cli) = trim(arg) 146 end if 147 end do 148 149 if (nfiles_cli == 0) then 150 ! Read names of files, operating mode (also check its value), 151 ! and short description of new database. 152 153 ! Read the name of the output ddb 154 write(std_out,*)' Give name for output derivative database : ' 155 read(std_in, '(a)' ) outname 156 write(std_out,'(2a)' )' ',trim(outname) 157 158 ! Read the description of the derivative database 159 write(std_out,*)' Give short description of the derivative database :' 160 read(std_in, '(a)' )dscrpt 161 write(std_out,'(2a)' )' ',trim(dscrpt) 162 163 ! Read the number of input ddbs, and check its value 164 ! MG NOTE: In the documentation of mrgddb_init I found: 165 ! 166 ! nddb = (=1 => will initialize the ddb, using an input GS file) 167 ! (>1 => will merge the whole set of ddbs listed) 168 ! if nddb==1, 169 ! (2) Formatted input file for the Corning ground-state code 170 ! 171 ! but the case nddb=1 with input file is not supported anymore! 172 173 write(std_out,*)' Give number of input ddbs' 174 read(std_in,*)nddb 175 write(std_out,*)nddb 176 ABI_MALLOC(filnam, (nddb)) 177 178 ! Read the file names 179 do iddb=1,nddb 180 !Added to catch error message if the number of input ddbs is greater than the 181 !actually number of ddb files entered by the user. 182 read(std_in, '(a)',IOSTAT =ios ) filnam(iddb) 183 if (ios < 0) then 184 write(msg, '(a,i0,4a)' )& 185 & 'The number of input ddb files: ',nddb,' exceeds the number ',& 186 & 'of ddb file names.', ch10, & 187 & 'Action: change the number of ddb files in the mrgddb input file.' 188 ABI_ERROR(msg) 189 else 190 write(std_out,*)' Give name for derivative database number',iddb,' : ' 191 write(std_out,'(2a)' )' ',trim(filnam(iddb)) 192 end if 193 end do 194 195 else 196 ! Command-line interface. 197 if (nfiles_cli == 1) then 198 ABI_ERROR("Need more than one argument") 199 end if 200 if (cannot_overwrite .and. file_exists(outname)) then 201 ABI_ERROR(sjoin("Cannot overwrite existing file:", outname)) 202 end if 203 nddb = nfiles_cli - 1 204 dscrpt = sjoin("Generated by mrgddb on:", asctime()) 205 end if 206 207 ! Call the main merging routine 208 call merge_ddb(nddb, filnam, outname, dscrpt, chkopt) 209 210 ABI_FREE(filnam) 211 212 !********************************************************************** 213 214 call timein(tcpu,twall) 215 216 tsec(1)=tcpu-tcpui 217 tsec(2)=twall-twalli 218 219 write(std_out, '(3a,f13.1,a,f13.1)' ) & 220 & '-',ch10,'- Proc. 0 individual time (sec): cpu=',tsec(1),' wall=',tsec(2) 221 call wrtout(std_out,'+mrgddb : the run completed successfully ','COLL', do_flush=.True.) 222 223 call abinit_doctor("__mrgddb") 224 225 100 call xmpi_end() 226 227 end program mrgddb