TABLE OF CONTENTS


ABINIT/gstateimg [ Functions ]

[ Top ] [ Functions ]

NAME

 gstateimg

FUNCTION

 Routine for conducting DFT calculations for a set of (dynamical) images

INPUTS

  codvsn=code version
  cpui=initial CPU time
  nimage=number of images of the cell (treated by current proc)
  === Optional arguments (needed when nimage>1) ===
    filnam(5)=character strings giving file names
    filstat=character strings giving name of status file
    idtset=index of the dataset
    jdtset(0:ndtset)=actual index of the datasets
    ndtset=number of datasets

OUTPUT

  etotal_img=total energy, for each image
  fcart_img(3,natom,nimage)=forces, in cartesian coordinates, for each image
  fred_img(3,natom,nimage)=forces, in reduced coordinates, for each image
  npwtot(nkpt) = total number of plane waves at each k point
  strten_img(6,nimage)=stress tensor, for each image

SIDE EFFECTS

  acell_img(3,nimage)=unit cell length scales (bohr), for each image
  amu_img(ntypat,nimage)=value of mass for each atomic type, for each image
  dtfil <type(datafiles_type)>=variables related to files
  dtset <type(dataset_type)>=all input variables in this dataset
   | mband =maximum number of bands (IN)
   | mgfft =maximum single fft dimension (IN)
   | mkmem =maximum number of k points which can fit in core memory (IN)
   | mpw   =maximum number of planewaves in basis sphere (large number) (IN)
   | natom =number of atoms in unit cell (IN)
   | nfft  =(effective) number of FFT grid points (for this processor) (IN)
   | nkpt  =number of k points (IN)
   | nspden=number of spin-density components (IN)
   | nsppol=number of channels for spin-polarization (1 or 2) (IN)
   | nsym  =number of symmetry elements in space group
  iexit= exit flag
  mixalch_img(npspalch,ntypalch,nimage)=value of alchemical mixing factors,for each image
  mpi_enreg=MPI-parallelisation information (some already initialized,
            some others to be initialized here)
  occ_img(mband*nkpt*nsppol,nimage) = occupation number for each band and k, for each image
  pawang <type(pawang_type)>=paw angular mesh and related data
  pawrad(ntypat*usepaw) <type(pawrad_type)>=paw radial mesh and related data
  pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
  psps <type(pseudopotential_type)>=variables related to pseudopotentials
   Before entering the first time in gstateimg, a significant part of
   psps has been initialized :
   the integers dimekb,lmnmax,lnmax,mpssang,mpssoang,mpsso,mgrid,
     ntypat,n1xccc,usepaw,useylm, and the arrays dimensioned to npsp
   All the remaining components of psps are to be initialized in the call
   to pspini .
   The next time the code enters gstateimg, psps might be identical to the
   one of the previous dtset, in which case, no reinitialisation is scheduled
   in pspini.f .
  rprim_img(3,3,nimage)=dimensionless real space primitive translations, for each image
  vel_cell_img(3,3,nimage)=value of cell parameters velocities, for each image
  vel_img(3,natom,nimage)=value of atomic velocities,for each image
  xred_img(3,natom,nimage) = reduced atomic coordinates, for each image

NOTES

 USE OF FFT GRIDS:
 =================
 In case of PAW:
 ---------------
    Two FFT grids are used:
    - A "coarse" FFT grid (defined by ecut)
      for the application of the Hamiltonian on the plane waves basis.
      It is defined by nfft, ngfft, mgfft, ...
      Hamiltonian, wave-functions, density related to WFs (rhor here), ...
      are expressed on this grid.
    - A "fine" FFT grid (defined) by ecutdg)
      for the computation of the density inside PAW spheres.
      It is defined by nfftf, ngfftf, mgfftf, ...
      Total density, potentials, ...
      are expressed on this grid.
 In case of norm-conserving:
 ---------------------------
    - Only the usual FFT grid (defined by ecut) is used.
      It is defined by nfft, ngfft, mgfft, ...
      For compatibility reasons, (nfftf,ngfftf,mgfftf)
      are set equal to (nfft,ngfft,mgfft) in that case.
 In case of wavelets:
 --------------------
    - Only the usual FFT grid (defined by wvl_crmult) is used.
      It is defined by nfft, ngfft, mgfft, ... This is strictly not
      an FFT grid since its dimensions are not suited for FFTs. They are
      defined by wvl_setngfft().
      For compatibility reasons, (nfftf,ngfftf,mgfftf)
      are set equal to (nfft,ngfft,mgfft) in that case.

TODO

 Not yet possible to use restartxf in parallel when localrdwf==0

PARENTS

      driver,gwls_sternheimer

CHILDREN

      abihist_bcast,abihist_copy,abihist_free,abihist_init,args_gs_free
      args_gs_init,copy_results_img,destroy_results_img,dtfil_init,fcart2fred
      ga_destroy,ga_init,gstate,init_results_img,libpaw_spmsg_mpisum
      localfilnam,localrdfile,localredirect,localwrfile,mep_destroy,mep_init
      mkradim,mkrdim,pimd_destroy,pimd_init,pimd_skip_qtb,predictimg,prtimg
      read_md_hist_img,scf_history_free,scf_history_nullify,specialmsg_mpisum
      status,timab,var2hist,vel2hist,write_md_hist_img,wrtout,xmpi_barrier
      xmpi_sum

SOURCE

