TABLE OF CONTENTS
ABINIT/aim [ Programs ]
NAME
aim
FUNCTION
Main routine for Bader Atom-In-Molecule analysis.
COPYRIGHT
Copyright (C) 2002-2024 ABINIT group (PCasek,FF,XG) 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) WARNING ABINIT rules are not yet followed in the present routine.
SOURCE
26 #if defined HAVE_CONFIG_H 27 #include "config.h" 28 #endif 29 30 #include "abi_common.h" 31 32 program aim 33 34 use defs_basis 35 use m_abicore 36 use m_xmpi 37 use m_errors 38 use m_nctk 39 #ifdef HAVE_NETCDF 40 use netcdf 41 #endif 42 43 use m_build_info, only : abinit_version 44 use m_time, only : timein 45 use m_io_tools, only : open_file, file_exists 46 use m_specialmsg, only : specialmsg_getcount, herald 47 use m_fstrings, only : int2char4 48 use m_bader !, only : adini, drvaim, inpar, defad, aim_shutdown 49 50 implicit none 51 52 !Arguments ----------------------------------- 53 54 !Local variables------------------------------- 55 integer,parameter :: master = 0 56 integer :: fin,ii,ios,iunt,ierr,nfcfile 57 ! Allow for maximum of 100 fc files 58 integer,parameter :: natm=500 59 integer :: lenstr,me,nproc,comm 60 real(dp) :: tcpu,tcpui,twall,twalli 61 real(dp) :: tsec(2) 62 character(len=fnlen) :: dnfile,hname,infile,ofile,tmpfilename 63 character(len=strlen) :: instr 64 character(len=10) :: procstr 65 character(len=24) :: codename 66 type(aim_dataset_type) :: aim_dtset 67 character(len=500) :: msg 68 character(len=fnlen) :: fcfile(natm) 69 70 !****************************************************************** 71 72 !Change communicator for I/O (mandatory!) 73 call abi_io_redirect(new_io_comm=xmpi_world) 74 75 !Initialize MPI 76 call xmpi_init() 77 78 comm = xmpi_world 79 80 me=xmpi_comm_rank(comm); nproc=xmpi_comm_size(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 !Initialize the code, master only : write heading, and read names of files. 92 if (me == master) then 93 read(*,'(a)') infile 94 infile = trim(infile) 95 read(*,'(a)') dnfile 96 dnfile = trim(dnfile) 97 read(*,'(a)') ofile 98 ofile = trim(ofile) 99 do ii=1,natm 100 iunt=unt+ii 101 read(*,'(a)',iostat=ios) fcfile(ii) 102 if (ios /=0) exit 103 end do 104 nfcfile=ii-1 105 end if 106 107 !Transfer file names to other procs 108 call xmpi_bcast (infile, master, comm, ierr) 109 call xmpi_bcast (dnfile, master, comm, ierr) 110 call xmpi_bcast (ofile, master, comm, ierr) 111 call xmpi_bcast (nfcfile, master, comm, ierr) 112 do ii=1,nfcfile 113 call xmpi_bcast (fcfile(ii), master, comm, ierr) 114 end do 115 116 !Prepare initialization of the log and output files 117 fin=len_trim(ofile) 118 hname(1:fin)=ofile(1:fin) 119 hname(fin+1:fin+1)='.' 120 untout=14 121 codename='AIM '//repeat(' ',18) 122 123 !Open main output file and main log file, then print herald at top of files 124 if(me==master)then 125 126 hname(fin+2:fin+4)='out' 127 if (open_file(hname(1:fin+4),msg,unit=untout,status='unknown',form='formatted') /= 0) then 128 ABI_ERROR(msg) 129 end if 130 rewind (unit=untout) 131 call herald(codename,abinit_version,untout) 132 call herald(codename,abinit_version,std_out) 133 end if 134 135 !Open log file for non-master procs, then print herald at top of files 136 if (me /= master) then 137 call int2char4(me, procstr) 138 tmpfilename = hname(1:fin) // "_LOG_P" // trim(procstr) 139 !close(std_out) 140 if (open_file(tmpfilename, msg, newunit=std_out, status='unknown',form='formatted') /= 0) then 141 ABI_ERROR(msg) 142 end if 143 rewind (unit=std_out) 144 call herald(codename,abinit_version,std_out) 145 end if 146 147 unt=21 ! WARNING : this number is used to define other unit numbers, in init.f 148 unt0=9 149 150 unto=6 ! XG020629 use standard IO unit => easier testing . So, should replace everywhere unto by std_out 151 152 untc=11 153 unts=12 154 untad=19 155 untd=17 156 untl=18 157 unta=15 158 untp=13 159 untg=20 160 161 !OPENING OF THE INPUT FILES 162 163 aim_iomode = IO_MODE_FORTRAN 164 if(me==master)then 165 if (open_file(infile,msg,unit=unt0,status='old',form='formatted') /= 0) then 166 ABI_ERROR(msg) 167 end if 168 169 if (file_exists(dnfile)) then 170 if (open_file(dnfile,msg,unit=untad,status='old',form='unformatted') /=0) then 171 ABI_ERROR(msg) 172 end if 173 else 174 if (file_exists(nctk_ncify(dnfile))) then 175 ! Use netcdf-io 176 write(std_out,"(3a)")"- File: ",trim(dnfile)," does not exist but found netcdf file with similar name." 177 dnfile = nctk_ncify(dnfile) 178 aim_iomode = IO_MODE_ETSF 179 #ifdef HAVE_NETCDF 180 NCF_CHECK(nctk_open_read(untad, dnfile, xmpi_comm_self)) 181 #else 182 ABI_ERROR("Cannot read netcdf file because netcdf support in Abinit is missing.") 183 #endif 184 else 185 ABI_ERROR('Missing data file: '//TRIM(dnfile)) 186 end if 187 end if 188 189 do ii=1,nfcfile 190 iunt=unt+ii 191 if (open_file(fcfile(ii),msg,unit=iunt,status='old',form='formatted') /= 0) then 192 ABI_ERROR(msg) 193 end if 194 end do 195 end if 196 call xmpi_bcast(aim_iomode, master, comm, ierr) 197 198 !call dump_config(std_out) 199 200 !READING OF THE MAIN INPUT FILE 201 202 !Setting the input variables to their default values 203 call defad(aim_dtset) 204 205 !Reading of the input file -> one string called instr 206 if(me==master)then 207 call inpar(instr,lenstr) 208 end if 209 call xmpi_bcast (lenstr, master, comm, ierr) 210 call xmpi_bcast (instr(1:lenstr), master, comm, ierr) 211 212 !Analysis of the input string, setting of input variables in aim_dtset 213 call adini(aim_dtset,instr,lenstr) 214 215 !OPENING OF THE OUTPUT FILES 216 if(me==master)then 217 218 if (aim_dtset%isurf/=0) hname(fin+2:fin+5)='surf' 219 220 if (aim_dtset%isurf==1) then 221 ierr = open_file(hname(1:fin+5),msg,unit=unts,status='unknown',form='formatted') 222 elseif (aim_dtset%isurf==-1) then 223 ierr = open_file(hname(1:fin+5),msg,unit=unts,status='old',action='read',form='formatted') 224 end if 225 if (ierr /= 0) then 226 ABI_ERROR(msg) 227 end if 228 229 if (aim_dtset%crit/=0) hname(fin+2:fin+5)='crit' 230 231 if (aim_dtset%crit>0) then 232 ierr = open_file(hname(1:fin+5),msg,unit=untc,status='unknown',form='formatted') 233 else if (aim_dtset%crit==-1) then 234 ierr = open_file(hname(1:fin+5),msg,unit=untc,status='old',action='read',form='formatted') 235 end if 236 if (ierr /= 0) then 237 ABI_ERROR(msg) 238 end if 239 240 if (aim_dtset%denout==1) then 241 hname(fin+2:fin+4)='dn1' 242 ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted') 243 elseif (aim_dtset%denout==2) then 244 hname(fin+2:fin+4)='dn2' 245 ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted') 246 elseif (aim_dtset%denout==3) then 247 hname(fin+2:fin+4)='dn3' 248 ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted') 249 elseif (aim_dtset%denout==-1) then 250 hname(fin+2:fin+4)='dna' 251 ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='unformatted') 252 end if 253 254 if (ierr /= 0) then 255 ABI_ERROR(msg) 256 end if 257 258 if (aim_dtset%lapout==1) then 259 hname(fin+2:fin+4)='lp1' 260 ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted') 261 elseif (aim_dtset%lapout==2) then 262 hname(fin+2:fin+4)='lp2' 263 ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted') 264 elseif (aim_dtset%lapout==3) then 265 hname(fin+2:fin+4)='lp3' 266 ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted') 267 elseif (aim_dtset%lapout==-1) then 268 hname(fin+2:fin+4)='lpa' 269 ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='unformatted') 270 end if 271 272 if (ierr /= 0) then 273 ABI_ERROR(msg) 274 end if 275 276 if (aim_dtset%gpsurf==1) then 277 hname(fin+2:fin+3)='gp' 278 if (open_file(hname(1:fin+3),msg,unit=untg,status='unknown',form='formatted') /= 0) then 279 ABI_ERROR(msg) 280 end if 281 end if 282 283 if (aim_dtset%plden==1) then 284 hname(fin+2:fin+4)='pld' 285 if (open_file(hname(1:fin+5),msg,unit=untp,status='unknown',form='formatted') /= 0) then 286 ABI_ERROR(msg) 287 end if 288 end if 289 290 end if 291 292 call timein(tcpu,twall) 293 tsec(1)=tcpu-tcpui 294 tsec(2)=twall-twalli 295 write(std_out, '(5a,f13.1,a,f13.1)' ) & 296 & '-',ch10,'- After reading the input file and opening the output files ',ch10,& 297 & '- Proc. 0 individual time (sec): cpu=',tsec(1),' wall=',tsec(2) 298 299 !MAIN DRIVER OF THE ANALYSIS 300 301 write(std_out,'(a,a,a,a,i4)' )char(10),& 302 & ' aim : read density file ',trim(dnfile),' from unit number ',untad 303 304 call drvaim(aim_dtset,tcpui,twalli) 305 306 !THE TOTAL TIME NEEDED 307 call timein(tcpu,twall) 308 tsec(1)=tcpu-tcpui 309 tsec(2)=twall-twalli 310 write(std_out, '(a,a,a,f13.1,a,f13.1)' ) & 311 & '-',ch10,'- Proc. 0 individual time (sec): cpu=',tsec(1),' wall=',tsec(2) 312 313 write(std_out,'(/," Time needed (seconds) - total, CP analyse, SURF determination:",/,/,"- ",3F16.8)')& 314 & tsec(2),ttcp,ttsrf 315 316 317 if(me==0)then 318 write(untout,*) 319 write(untout,*) "TIME ANALYSIS" 320 write(untout,*) "============" 321 write(untout,'(/," Time needed (seconds) - total, CP analyse, SURF determination:",/,/,"- ",3F16.8)') & 322 & tsec(1),ttcp,ttsrf 323 write(untout,'(a,a,f11.3,a,f11.3,a,a,a,a)') char(10),& 324 & '+Total cpu time',tsec(1),& 325 & ' and wall time',tsec(2),' sec',char(10),char(10),& 326 & ' aim : the run completed succesfully.' 327 end if 328 329 !CLOSING OF THE FILES 330 if(me==0)then 331 close(unt0) 332 close(untout) 333 if (aim_dtset%isurf/=0) close(unts) 334 if (aim_dtset%crit/=0) close(untc) 335 if (aim_dtset%denout/=0) close(untd) 336 if (aim_dtset%lapout/=0) close(untl) 337 if (aim_dtset%gpsurf/=0) close(untg) 338 if (aim_dtset%plden/=0) close(untp) 339 end if 340 341 call aim_shutdown() 342 343 !call abinit_doctor("__aim") 344 345 !Eventual cleaning of MPI run 346 call xmpi_end() 347 348 end program aim