TABLE OF CONTENTS


ABINIT/aim [ Programs ]

[ Top ] [ Programs ]

NAME

 aim

FUNCTION

 Main routine for Bader Atom-In-Molecule analysis.

COPYRIGHT

 Copyright (C) 2002-2018 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.

PARENTS

CHILDREN

      abi_io_redirect,abimem_init,adini,aim_shutdown,defad,drvaim,herald
      inpar,int2char4,timein,xmpi_bcast,xmpi_end,xmpi_init

SOURCE

 33 #if defined HAVE_CONFIG_H
 34 #include "config.h"
 35 #endif
 36 
 37 #include "abi_common.h"
 38 
 39 program aim
 40 
 41  use defs_basis
 42  use defs_abitypes
 43  use m_abicore
 44  use m_xmpi
 45  use m_build_info
 46  use m_errors
 47  use m_nctk
 48 #ifdef HAVE_NETCDF
 49  use netcdf
 50 #endif
 51 
 52  use m_time,     only : timein
 53  use m_io_tools, only : open_file, file_exists
 54  use m_specialmsg,  only : specialmsg_getcount, herald
 55  use m_fstrings, only : int2char4
 56  use m_bader !,    only : adini, drvaim, inpar, defad, aim_shutdown
 57 
 58 !This section has been created automatically by the script Abilint (TD).
 59 !Do not modify the following lines by hand.
 60 #undef ABI_FUNC
 61 #define ABI_FUNC 'aim'
 62 !End of the abilint section
 63 
 64  implicit none
 65 
 66 !Arguments -----------------------------------
 67 
 68 !Local variables-------------------------------
 69  integer,parameter :: master = 0
 70  integer :: fin,ii,ios,iunt,ierr,nfcfile
 71 ! Allow for maximum of 100 fc files
 72  integer,parameter :: natm=500
 73  integer :: lenstr,me,nproc,comm
 74  real(dp) :: tcpu,tcpui,twall,twalli
 75  real(dp) :: tsec(2)
 76  character(len=fnlen) :: dnfile,hname,infile,ofile
 77  character(len=strlen) :: instr
 78  character(len=fnlen) :: tmpfilename
 79  character(len=10) :: procstr
 80  character(len=24) :: codename
 81  type(aim_dataset_type) :: aim_dtset
 82  character(len=500) :: msg
 83  character(len=fnlen) :: fcfile(natm)
 84 
 85 !******************************************************************
 86 
 87 !Change communicator for I/O (mandatory!)
 88  call abi_io_redirect(new_io_comm=xmpi_world)
 89 
 90 !Initialize MPI
 91  call xmpi_init()
 92 
 93  comm = xmpi_world
 94 
 95  me=xmpi_comm_rank(comm); nproc=xmpi_comm_size(comm)
 96 
 97 !Initialize memory profiling if it is activated
 98 !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
 99 !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