189 subroutine gstateimg(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
190 &                    fred_img,iexit,mixalch_img,mpi_enreg,nimage,npwtot,occ_img,&
191 &                    pawang,pawrad,pawtab,psps,&
192 &                    rprim_img,strten_img,vel_cell_img,vel_img,wvl,xred_img,&
193 &                    filnam,filstat,idtset,jdtset,ndtset) ! optional arguments
194 
195 
196 !This section has been created automatically by the script Abilint (TD).
197 !Do not modify the following lines by hand.
198 #undef ABI_FUNC
199 #define ABI_FUNC 'gstateimg'
200 !End of the abilint section
201 
202  implicit none
203 
204 !Arguments ------------------------------------
205 !scalars
206  integer,intent(in) :: nimage
207  integer,optional,intent(in) :: idtset,ndtset
208  integer,intent(inout) :: iexit
209  real(dp),intent(in) :: cpui
210  character(len=6),intent(in) :: codvsn
211  character(len=fnlen),optional,intent(in) :: filstat
212  type(MPI_type),intent(inout) :: mpi_enreg
213  type(datafiles_type),target,intent(inout) :: dtfil
214  type(dataset_type),target,intent(inout) :: dtset
215  type(pawang_type),intent(inout) :: pawang
216  type(pseudopotential_type),intent(inout) :: psps
217  type(wvl_data),intent(inout) :: wvl
218 !arrays
219  integer,optional,intent(in) :: jdtset(:)
220  integer,intent(out) :: npwtot(dtset%nkpt)
221  character(len=fnlen),optional,intent(in) :: filnam(:)
222  real(dp), intent(out) :: etotal_img(nimage),fcart_img(3,dtset%natom,nimage)
223  real(dp), intent(out) :: fred_img(3,dtset%natom,nimage),strten_img(6,nimage)
224  real(dp),intent(inout) :: acell_img(3,nimage),amu_img(dtset%ntypat,nimage)
225  real(dp),intent(inout) :: mixalch_img(dtset%npspalch,dtset%ntypalch,nimage)
226  real(dp),intent(inout) :: occ_img(dtset%mband*dtset%nkpt*dtset%nsppol,nimage)
227  real(dp),intent(inout) :: rprim_img(3,3,nimage),vel_cell_img(3,3,nimage),vel_img(3,dtset%natom,nimage)
228  real(dp),intent(inout) :: xred_img(3,dtset%natom,nimage)
229  type(pawrad_type),intent(inout) :: pawrad(psps%ntypat*psps%usepaw)
230  type(pawtab_type),intent(inout) :: pawtab(psps%ntypat*psps%usepaw)
231 
232 !Local variables-------------------------------
233 !Define file format for different type of files. Presently,
234 !only one file format is supported for each type of files, but this might
235 !change soon ...
236 !2   for wavefunction file, new format (version 2.0 and after)    (fform)   NOT USED
237 !52  for density rho(r)       (fformr)
238 !102 for potential V(r) file. (fformv)  NOT USED
239 !scalars
240  integer,parameter :: formeig=0,level=100,ndtpawuj=0,response=0
241  integer :: history_size,idelta,idynimage,ierr,ifirst
242  integer :: ii,iimage,ih,itimimage,itimimage_eff,itimimage_prev,ndynimage,nocc
243  integer :: ntimimage,ntimimage_stored,ntimimage_max
244  logical :: check_conv,compute_all_images,compute_static_images
245  logical :: isVused,isARused,is_master,is_mep,is_pimd
246  logical :: call_predictor,use_hist,use_hist_prev
247  real(dp) :: delta_energy,dtion
248  character(len=500) :: hist_filename,msg
249  type(args_gs_type) :: args_gs
250  type(mep_type) :: mep_param
251  type(ga_type) :: ga_param
252  type(m1geo_type) :: m1geo_param
253  type(pimd_type) :: pimd_param
254 !arrays
255  integer,allocatable :: list_dynimage(:),scf_initialized(:)
256  character(len=60),parameter :: imagealgo_str(0:13)=(/ &
257 &   'IMAGE COPY                                                  ',& ! 0
258 &   'IMAGE STEEPEST DESCENT                                      ',& ! 1
259 &   'STRING METHOD                                               ',& ! 2
260 &   'METADYNAMICS                                                ',& ! 3
261 &   'GENETIC ALGORITHM                                           ',& ! 4
262 &   'NUDGED ELASTIC BAND                                         ',& ! 5
263 &   'LINEAR COMBINATION OF CONSTRAINED DFT ENERGIES              ',& ! 6
264 &   '                                                            ',& ! 7
265 &   '                                                            ',& ! 8
266 &   'PATH-INTEGRAL MOLECULAR DYNAMICS (LANGEVIN)                 ',& ! 9
267 &   'PATH-INTEGRAL MOLECULAR DYNAMICS (QUANTUM THERMAL BATH)     ',& ! 10
268 &   '                                                            ',& ! 11
269 &   '                                                            ',& ! 12
270 &   'PATH-INTEGRAL MOLECULAR DYNAMICS (CHAIN OF THERMOSTATS)     '/) ! 13
271  real(dp) :: acell(3),rprim(3,3),rprimd(3,3),tsec(2),vel_cell(3,3)
272  real(dp),allocatable :: amass(:,:),occ(:),vel(:,:),xred(:,:)
273  type(abihist),allocatable :: hist(:),hist_prev(:)
274  type(results_img_type),pointer :: results_img(:,:),res_img(:)
275  type(scf_history_type),allocatable :: scf_history(:)
276  type(abiforstr) :: preconforstr ! Preconditioned forces and stress ... Only needed to deallocate an internal matrix in prec_simple
277 
278 ! ***********************************************************************
279 
280  DBG_ENTER("COLL")
281 
282  call timab(700,1,tsec)
283  call timab(703,3,tsec)
284 
285  call status(0,dtfil%filstat,iexit,level,'enter         ')
286 
287 !Arguments check
288  if (dtset%nimage>1) then
289    if ((.not.present(filnam)).or.(.not.present(filnam)).or.(.not.present(idtset)).or.&
290 &   (.not.present(ndtset)).or.(.not.present(jdtset))) then
291      write(msg,'(3a)') &
292 &     'When nimage>1, all the following argument should be present:',ch10,&
293 &     'filnam, filstat, idtset, ndtset, jdtset  !'
294      MSG_BUG(msg)
295    end if
296  end if
297 
298 !Set flag for the effective computation (once) of static images
299 !For the time being only set when parallelization is activated
300 !Note: if you modify this flag, do not forget to change it in outvars and outvar1
301  compute_static_images=(dtset%istatimg>0)
302 
303 !Prepare the calculation, by computing flags and dimensions
304  is_pimd=(dtset%imgmov==9.or.dtset%imgmov==10.or.dtset%imgmov==13)
305  is_mep =(dtset%imgmov==1.or.dtset%imgmov== 2.or.dtset%imgmov== 5)
306  ntimimage=dtset%ntimimage
307  ntimimage_stored=ntimimage;if(is_pimd)ntimimage_stored=2
308  nocc=dtset%mband*dtset%nkpt*dtset%nsppol
309  is_master=(mpi_enreg%me_cell==0.and.mpi_enreg%me_img==0)
310  delta_energy=zero
311 
312 !Management of dynamics/relaxation history (positions, forces, stresses, ...)
313  use_hist=(dtset%imgmov/=0.and.nimage>0) ; use_hist_prev=.false.
314  isVused=is_pimd;isARused=(dtset%optcell/=0)
315  if (use_hist) then
316    !Read history from file (and broadcast if MPI)
317 #if defined HAVE_NETCDF
318    use_hist_prev=(dtset%restartxf==-1.and.nimage>0)
319 #endif
320    hist_filename=trim(dtfil%filnam_ds(4))//'_HIST.nc'
321    if (use_hist_prev)then
322      ABI_DATATYPE_ALLOCATE(hist_prev,(nimage))
323      if (mpi_enreg%me_cell==0) then
324        call read_md_hist_img(hist_filename,hist_prev,isVused,isARused,&
325 &       imgtab=mpi_enreg%my_imgtab)
326      end if
327      call abihist_bcast(hist_prev,0,mpi_enreg%comm_cell)
328      if (nimage>0) then
329        if (any(hist_prev(:)%mxhist/=hist_prev(1)%mxhist)) then
330          msg='History problem: all images should have the same number of time steps!'
331          MSG_ERROR(msg)
332        end if
333        use_hist_prev=(hist_prev(1)%mxhist>0)
334        if (use_hist_prev) ntimimage=ntimimage+hist_prev(1)%mxhist
335      end if
336      if (.not.use_hist_prev) then
337        call abihist_free(hist_prev)
338        ABI_DATATYPE_DEALLOCATE(hist_prev)
339      end if
340    end if
341    !Initialize a variable to write the history
342    ABI_DATATYPE_ALLOCATE(hist,(nimage))
343    call abihist_init(hist,dtset%natom,ntimimage,isVused,isARused)
344  end if ! imgmov/=0
345 
346 !Allocations
347  ABI_ALLOCATE(occ,(nocc))
348  ABI_ALLOCATE(vel,(3,dtset%natom))
349  ABI_ALLOCATE(xred,(3,dtset%natom))
350  ABI_DATATYPE_ALLOCATE(results_img,(nimage,ntimimage_stored))
351  ABI_ALLOCATE(list_dynimage,(dtset%ndynimage))
352  do itimimage=1,ntimimage_stored
353    res_img => results_img(:,itimimage)
354    call init_results_img(dtset%natom,dtset%npspalch,dtset%nsppol,dtset%ntypalch,&
355 &   dtset%ntypat,res_img)
356    do iimage=1,nimage
357      res_img(iimage)%acell(:)     =acell_img(:,iimage)
358      res_img(iimage)%amu(:)       =amu_img(:,iimage)
359      res_img(iimage)%mixalch(:,:) =mixalch_img(:,:,iimage)
360      res_img(iimage)%rprim(:,:)   =rprim_img(:,:,iimage)
361      res_img(iimage)%xred(:,:)    =xred_img(:,:,iimage)
362      res_img(iimage)%vel(:,:)     =vel_img(:,:,iimage)
363      res_img(iimage)%vel_cell(:,:)=vel_cell_img(:,:,iimage)
364    end do
365  end do
366  ndynimage=0
367  do iimage=1,nimage
368    ii=mpi_enreg%my_imgtab(iimage)
369    if (dtset%dynimage(ii)==1) then
370      ndynimage=ndynimage+1
371      list_dynimage(ndynimage)=iimage
372    end if
373  end do
374 
375 !Management of SCF history (density/WF predictions from one time step to another)
376  ABI_DATATYPE_ALLOCATE(scf_history,(nimage))
377  ABI_ALLOCATE(scf_initialized,(nimage))
378  scf_initialized=0
379  history_size=-1
380  if (dtset%ntimimage<=1) then
381    if (dtset%usewvl==0.and.dtset%ionmov>0.and. &
382 &   (abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) history_size=2
383  else
384    if (abs(dtset%densfor_pred)==2.or.abs(dtset%densfor_pred)==3) history_size=0
385    if (dtset%imgwfstor==1) history_size=1
386    if (dtset%usewvl==0.and.(abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) history_size=2
387  end if
388  do iimage=1,nimage
389    call scf_history_nullify(scf_history(iimage))
390    scf_history(iimage)%history_size=history_size
391  end do
392 
393 !In some cases, need amass variable
394  if (use_hist) then
395    ABI_ALLOCATE(amass,(dtset%natom,nimage))
396    do iimage=1,nimage
397      if (any(amu_img(:,iimage)/=amu_img(:,1))) then
398        msg='HIST file is not compatible with variable masses!'
399        MSG_ERROR(msg)
400      end if
401      amass(:,iimage)=amu_emass*amu_img(dtset%typat(:),iimage)
402    end do
403  end if
404 
405 !In the case of the 4th-order Runge-Kutta solver,
406 !one must have a number of step multiple of 4.
407  ntimimage_max=ntimimage;idelta=1
408  if (dtset%imgmov==2.and.mep_param%mep_solver==4) then
409    ntimimage_max=4*(ntimimage_max/4)
410    idelta=4
411  end if
412 
413 !MEP search: fill in eventually the data structure mep_param
414  call mep_init(dtset,mep_param)
415 
416 !GA search: fill in eventually the data structure ga_param
417  call ga_init(dtset,ga_param)
418 
419 !Move 1GEO approach: fill the data structure m1geo_param
420  call m1geo_init(dtfil,dtset,m1geo_param)
421 
422 !PIMD: fill in eventually the data structure pimd_param
423  call pimd_init(dtset,pimd_param,is_master)
424  dtion=one;if (is_pimd) dtion=pimd_param%dtion
425 
426  call timab(703,2,tsec)
427 
428 !-----------------------------------------------------------------------------------------
429 !Big loop on the propagation of all images
430  itimimage_eff=1
431  do itimimage=1,ntimimage
432 
433    res_img => results_img(:,itimimage_eff)
434    call_predictor=(ntimimage>1)
435 
436 !  If history is activated and if current image is inside it: do not compute anything
437    if (use_hist_prev) then
438      if (all(hist_prev(:)%ihist<=hist_prev(:)%mxhist)) then
439        do iimage=1,nimage
440          ih=hist_prev(iimage)%ihist
441          call abihist_copy(hist_prev(iimage),hist(iimage))
442          call mkradim(hist_prev(iimage)%acell(:,ih),rprim,hist_prev(iimage)%rprimd(:,:,ih))
443          res_img(iimage)%acell(:)=hist_prev(iimage)%acell(:,ih)
444          res_img(iimage)%rprim(:,:)=rprim
445          res_img(iimage)%xred(:,:)=hist_prev(iimage)%xred(:,:,ih)
446          res_img(iimage)%vel(:,:)=hist_prev(iimage)%vel(:,:,ih)
447          res_img(iimage)%vel_cell(:,:)=hist_prev(iimage)%vel_cell(:,:,ih)
448          res_img(iimage)%results_gs%fcart(:,:)=hist_prev(iimage)%fcart(:,:,ih)
449          res_img(iimage)%results_gs%strten(:)=hist_prev(iimage)%strten(:,ih)
450          res_img(iimage)%results_gs%etotal=hist_prev(iimage)%etot(ih)
451          res_img(iimage)%results_gs%energies%entropy=hist_prev(iimage)%entropy(ih)
452          call fcart2fred(res_img(iimage)%results_gs%fcart,res_img(iimage)%results_gs%fred,&
453 &         hist_prev(iimage)%rprimd(:,:,ih),dtset%natom)
454          hist_prev(iimage)%ihist=hist_prev(iimage)%ihist+1
455        end do
456        !PI-QTB: skip a record in random force file
457        if (pimd_param%use_qtb==1) call pimd_skip_qtb(pimd_param)
458        !call_predictor=.false.
459        goto 110 ! This is temporary
460      end if
461    end if
462 
463    call timab(704,1,tsec)
464    call localfilnam(mpi_enreg%comm_img,mpi_enreg%comm_cell,mpi_enreg%comm_world,filnam,'_IMG',dtset%nimage)
465    compute_all_images=(compute_static_images.and.itimimage==1)
466 
467 !  Print title for time step
468    if (dtset%nimage>1.or.dtset%ntimimage>1) then
469      if (dtset%prtvolimg<2) then
470        msg=ch10;if (itimimage >1) write(msg,'(2a)') ch10,ch10
471        write(msg,'(5a)') trim(msg),&
472 &       '================================================================================',&
473 &       ch10,' ',trim(imagealgo_str(dtset%imgmov))
474      else
475        msg='';if (itimimage >1) msg=ch10
476        write(msg,'(5a)') trim(msg),&
477 &       '--------------------------------------------------------------------------------',&
478 &       ch10,' ',trim(imagealgo_str(dtset%imgmov))
479      end if
480      if (dtset%ntimimage==1) write(msg,'(2a)')    trim(msg),' FOR 1 TIME STEP'
481      if (dtset%ntimimage >1) write(msg,'(2a,i5)') trim(msg),' - TIME STEP ',itimimage
482      if (dtset%prtvolimg<2) then
483        write(msg,'(3a)') trim(msg),ch10,&
484 &       '================================================================================'
485      end if
486      call wrtout(ab_out ,msg,'COLL')
487      call wrtout(std_out,msg,'PERS')
488    end if
489 
490    call timab(704,2,tsec)
491 
492 !  Loop on the dynamical images
493    idynimage=0
494    do iimage=1,nimage
495 
496      call timab(705,1,tsec)
497 
498      ii=mpi_enreg%my_imgtab(iimage)
499      if (dtset%dynimage(ii)==1) idynimage=idynimage+1
500 
501 !    Compute static image only at first time step
502      if (dtset%dynimage(ii)==1.or.compute_all_images) then
503 
504 !      Change file names according to image index (if nimage>1)
505        if (dtset%nimage>1) then
506          call dtfil_init(dtfil,dtset,filnam,filstat,idtset,jdtset,mpi_enreg,ndtset,&
507 &         image_index=ii)
508          if (itimimage>1) then
509            dtfil%ireadwf=0;dtfil%ireadden=0;dtfil%ireadkden=0
510          end if
511        end if
512 
513 !      Redefine output units
514        call localwrfile(mpi_enreg%comm_cell,ii,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg)
515 
516 !      Print title for image
517        if (dtset%nimage>1.and.(dtset%prtvolimg==0.or.do_write_log)) then
518          if (ii==1) write(msg,'(a)' ) ch10
519          if (ii >1) write(msg,'(2a)') ch10,ch10
520          write(msg,'(6a,i4,a,i4,3a)') trim(msg),&
521 &         '--------------------------------------------------------------------------------',ch10,&
522 &         ' ',trim(imagealgo_str(dtset%imgmov)),' - CELL # ',ii,'/',dtset%nimage,ch10,&
523 &         '--------------------------------------------------------------------------------',ch10
524          if (dtset%prtvolimg==0) then
525            call wrtout(ab_out ,msg,'COLL')
526          end if
527          if (do_write_log) then
528            call wrtout(std_out,msg,'PERS')
529          end if
530        end if
531 
532        acell(:)     =res_img(iimage)%acell(:)
533        rprim(:,:)   =res_img(iimage)%rprim(:,:)
534        vel(:,:)     =res_img(iimage)%vel(:,:)
535        vel_cell(:,:)=res_img(iimage)%vel_cell(:,:)
536        xred(:,:)    =res_img(iimage)%xred(:,:)
537        occ(:)       =occ_img(:,iimage)
538 
539        call args_gs_init(args_gs, &
540 &       res_img(iimage)%amu(:),res_img(iimage)%mixalch(:,:),&
541 &       dtset%dmatpawu(:,:,:,:,ii),dtset%upawu(:,ii),dtset%jpawu(:,ii),&
542 &       dtset%rprimd_orig(:,:,ii))
543 
544        call timab(705,2,tsec)
545 
546        call status(idynimage+100*itimimage,dtfil%filstat,iexit,level,'call gstate   ')
547        call gstate(args_gs,acell,codvsn,cpui,dtfil,dtset,iexit,scf_initialized(iimage),&
548 &       mpi_enreg,npwtot,occ,pawang,pawrad,pawtab,psps,&
549 &       res_img(iimage)%results_gs,&
550 &       rprim,scf_history(iimage),vel,vel_cell,wvl,xred)
551 
552        call timab(706,1,tsec)
553 
554        call args_gs_free(args_gs)
555 
556        if (dtset%dynimage(ii)==1) then
557          res_img(iimage)%acell(:)     =acell(:)
558          res_img(iimage)%rprim(:,:)   =rprim(:,:)
559          res_img(iimage)%vel(:,:)     =vel(:,:)
560          res_img(iimage)%vel_cell(:,:)=vel_cell(:,:)
561          res_img(iimage)%xred(:,:)    =xred(:,:)
562          occ_img(:,iimage)            =occ(:)
563        end if
564 
565 !    check change of rprim and rewriting in hist
566 !    check change of xred and rewriting in hist
567 
568 !      Close output units ; restore defaults
569        call localredirect(mpi_enreg%comm_cell,mpi_enreg%comm_world,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg)
570        call timab(706,2,tsec)
571 
572      else if (itimimage>1) then ! For static images, simply copy one time step to the other
573        itimimage_prev=itimimage_eff-1
574        if (itimimage_prev<1) itimimage_prev=ntimimage_stored
575        call copy_results_img(results_img(iimage,itimimage_prev), &
576 &       results_img(iimage,itimimage_eff ))
577      end if
578 
579 !    Store results in hist datastructure
580      if (use_hist) then
581        ih=hist(iimage)%ihist
582        call mkrdim(res_img(iimage)%acell(:),res_img(iimage)%rprim(:,:),rprimd)
583        call var2hist(res_img(iimage)%acell(:),hist(iimage),dtset%natom,&
584 &       rprimd,res_img(iimage)%xred(:,:),.FALSE.)
585        call vel2hist(amass(:,iimage),hist(iimage),res_img(iimage)%vel(:,:),&
586 &       res_img(iimage)%vel_cell(:,:))
587        hist(iimage)%fcart(:,:,ih)=res_img(iimage)%results_gs%fcart(:,:)
588        hist(iimage)%strten(:,ih)=res_img(iimage)%results_gs%strten(:)
589        hist(iimage)%etot(ih)=res_img(iimage)%results_gs%etotal
590        hist(iimage)%entropy(ih)=res_img(iimage)%results_gs%energies%entropy
591        hist(iimage)%time(ih)=real(itimimage,kind=dp)*dtion
592      end if
593 
594    end do ! iimage
595 
596    if(mpi_enreg%paral_img==1)then
597      call timab(702,1,tsec)
598      call xmpi_barrier(mpi_enreg%comm_img)
599      call timab(702,2,tsec)
600    end if
601 
602    call timab(707,1,tsec)
603 
604 !  Output when images are used
605    if (dtset%nimage>1) then
606 !    === 1st option: reduced outputs ===
607      if (dtset%prtvolimg>0) then
608        call prtimg(dtset%dynimage,imagealgo_str(dtset%imgmov),dtset%imgmov,ab_out,&
609 &       mpi_enreg,nimage,dtset%nimage,compute_all_images,dtset%prtvolimg,res_img)
610      end if
611    end if
612 
613 !  Manage log files when images are used
614    call localrdfile(mpi_enreg%comm_img,mpi_enreg%comm_world,compute_all_images,&
615 &   dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg,dyn=dtset%dynimage)
616 
617 !  Write hist datastructure in HIST file
618 #if defined HAVE_NETCDF
619    if (use_hist.and.mpi_enreg%me_cell==0) then
620      ifirst=merge(0,1,itimimage>1)
621      call write_md_hist_img(hist,hist_filename,ifirst,itimimage,dtset%natom,dtset%ntypat,&
622 &     dtset%typat,amu_img(:,1),dtset%znucl,dtion,&
623 &     nimage=dtset%nimage,imgmov=dtset%imgmov,mdtemp=dtset%mdtemp,comm_img=mpi_enreg%comm_img,&
624 &     imgtab=mpi_enreg%my_imgtab)
625    end if
626 #endif
627 
628 !  TESTS WHETHER ONE CONTINUES THE LOOP
629 !  Here we calculate the change in energy, and exit if delta_energy < tolimg
630    delta_energy=zero
631 !  Doesn't check convergence in case of PIMD
632    check_conv=((.not.is_pimd).and.itimimage>1)
633 !  In case of 4th-order Runge-Kutta, does check convergence every 4 steps
634    if (dtset%imgmov==2.and.mep_param%mep_solver==4) then
635      check_conv=(mod(itimimage,4)==0.and.itimimage>4)
636    end if
637    if (check_conv) then
638      do idynimage=1,ndynimage
639        iimage=list_dynimage(idynimage)
640        delta_energy=delta_energy &
641 &       +abs(results_img(iimage,itimimage)%results_gs%etotal &
642 &       -results_img(iimage,itimimage-idelta)%results_gs%etotal)
643      end do
644      if (mpi_enreg%paral_img==1) then
645        call xmpi_sum(delta_energy,mpi_enreg%comm_img,ierr)
646      end if
647      delta_energy=delta_energy/dtset%ndynimage
648      if (delta_energy<dtset%tolimg) then
649        if (dtset%prtvolimg<2) then
650          write(msg,'(5a,i5,6a,es11.3,a,es11.3,2a)') ch10,ch10,&
651 &         '================================================================================',ch10,&
652 &         ' At time step ',itimimage,ch10,&
653 &         ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,&
654 &         ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg,ch10,&
655 &         '================================================================================'
656        else
657          write(msg,'(4a,i5,6a,es11.3,a,es11.3)') ch10,&
658 &         '--------------------------------------------------------------------------------',ch10,&
659 &         ' At time step ',itimimage,ch10,&
660 &         ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,&
661 &         ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg
662        end if
663        call wrtout(ab_out ,msg,'COLL')
664        call wrtout(std_out,msg,'COLL')
665        call timab(707,2,tsec)
666        exit   ! exit itimimage
667      end if
668    end if
669 
670 !Temporary statement
671    110 continue
672 
673 !  Dont call the predictor at last time step
674    if (itimimage>=ntimimage_max) call_predictor=(call_predictor.and.is_pimd)
675 
676 !  Predict the next value of the images
677    if (call_predictor) then
678      call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,&
679 &     itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,&
680 &     nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img)
681    end if
682 
683 !  Increment indexes
684    if (itimimage>=ntimimage_max) exit
685    itimimage_eff=itimimage_eff+1;if (itimimage_eff>ntimimage_stored) itimimage_eff=1
686    if (use_hist) then
687      do iimage=1,nimage
688        hist(iimage)%ihist=hist(iimage)%ihist+1
689      end do
690    end if
691 
692    call timab(707,2,tsec)
693 
694  end do ! itimimage
695 !-----------------------------------------------------------------------------------------
696 
697  call timab(708,1,tsec)
698 
699 !Copy the results of the computation in the appropriate arguments of the routine
700  do iimage=1,nimage
701    ii=mpi_enreg%my_imgtab(iimage)
702    if (dtset%dynimage(ii)==1) then
703      acell_img(:,iimage)     =results_img(iimage,itimimage_eff)%acell(:)
704      amu_img(:,iimage)       =results_img(iimage,itimimage_eff)%amu(:)
705      mixalch_img(:,:,iimage) =results_img(iimage,itimimage_eff)%mixalch(:,:)
706      rprim_img(:,:,iimage)   =results_img(iimage,itimimage_eff)%rprim(:,:)
707      vel_img(:,:,iimage)     =results_img(iimage,itimimage_eff)%vel(:,:)
708      vel_cell_img(:,:,iimage)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
709      xred_img(:,:,iimage)    =results_img(iimage,itimimage_eff)%xred(:,:)
710      etotal_img(iimage)      =results_img(iimage,itimimage_eff)%results_gs%etotal
711      fcart_img(:,:,iimage)   =results_img(iimage,itimimage_eff)%results_gs%fcart(:,:)
712      fred_img(:,:,iimage)    =results_img(iimage,itimimage_eff)%results_gs%fred(:,:)
713      strten_img(:,iimage)    =results_img(iimage,itimimage_eff)%results_gs%strten(:)
714    else if (compute_static_images) then
715      etotal_img(iimage)    =results_img(iimage,1)%results_gs%etotal
716      fcart_img(:,:,iimage) =results_img(iimage,1)%results_gs%fcart(:,:)
717      fred_img(:,:,iimage)  =results_img(iimage,1)%results_gs%fred(:,:)
718      strten_img(:,iimage)  =results_img(iimage,1)%results_gs%strten(:)
719    end if
720  end do
721 
722 !When parallelizattion over images is activated, has to sum number of warnings
723 !and comments written in log file
724  if (mpi_enreg%paral_img==1.and.mpi_enreg%comm_cell==0) then
725    call specialmsg_mpisum(mpi_enreg%comm_img)
726    call libpaw_spmsg_mpisum(mpi_enreg%comm_img)
727  end if
728 
729 
730 !Final deallocations
731 
732 !This call is needed to free internal storages in different routines (prec_simple, pred_bfgs ...)
733  if(dtset%imgmov==6)then
734    m1geo_param%iexit=1
735    call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,&
736 &   itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,&
737 &   nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img)
738  endif
739 
740  ABI_DEALLOCATE(occ)
741  ABI_DEALLOCATE(vel)
742  ABI_DEALLOCATE(xred)
743  ABI_DEALLOCATE(list_dynimage)
744 
745  if (allocated(amass)) then
746    ABI_DEALLOCATE(amass)
747  end if
748 
749  do itimimage=1,ntimimage_stored
750    call destroy_results_img(results_img(:,itimimage))
751  end do
752  ABI_DATATYPE_DEALLOCATE(results_img)
753  do iimage=1,nimage
754    call scf_history_free(scf_history(iimage))
755  end do
756  ABI_DATATYPE_DEALLOCATE(scf_history)
757  ABI_DEALLOCATE(scf_initialized)
758  if (allocated(hist_prev)) then
759    call abihist_free(hist_prev)
760    ABI_DATATYPE_DEALLOCATE(hist_prev)
761  end if
762  if (allocated(hist)) then
763    call abihist_free(hist)
764    ABI_DATATYPE_DEALLOCATE(hist)
765  end if
766 
767  call mep_destroy(mep_param)
768  call ga_destroy(ga_param)
769  call m1geo_destroy(m1geo_param)
770  call pimd_destroy(pimd_param)
771 
772  call status(0,dtfil%filstat,iexit,level,'exit          ')
773 
774  call timab(708,2,tsec)
775  call timab(700,2,tsec)
776 
777  DBG_EXIT("COLL")
778 
779 end subroutine gstateimg

ABINIT/m_gstateimg [ Modules ]

[ Top ] [ Modules ]

NAME

  m_gstateimg

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (XG, AR, GG, MT)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_gstateimg
27 
28  use defs_basis
29  use defs_datatypes
30  use defs_abitypes
31  use defs_wvltypes
32  use defs_rectypes
33  use m_abicore
34  use m_abihist
35  use m_mep
36  use m_ga
37  use m_use_ga
38  use m_pimd
39  use m_xmpi
40  use m_errors
41  use m_rec
42  use m_args_gs
43  use m_results_img
44  use m_scf_history
45  use m_io_redirect
46  use m_m1geo
47  use m_abimover
48 
49  use m_time,         only : timab
50  use m_geometry,     only : mkradim, mkrdim, fcart2fred, xred2xcart, metric
51  use m_specialmsg,   only : specialmsg_mpisum
52  use m_libpaw_tools, only : libpaw_spmsg_mpisum
53  use m_pawang,       only : pawang_type
54  use m_pawrad,       only : pawrad_type
55  use m_pawtab,       only : pawtab_type
56  use m_dtfil,        only : dtfil_init, status
57  use m_gstate,       only : gstate
58  use m_predtk,       only : prtxvf
59  use m_precpred_1geo, only : precpred_1geo
60  use m_pred_simple,  only : prec_simple
61 
62 #if defined  HAVE_BIGDFT
63  use BigDFT_API, only: mpi_environment_set
64 #endif
65 
66  implicit none
67 
68  private