100 #ifdef HAVE_MEM_PROFILING
101  call abimem_init(0)
102 #endif
103 
104  call timein(tcpui,twalli)
105 
106 !Initialize the code, master only : write heading, and read names of files.
107  if (me == master) then
108    read(*,'(a)') infile
109    infile = trim(infile)
110    read(*,'(a)') dnfile
111    dnfile = trim(dnfile)
112    read(*,'(a)') ofile
113    ofile = trim(ofile)
114    do ii=1,natm
115      iunt=unt+ii
116      read(*,'(a)',iostat=ios) fcfile(ii)
117      if (ios /=0) exit
118    end do
119    nfcfile=ii-1
120  end if
121 
122 !Transfer file names to other procs
123  call xmpi_bcast (infile, master, comm, ierr)
124  call xmpi_bcast (dnfile, master, comm, ierr)
125  call xmpi_bcast (ofile, master, comm, ierr)
126  call xmpi_bcast (nfcfile, master, comm, ierr)
127  do ii=1,nfcfile
128    call xmpi_bcast (fcfile(ii), master, comm, ierr)
129  end do
130 
131 !Prepare initialization of the log and output files
132  fin=len_trim(ofile)
133  hname(1:fin)=ofile(1:fin)
134  hname(fin+1:fin+1)='.'
135  untout=14
136  codename='AIM   '//repeat(' ',18)
137 
138 !Open main output file and main log file, then print herald at top of files
139  if(me==master)then
140 
141    hname(fin+2:fin+4)='out'
142    if (open_file(hname(1:fin+4),msg,unit=untout,status='unknown',form='formatted') /= 0) then
143      MSG_ERROR(msg)
144    end if
145    rewind (unit=untout)
146    call herald(codename,abinit_version,untout)
147    call herald(codename,abinit_version,std_out)
148  end if
149 
150 !Open log file for non-master procs, then print herald at top of files
151  if (me /= master) then
152    call int2char4(me, procstr)
153    tmpfilename = hname(1:fin) // "_LOG_P" // trim(procstr)
154    !close(std_out)
155    if (open_file(tmpfilename, msg, newunit=std_out, status='unknown',form='formatted') /= 0) then
156      MSG_ERROR(msg)
157    end if
158    rewind (unit=std_out)
159    call herald(codename,abinit_version,std_out)
160  end if
161 
162 !DEBUG
163 !stop   ! Works in paralel up to here
164 !ENDDEBUG
165 
166  unt=21 ! WARNING : this number is used to define other unit numbers, in init.f
167  unt0=9
168 
169  unto=6  ! XG020629 use standard IO unit => easier testing . So, should replace everywhere unto by std_out
170 
171  untc=11
172  unts=12
173  untad=19
174  untd=17
175  untl=18
176  unta=15
177  untp=13
178  untg=20
179 
180 !OPENING OF THE INPUT FILES
181 
182 !DEBUG
183 !write(std_out,*) ' me,master=',me,master
184 !call flush(std_out)
185 !call xmpi_end()
186 !stop  ! OK until here
187 !ENDDEBUG
188 
189  aim_iomode = IO_MODE_FORTRAN
190  if(me==master)then
191    if (open_file(infile,msg,unit=unt0,status='old',form='formatted') /= 0) then
192      MSG_ERROR(msg)
193    end if
194 
195    if (file_exists(dnfile)) then
196      if (open_file(dnfile,msg,unit=untad,status='old',form='unformatted') /=0) then
197        MSG_ERROR(msg)
198      end if
199    else
200      if (file_exists(nctk_ncify(dnfile))) then
201        ! Use netcdf-io
202        write(std_out,"(3a)")"- File: ",trim(dnfile)," does not exist but found netcdf file with similar name."
203        dnfile = nctk_ncify(dnfile)
204        aim_iomode = IO_MODE_ETSF
205 #ifdef HAVE_NETCDF
206        NCF_CHECK(nctk_open_read(untad, dnfile, xmpi_comm_self))
207 #else
208        MSG_ERROR("Cannot read netcdf file because netcdf support in Abinit is missing.")
209 #endif
210      else
211        MSG_ERROR('Missing data file: '//TRIM(dnfile))
212      end if
213    end if
214 
215    do ii=1,nfcfile
216      iunt=unt+ii
217      if (open_file(fcfile(ii),msg,unit=iunt,status='old',form='formatted') /= 0) then
218        MSG_ERROR(msg)
219      end if
220    end do
221  end if
222  call xmpi_bcast(aim_iomode, master, comm, ierr)
223 
224 !call dump_config(std_out)
225 
226 !READING OF THE MAIN INPUT FILE
227 
228 !DEBUG
229 !write(std_out,*) ' me,master=',me,master
230 !call flush(std_out)
231 !call xmpi_end()
232 !stop ! OK until here
233 !ENDDEBUG
234 
235 
236 !Setting the input variables to their default values
237  call defad(aim_dtset)
238 
239 !Reading of the input file -> one string called instr
240  if(me==master)then
241    call inpar(instr,lenstr)
242  end if
243  call xmpi_bcast (lenstr, master, comm, ierr)
244  call xmpi_bcast (instr(1:lenstr), master, comm, ierr)
245 
246 !DEBUG
247 !write(std_out,*) ' after inpar, infile= ',infile
248 !write(std_out,*) ' after inpar, dnfile= ',dnfile
249 !write(std_out,*) ' after inpar, ofile= ',ofile
250 !do ii=1,nfcfile
251 !write(std_out,*) ' after inpar, fcfile(',ii,')= ',fcfile(ii)
252 !enddo
253 !write(std_out,*) ' after inpar, instr= ',instr(1:lenstr)
254 !write(std_out,*) ' after inpar, lenstr= ',lenstr
255 !call flush(std_out)
256 !call xmpi_end()
257 !stop   ! OK until here
258 !ENDDEBUG
259 
260 !Analysis of the input string, setting of input variables in aim_dtset
261  call adini(aim_dtset,instr,lenstr)
262 
263 !OPENING OF THE OUTPUT FILES
264  if(me==master)then
265 
266    if (aim_dtset%isurf/=0) hname(fin+2:fin+5)='surf'
267 
268    if (aim_dtset%isurf==1) then
269      ierr = open_file(hname(1:fin+5),msg,unit=unts,status='unknown',form='formatted')
270    elseif (aim_dtset%isurf==-1) then
271      ierr = open_file(hname(1:fin+5),msg,unit=unts,status='old',action='read',form='formatted')
272    end if
273    if (ierr /= 0) then
274      MSG_ERROR(msg)
275    end if
276 
277    if (aim_dtset%crit/=0) hname(fin+2:fin+5)='crit'
278 
279    if (aim_dtset%crit>0) then
280      ierr = open_file(hname(1:fin+5),msg,unit=untc,status='unknown',form='formatted')
281    else if (aim_dtset%crit==-1) then
282      ierr = open_file(hname(1:fin+5),msg,unit=untc,status='old',action='read',form='formatted')
283    end if
284    if (ierr /= 0) then
285      MSG_ERROR(msg)
286    end if
287 
288    if (aim_dtset%denout==1) then
289      hname(fin+2:fin+4)='dn1'
290      ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted')
291    elseif (aim_dtset%denout==2) then
292      hname(fin+2:fin+4)='dn2'
293      ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted')
294    elseif (aim_dtset%denout==3) then
295      hname(fin+2:fin+4)='dn3'
296      ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='formatted')
297    elseif (aim_dtset%denout==-1) then
298      hname(fin+2:fin+4)='dna'
299      ierr = open_file(hname(1:fin+4),msg,unit=untd,status='unknown',form='unformatted')
300    end if
301 
302    if (ierr /= 0) then
303      MSG_ERROR(msg)
304    end if
305 
306    if (aim_dtset%lapout==1) then
307      hname(fin+2:fin+4)='lp1'
308      ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted')
309    elseif (aim_dtset%lapout==2) then
310      hname(fin+2:fin+4)='lp2'
311      ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted')
312    elseif (aim_dtset%lapout==3) then
313      hname(fin+2:fin+4)='lp3'
314      ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='formatted')
315    elseif (aim_dtset%lapout==-1) then
316      hname(fin+2:fin+4)='lpa'
317      ierr = open_file(hname(1:fin+4),msg,unit=untl,status='unknown',form='unformatted')
318    end if
319 
320    if (ierr /= 0) then
321      MSG_ERROR(msg)
322    end if
323 
324    if (aim_dtset%gpsurf==1) then
325      hname(fin+2:fin+3)='gp'
326      if (open_file(hname(1:fin+3),msg,unit=untg,status='unknown',form='formatted') /= 0) then
327        MSG_ERROR(msg)
328      end if
329    end if
330 
331    if (aim_dtset%plden==1) then
332      hname(fin+2:fin+4)='pld'
333      if (open_file(hname(1:fin+5),msg,unit=untp,status='unknown',form='formatted') /= 0) then
334        MSG_ERROR(msg)
335      end if
336    end if
337 
338  end if
339 
340  call timein(tcpu,twall)
341  tsec(1)=tcpu-tcpui
342  tsec(2)=twall-twalli
343  write(std_out, '(5a,f13.1,a,f13.1)' ) &
344 & '-',ch10,'- After reading the input file and opening the output files ',ch10,&
345 & '- Proc.   0 individual time (sec): cpu=',tsec(1),'  wall=',tsec(2)
346 
347 !MAIN DRIVER OF THE ANALYSIS
348 
349  write(std_out,'(a,a,a,a,i4)' )char(10),&
350 & ' aim : read density file ',trim(dnfile),' from unit number ',untad
351 
352  call drvaim(aim_dtset,tcpui,twalli)
353 
354 !THE TOTAL TIME NEEDED
355  call timein(tcpu,twall)
356  tsec(1)=tcpu-tcpui
357  tsec(2)=twall-twalli
358  write(std_out, '(a,a,a,f13.1,a,f13.1)' ) &
359 & '-',ch10,'- Proc.   0 individual time (sec): cpu=',tsec(1),'  wall=',tsec(2)
360 
361  write(std_out,'(/," Time needed (seconds) - total, CP analyse, SURF determination:",/,/,"-         ",3F16.8)')&
362 & tsec(2),ttcp,ttsrf
363 
364 
365  if(me==0)then
366    write(untout,*)
367    write(untout,*) "TIME ANALYSIS"
368    write(untout,*) "============"
369    write(untout,'(/," Time needed (seconds) - total, CP analyse, SURF determination:",/,/,"-         ",3F16.8)') &
370 &   tsec(1),ttcp,ttsrf
371    write(untout,'(a,a,f11.3,a,f11.3,a,a,a,a)') char(10),&
372 &   '+Total cpu time',tsec(1),&
373 &   '  and wall time',tsec(2),' sec',char(10),char(10),&
374 &   ' aim : the run completed succesfully.'
375  end if
376 
377 !CLOSING OF THE FILES
378  if(me==0)then
379    close(unt0)
380    close(untout)
381    if (aim_dtset%isurf/=0) close(unts)
382    if (aim_dtset%crit/=0) close(untc)
383    if (aim_dtset%denout/=0) close(untd)
384    if (aim_dtset%lapout/=0) close(untl)
385    if (aim_dtset%gpsurf/=0) close(untg)
386    if (aim_dtset%plden/=0) close(untp)
387  end if
388 
389  call aim_shutdown()
390 
391  !call abinit_doctor("__aim")
392 
393 !Eventual cleaning of MPI run
394  call xmpi_end()
395 
396  end program aim