ABINIT/move_1geo [ Functions ]

[ Top ] [ Functions ]

NAME

 move_1geo

FUNCTION

 This subroutine uses the forces, stresses and other results obtained for several images with one, common, geometry, 
 weight them to deliver averaged forces, stresses, etc, and uses these to predict the next common geometry.
 All images must be dynamical.
 WARNING : at present, only forces are used, to change atomic positions. No change of cell geometry.
 Since this is not the PIMD algorithm, suppose ntimimage_stored=ntimimage, and itimimage=itimimage_eff.

INPUTS

 itimimage_eff=time index in the history
 nimage=number of images
 ntimimage_stored=number of time steps stored in the history
  mpi_enreg=MPI-parallelisation information
 m1geo_param=parameters for the 1geo algorithms

OUTPUT

SIDE EFFECTS

 results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
   results_img(:,:)%acell(3)
    at input, history of the values of acell for all images
    at output, the predicted values of acell for all images
   results_img(:,:)%results_gs
    at input, history of the values of energies and forces for all images
   results_img(:,:)%rprim(3,3)
    at input, history of the values of rprim for all images
    at output, the predicted values of rprim for all images
   results_img(:,:)%vel(3,natom)
    at input, history of the values of vel for all images
    at output, the predicted values of vel for all images
   results_img(:,:)%vel_cell(3,3)
    at input, history of the values of vel_cell for all images
    at output, the predicted values of vel_cell for all images
   results_img(:,:)%xred(3,natom)
    at input, history of the values of xred for all images
    at output, the predicted values of xred for all images

PARENTS

      predictimg

CHILDREN

SOURCE

1279 subroutine move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,ntimimage_stored,results_img)
1280 
1281 
1282 !This section has been created automatically by the script Abilint (TD).
1283 !Do not modify the following lines by hand.
1284 #undef ABI_FUNC
1285 #define ABI_FUNC 'move_1geo'
1286 !End of the abilint section
1287 
1288  implicit none
1289 
1290 !Arguments ------------------------------------
1291 !scalars
1292  integer,intent(in) :: itimimage_eff,nimage,ntimimage_stored
1293  type(MPI_type),intent(in) :: mpi_enreg
1294  type(m1geo_type),intent(inout) :: m1geo_param
1295 !arrays
1296  type(results_img_type),intent(inout) :: results_img(nimage,ntimimage_stored)
1297 
1298 !Local variables-------------------------------
1299 !scalars
1300  integer :: ihist,iimage,natom,next_itimimage
1301  real(dp) :: acell(3),rprim(3,3),rprimd(3,3),strten(6),vel_cell(3,3)
1302  real(dp),allocatable :: fcart(:,:),vel(:,:),xred(:,:)
1303  logical :: DEBUG=.FALSE.
1304 
1305 ! *************************************************************************
1306 
1307  natom=m1geo_param%ab_mover%natom
1308  ihist=m1geo_param%hist_1geo%ihist
1309 
1310  ABI_ALLOCATE(fcart,(3,natom))
1311  ABI_ALLOCATE(vel,(3,natom))
1312  ABI_ALLOCATE(xred,(3,natom))
1313  
1314 !Of course, assume that the geometry parameters are the same for all images, so take them from the first one.
1315  xred(:,:)    =results_img(1,itimimage_eff)%xred(:,:)
1316  acell(:)     =results_img(1,itimimage_eff)%acell(:)
1317  rprim(:,:)   =results_img(1,itimimage_eff)%rprim(:,:)
1318  vel(:,:)     =results_img(1,itimimage_eff)%vel(:,:)
1319  vel_cell(:,:)=results_img(1,itimimage_eff)%vel_cell(:,:)
1320 
1321  call mkrdim(acell,rprim,rprimd)
1322 
1323 !Fill history with the values of xred, acell and rprimd
1324  call var2hist(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG)
1325 
1326 !Fill history with velocities and ionic kinetic energy
1327  call vel2hist(m1geo_param%ab_mover%amass,m1geo_param%hist_1geo,vel,vel_cell)
1328  m1geo_param%hist_1geo%time(ihist)=zero
1329 
1330 !Compute forces and stresses for the 1geo : take the weighted average.
1331  fcart(:,:)=zero
1332  strten(:)=zero
1333  do iimage=1,nimage
1334    fcart(:,:)=fcart(:,:)+results_img(iimage,itimimage_eff)%results_gs%fcart(:,:)*m1geo_param%mixesimgf(iimage) 
1335    strten(:) =strten(:) +results_img(iimage,itimimage_eff)%results_gs%strten(:)*m1geo_param%mixesimgf(iimage) 
1336  enddo
1337 
1338 !Store them in hist_1geo 
1339  m1geo_param%hist_1geo%fcart(:,:,ihist)=fcart(:,:)
1340  m1geo_param%hist_1geo%strten(:,ihist) =strten(:)
1341 
1342 !Store them in ab_xfh
1343 !THIS IS TO BE DONE !
1344 
1345 !Compute new atomic positions and cell characteristics in the single geometry
1346  call precpred_1geo(m1geo_param%ab_mover,&
1347 & m1geo_param%ab_xfh_1geo,&
1348 & m1geo_param%ab_mover%amu_curr,&
1349 & m1geo_param%deloc,&
1350 & m1geo_param%dt_chkdilatmx,&
1351 & mpi_enreg%comm_cell,&
1352 & m1geo_param%dilatmx,&
1353 & m1geo_param%filnam_ds4,&
1354 & m1geo_param%hist_1geo,&
1355 & m1geo_param%hmctt,&
1356 & m1geo_param%icycle,&
1357 & m1geo_param%iexit,&
1358 !& m1geo_param%itime,&
1359   itimimage_eff,&       ! m1geo_param%itime should be eliminated, no need for it 
1360 & m1geo_param%mttk_vars,&
1361 & m1geo_param%nctime,&
1362 & m1geo_param%ncycle,&
1363 & m1geo_param%nerr_dilatmx,&
1364 & m1geo_param%npsp,&
1365 & m1geo_param%ntime,&
1366 & m1geo_param%rprimd_orig,&
1367 & m1geo_param%skipcycle,&
1368 & m1geo_param%usewvl)
1369 
1370 !Retrieve the new positions, cell parameters [and velocities ?!]
1371  call hist2var(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG)
1372 
1373 !Store acell, rprim, xred and vel for the new iteration
1374  next_itimimage=itimimage_eff+1
1375  if (next_itimimage>ntimimage_stored)then
1376    MSG_ERROR('next_itimimage>ntimimage_stored')
1377  endif
1378 
1379  do iimage=1,nimage
1380    results_img(iimage,next_itimimage)%xred(:,:)    =xred(:,:)
1381    results_img(iimage,next_itimimage)%acell(:)     =acell(:)
1382    results_img(iimage,next_itimimage)%rprim(:,:)   =rprim(:,:)
1383 !  WARNING : Should also store vel and vel_cell of course ... 
1384 !  results_img(iimage,next_itimimage)%vel(:,:)     =vel(:,:)
1385 !  results_img(iimage,next_itimimage)%vel_cell(:,:)=vel_cell(:,:)
1386  end do
1387 
1388  ABI_DEALLOCATE(fcart)
1389  ABI_DEALLOCATE(vel)
1390  ABI_DEALLOCATE(xred)
1391 
1392 end subroutine move_1geo

ABINIT/predict_copy [ Functions ]

[ Top ] [ Functions ]

NAME

 predict_copy

FUNCTION

 Given the past history of images, predict the new set of images.
 Here, simple copy of the previous image.

INPUTS

 itimimage_eff=time index in the history
 list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
       Example : in the NEB of string method, one expect the two end images to be fixed.
 ndynimage=number of dynamical images
 nimage=number of images
 ntimimage_stored=number of time steps stored in the history

OUTPUT

SIDE EFFECTS

 results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
   results_img(:,:)%acell(3)
    at input, history of the values of acell for all images
    at output, the predicted values of acell for all images
   results_img(:,:)%results_gs
    at input, history of the values of energies and forces for all images
   results_img(:,:)%rprim(3,3)
    at input, history of the values of rprim for all images
    at output, the predicted values of rprim for all images
   results_img(:,:)%vel(3,natom)
    at input, history of the values of vel for all images
    at output, the predicted values of vel for all images
   results_img(:,:)%vel_cell(3,3)
    at input, history of the values of vel_cell for all images
    at output, the predicted values of vel_cell for all images
   results_img(:,:)%xred(3,natom)
    at input, history of the values of xred for all images
    at output, the predicted values of xred for all images

PARENTS

      predictimg

CHILDREN

SOURCE

1189 subroutine predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,&
1190 &                       ntimimage_stored,results_img)
1191 
1192 
1193 !This section has been created automatically by the script Abilint (TD).
1194 !Do not modify the following lines by hand.
1195 #undef ABI_FUNC
1196 #define ABI_FUNC 'predict_copy'
1197 !End of the abilint section
1198 
1199  implicit none
1200 
1201 !Arguments ------------------------------------
1202 !scalars
1203  integer,intent(in) :: itimimage_eff,ndynimage,nimage,ntimimage_stored
1204 !arrays
1205  integer,intent(in) :: list_dynimage(ndynimage)
1206  type(results_img_type),intent(inout) :: results_img(nimage,ntimimage_stored)
1207 
1208 !Local variables-------------------------------
1209 !scalars
1210  integer :: idynimage,iimage,next_itimimage
1211 
1212 ! *************************************************************************
1213 
1214  next_itimimage=itimimage_eff+1
1215  if (next_itimimage>ntimimage_stored) next_itimimage=1
1216 
1217  do idynimage=1,ndynimage
1218 
1219    iimage=list_dynimage(idynimage)
1220 
1221    results_img(iimage,next_itimimage)%acell(:)     =results_img(iimage,itimimage_eff)%acell(:)
1222    results_img(iimage,next_itimimage)%rprim(:,:)   =results_img(iimage,itimimage_eff)%rprim(:,:)
1223    results_img(iimage,next_itimimage)%vel(:,:)     =results_img(iimage,itimimage_eff)%vel(:,:)
1224    results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
1225    results_img(iimage,next_itimimage)%xred(:,:)    =results_img(iimage,itimimage_eff)%xred(:,:)
1226 
1227  end do  ! idynimage
1228 
1229 end subroutine predict_copy

ABINIT/predictimg [ Functions ]

[ Top ] [ Functions ]

NAME

 predictimg

FUNCTION

 Given the past history of images, predict the new set of images

INPUTS

 deltae=averaged energy difference used to control convergence over images
 imagealgo_str=name of the algorithm (with images) used
 imgmov=gives the algorithm to be used for prediction of new set of images
 itimimage=time index for image propagation (itimimage+1 is to be predicted here)
 itimimage_eff=time index in the history
 list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
       Example : in the NEB method, or in the string method, one expect the two end images to be fixed.
 mep_param=several parameters for Minimal Energy Path (MEP) search
 mpi_enreg=MPI-parallelisation information
 natom= number of atoms
 ndynimage=number of dynamical images
 nimage=number of images (treated by current proc)
 nimage_tot=total number of images
 ntimimage_stored=number of time steps stored in the history
 pimd_param=several parameters for Path-Integral MD
 prtvolimg=printing volume

OUTPUT

SIDE EFFECTS

 results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
   results_img(:,:)%acell(3)
    at input, history of the values of acell for all images
    at output, the predicted values of acell for all images
   results_img(:,:)%results_gs
    at input, history of the values of energies and forces for all images
   results_img(:,:)%rprim(3,3)
    at input, history of the values of rprim for all images
    at output, the predicted values of rprim for all images
   results_img(:,:)%vel(3,natom)
    at input, history of the values of vel for all images
    at output, the predicted values of vel for all images
   results_img(:,:)%vel_cell(3,3)
    at input, history of the values of vel_cell for all images
    at output, the predicted values of vel_cell for all images
   results_img(:,:)%xred(3,natom)
    at input, history of the values of xred for all images
    at output, the predicted values of xred for all images

PARENTS

      gstateimg

CHILDREN

      predict_copy,predict_ga,predict_neb,predict_pimd,predict_steepest
      predict_string,wrtout

SOURCE

1026 subroutine predictimg(deltae,imagealgo_str,imgmov,itimimage,itimimage_eff,list_dynimage,&
1027 &                     ga_param,mep_param,mpi_enreg,m1geo_param,natom,ndynimage,nimage,nimage_tot,&
1028 &                     ntimimage_stored,pimd_param,prtvolimg,results_img)
1029 
1030  use m_results_gs , only : results_gs_type
1031  use m_predict_neb, only : predict_neb
1032  use m_predict_steepest, only : predict_steepest
1033  use m_predict_pimd,    only : predict_pimd
1034  use m_predict_string, only : predict_string
1035 
1036 !This section has been created automatically by the script Abilint (TD).
1037 !Do not modify the following lines by hand.
1038 #undef ABI_FUNC
1039 #define ABI_FUNC 'predictimg'
1040 !End of the abilint section
1041 
1042  implicit none
1043 
1044 !Arguments ------------------------------------
1045 !scalars
1046  integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom,ndynimage
1047  integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg
1048  character(len=60),intent(in) :: imagealgo_str
1049  real(dp),intent(in) :: deltae
1050  type(mep_type),intent(inout) :: mep_param
1051  type(m1geo_type),intent(inout) :: m1geo_param
1052  type(ga_type),intent(inout) :: ga_param
1053  type(pimd_type),intent(inout) :: pimd_param
1054  type(MPI_type),intent(in) :: mpi_enreg
1055 !arrays
1056  integer,intent(in) :: list_dynimage(ndynimage)
1057  type(results_img_type) :: results_img(nimage,ntimimage_stored)
1058 
1059 !Local variables-------------------------------
1060 !scalars
1061  integer,save :: idum=5
1062  logical :: is_pimd
1063  character(len=500) :: msg
1064 !arrays
1065 
1066 ! *************************************************************************
1067 
1068  is_pimd=(imgmov==9.or.imgmov==10.or.imgmov==13)
1069 
1070 !Write convergence info
1071  write(msg,'(3a)') ch10,&
1072 & '------------------------------------------------------------',ch10
1073  if (prtvolimg<2) write(msg,'(5a)') trim(msg),' ',trim(imagealgo_str),':',ch10
1074 
1075 !Specific case of 4th-order RK algorithm
1076  if (mep_param%mep_solver==4) then
1077    if (mod(itimimage,4)==0) then
1078      write(msg,'(2a,i1,2a)') trim(msg),&
1079 &     ' Fourth-order Runge-Kutta algorithm - final step',ch10
1080      if (itimimage>4) write(msg,'(2a,es11.3,2a)') trim(msg),&
1081 &     ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10
1082      write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1083    else
1084      write(msg,'(2a,i1,2a)') trim(msg),&
1085 &     ' Fourth-order Runge-Kutta algorithm - intermediate step ',mod(itimimage,4),ch10
1086      write(msg,'(2a)') trim(msg),' Computing new intermediate positions...'
1087    end if
1088  else if (is_pimd) then
1089 
1090 !  PIMD
1091    write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1092  else
1093 
1094 !  Other cases
1095    if (itimimage>1) write(msg,'(2a,es11.3,2a)') trim(msg),&
1096 &   ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10
1097    write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1098 
1099  end if
1100 
1101 !Prevent writing if iexit==1, which at present only happens for imgmov==6 algo
1102  if(imgmov/=6 .or. m1geo_param%iexit==0)then
1103    call wrtout(ab_out ,msg,'COLL')
1104    call wrtout(std_out,msg,'COLL')
1105  endif
1106 
1107  select case(imgmov)
1108 
1109  case(0)
1110    call predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,&
1111 &   ntimimage_stored,results_img)
1112 
1113  case(1)
1114    call predict_steepest(itimimage,itimimage_eff,list_dynimage,mep_param,natom,ndynimage,nimage,&
1115 &   ntimimage_stored,results_img)
1116 
1117  case(2)
1118    call predict_string(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
1119 &   ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
1120 
1121  case(4)
1122    call predict_ga(itimimage_eff,idum,ga_param,natom,nimage,&
1123 &   ntimimage_stored,results_img)
1124 
1125  case(5)
1126    call predict_neb(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
1127 &   ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
1128 
1129  case(6)
1130    call move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,ntimimage_stored,results_img)
1131 
1132  case(9, 10, 13)
1133 !    Path Integral Molecular Dynamics
1134    call predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,&
1135 &   ntimimage_stored,pimd_param,prtvolimg,results_img)
1136 
1137  case default
1138 
1139  end select
1140 
1141 end subroutine predictimg

ABINIT/prtimg [ Functions ]

[ Top ] [ Functions ]

NAME

 prtimg

FUNCTION

 Print out results obtained by as ground-state calculation of
 an image of the cell. The printing format is condensed in order
 to facilitate the reading.

INPUTS

  dynimage(nimagetot)=flags defining static/dynamic state of images
  imagealgo_str=name of the algorithm (with images) used
  imgmov=index of algorithm (with images) used
  iout=unit number for output
  mpi_enreg=MPI-parallelisation information
  nimage=number of images stored on current proc
  nimage_tot=total number of images (should be dtset%nimage)
  prt_all_images=true if all images have to be printed out (ignoring dynimage)
  prtvolimg=printing volume for each image
           <0 : nothing
            0 : only a title
            1 : energy, residuals, forces, stresses, velocities, atomic positions
            2 : energy, residuals
  resimg(nimage) <type(results_img_type)>=results of the ground-state computations
                                          for all images treated by current proc

OUTPUT

  (data written to unit iout)

PARENTS

      gstateimg

CHILDREN

      destroy_results_img,gather_results_img,metric,prtxvf,wrtout,xred2xcart

SOURCE

819 subroutine prtimg(dynimage,imagealgo_str,imgmov,iout,mpi_enreg,nimage,nimage_tot,&
820 &                 prt_all_images,prtvolimg,resimg)
821 
822 
823 !This section has been created automatically by the script Abilint (TD).
824 !Do not modify the following lines by hand.
825 #undef ABI_FUNC
826 #define ABI_FUNC 'prtimg'
827 !End of the abilint section
828 
829  implicit none
830 
831 !Arguments ------------------------------------
832 !scalars
833  integer,intent(in) :: nimage_tot,dynimage(nimage_tot),imgmov,iout,nimage,prtvolimg !vz_d
834  logical,intent(in) :: prt_all_images
835  character(len=60),intent(in) :: imagealgo_str
836  type(MPI_type),intent(in) :: mpi_enreg
837 !arrays
838  type(results_img_type),target,intent(inout) :: resimg(nimage)
839 
840 !Local variables-------------------------------
841 !scalars
842  integer :: ii,prtvel
843  logical :: test_img
844  real(dp) :: ucvol_img
845  character(len=500) :: msg
846 !arrays
847  integer,allocatable :: iatfix_img(:,:)
848  real(dp),allocatable :: gmet_img(:,:),gprimd_img(:,:),rmet_img(:,:),xcart_img(:,:)
849  type(results_img_type),pointer :: resimg_all(:)
850 
851 ! ****************************************************************
852 
853  DBG_ENTER('COLL')
854 
855  if (prtvolimg<=0) return
856  if (mpi_enreg%me_cell/=0) return
857 
858 !Gather data
859  if (prtvolimg==1.or.prtvolimg==2) then
860    test_img=(nimage_tot/=1.and.mpi_enreg%paral_img==1)
861    if (test_img) then
862      if (mpi_enreg%me==0)  then
863        ABI_DATATYPE_ALLOCATE(resimg_all,(nimage_tot))
864      end if
865      call gather_results_img(mpi_enreg,resimg,resimg_all,master=0,&
866 &     allgather=.false.,only_one_per_img=.true.)
867    else
868      resimg_all => resimg
869    end if
870  end if
871 
872 !===== First option for the printing volume ===
873  if (prtvolimg==1.and.mpi_enreg%me==0) then
874 
875    prtvel=0;if (imgmov==0.or.imgmov==9.or.imgmov==10.or.imgmov==13) prtvel=1
876 
877    do ii=1,nimage_tot
878      if (dynimage(ii)==1.or.prt_all_images) then
879 
880 !      Title
881        write(msg,'(6a,i4,a,i4,2a)') ch10,&
882 &       '----------------------------------------------------------------------',ch10,&
883 &       ' ',trim(imagealgo_str),' - CELL # ',ii,'/',nimage_tot,ch10,&
884 &       '----------------------------------------------------------------------'
885        call wrtout(iout,msg,'COLL')
886 
887 !      Total energy
888        write(msg,'(2a,es20.12)') ch10,' Total energy for the cell [Ha]: ',resimg_all(ii)%results_gs%etotal
889        call wrtout(iout,msg,'COLL')
890 
891 !      Residuals of the SCF cycle
892        write(msg,'(3a,4(a,es16.8,a))') ch10,&
893 &       ' Residuals from SCF cycle: ',ch10,&
894 &       '    Total energy difference        =',resimg_all(ii)%results_gs%deltae,ch10,&
895 &       '    Maximal forces difference      =',resimg_all(ii)%results_gs%diffor,ch10,&
896 &       '    Max. residual of wave-functions=',resimg_all(ii)%results_gs%residm,ch10,&
897 &       '    Density/potential residual (^2)=',resimg_all(ii)%results_gs%res2,ch10
898        call wrtout(iout,msg,'COLL')
899 
900 !      Cell parameters
901        ABI_ALLOCATE(rmet_img,(3,3))
902        ABI_ALLOCATE(gmet_img,(3,3))
903        ABI_ALLOCATE(gprimd_img,(3,3))
904        call metric(gmet_img,gprimd_img,iout,rmet_img,resimg_all(ii)%rprim,ucvol_img)
905        ABI_DEALLOCATE(rmet_img)
906        ABI_DEALLOCATE(gmet_img)
907        ABI_DEALLOCATE(gprimd_img)
908 
909 !      Positions, forces and velocities
910        ABI_ALLOCATE(iatfix_img,(3,resimg_all(ii)%natom))
911        ABI_ALLOCATE(xcart_img,(3,resimg_all(ii)%natom))
912        iatfix_img=0
913        call xred2xcart(resimg_all(ii)%natom,resimg_all(ii)%rprim,xcart_img,resimg_all(ii)%xred)
914        call prtxvf(resimg_all(ii)%results_gs%fcart,resimg_all(ii)%results_gs%fred,&
915 &       iatfix_img,iout,resimg_all(ii)%natom,prtvel,&
916 &       resimg_all(ii)%vel,xcart_img,resimg_all(ii)%xred)
917        ABI_DEALLOCATE(iatfix_img)
918        ABI_DEALLOCATE(xcart_img)
919 
920 !      Stress tensor
921        write(msg, '(a,es12.4,a)' ) &
922 &       '-Cartesian components of stress tensor (GPa)         [Pressure=',&
923 &       -(resimg_all(ii)%results_gs%strten(1)+resimg_all(ii)%results_gs%strten(2) &
924 &       +resimg_all(ii)%results_gs%strten(3))*HaBohr3_GPa/three,' GPa]'
925        call wrtout(iout,msg,'COLL')
926        write(msg, '(2(a,1p,e16.8))' ) '- sigma(1 1)=',resimg_all(ii)%results_gs%strten(1)*HaBohr3_GPa,&
927 &       '  sigma(3 2)=',resimg_all(ii)%results_gs%strten(4)*HaBohr3_GPa
928        call wrtout(iout,msg,'COLL')
929        write(msg, '(2(a,1p,e16.8))' ) '- sigma(2 2)=',resimg_all(ii)%results_gs%strten(2)*HaBohr3_GPa,&
930 &       '  sigma(3 1)=',resimg_all(ii)%results_gs%strten(5)*HaBohr3_GPa
931        call wrtout(iout,msg,'COLL')
932        write(msg, '(2(a,1p,e16.8))' ) '- sigma(3 3)=',resimg_all(ii)%results_gs%strten(3)*HaBohr3_GPa,&
933 &       '  sigma(2 1)=',resimg_all(ii)%results_gs%strten(6)*HaBohr3_GPa
934        call wrtout(iout,msg,'COLL')
935      end if
936    end do
937  end if
938 
939 
940 !===== 2nd option for the printing volume ===
941  if (prtvolimg==2.and.mpi_enreg%me==0) then
942    write(msg,'(a,1x,a)') ch10,&
943 &   'Cell   Total_energy[Ha]     deltae       diffor       residm         res2'
944    call wrtout(iout,msg,'COLL')
945    do ii=1,nimage_tot
946      if (dynimage(ii)==1.or.prt_all_images) then
947        write(msg,'(1x,i4,2x,es16.8,4(1x,es13.5))') &
948 &       ii,resimg_all(ii)%results_gs%etotal,resimg_all(ii)%results_gs%deltae,&
949 &       resimg_all(ii)%results_gs%diffor,resimg_all(ii)%results_gs%residm,&
950 &       resimg_all(ii)%results_gs%res2
951        call wrtout(iout,msg,'COLL')
952      end if
953    end do
954  end if
955 
956 !=====
957  if (prtvolimg==1.or.prtvolimg==2) then
958    if (test_img.and.mpi_enreg%me==0) then
959      call destroy_results_img(resimg_all)
960      ABI_DATATYPE_DEALLOCATE(resimg_all)
961    end if
962    nullify(resimg_all)
963  end if
964 
965  DBG_EXIT('COLL')
966 
967 end subroutine prtimg