TABLE OF CONTENTS


ABINIT/m_precpred_1geo [ Modules ]

[ Top ] [ Modules ]

NAME

  m_precpred_1geo

FUNCTION

 Single geometry : apply force and stress preconditioner followed by geometry predictor. 
 Choose among the whole set of geometry predictors defined by iomov.

COPYRIGHT

  Copyright (C) 2018 ABINIT group (DCA, XG, GMR, SE)
  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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_precpred_1geo
28 
29  use defs_basis
30  use defs_abitypes
31  use m_errors
32  use m_abimover
33  use m_abihist
34  use m_xmpi
35  use m_nctk
36 #ifdef HAVE_NETCDF
37  use netcdf
38 #endif
39 #if defined HAVE_LOTF
40  use lotfpath
41  use m_pred_lotf
42 #endif
43 
44  use m_fstrings,           only : strcat
45  use m_geometry,           only : chkdilatmx
46  use m_crystal,            only : crystal_init, crystal_free, crystal_t
47  use m_crystal_io,         only : crystal_ncwrite_path
48  use m_pred_bfgs,          only : pred_bfgs, pred_lbfgs
49  use m_pred_delocint,      only : pred_delocint
50  use m_pred_fire,          only : pred_fire
51  use m_pred_isokinetic,    only : pred_isokinetic
52  use m_pred_diisrelax,     only : pred_diisrelax
53  use m_pred_nose,          only : pred_nose
54  use m_pred_srkhna14,      only : pred_srkna14
55  use m_pred_isothermal,    only : pred_isothermal
56  use m_pred_verlet,        only : pred_verlet
57  use m_pred_velverlet,     only : pred_velverlet
58  use m_pred_moldyn,        only : pred_moldyn
59  use m_pred_langevin,      only : pred_langevin
60  use m_pred_steepdesc,     only : pred_steepdesc
61  use m_pred_simple,        only : pred_simple, prec_simple
62  use m_pred_hmc,           only : pred_hmc
63  use m_generate_training_set, only : generate_training_set
64 
65  implicit none
66 
67  private

ABINIT/precpred_1geo [ Functions ]

[ Top ] [ Functions ]

NAME

 mover

FUNCTION

 Single geometry : apply force and stress preconditioner followed by geometry predictor.
 Choose among the whole set of geometry predictors defined by iomov.

INPUTS

  dtfil <type(datafiles_type)>=variables related to files
  dtset <type(dataset_type)>=all input variables for this dataset
   | mband=maximum number of bands
   | mgfft=maximum size of 1D FFTs
   | mkmem =number of k points treated by this node
   |  angular momentum for nonlocal pseudopotential
   | mpw=maximum dimensioned size of npw.
   | natom=number of atoms in unit cell
   |  except on first call (hartree/bohr); updated on output
   | nfft=(effective) number of FFT grid points (for this processor)
   |      for the "coarse" grid (see NOTES below)
   | nkpt=number of k points.
   | nspden=number of spin-density components
   | nsppol=1 for unpolarized, 2 for spin-polarized
   | nsym=number of symmetry elements in space group
  mpi_enreg=informations about MPI parallelization
  rprimd(3,3)=dimensional primitive translations (bohr)

OUTPUT

SIDE EFFECTS

 Rest of i/o is related to lda
  xred(3,natom)=reduced dimensionless atomic coordinates; updated on output
  write_HIST = optional, default is true, flag to disble the write of the HIST file

NOTES

 This subroutine uses the arguments natom, xred, 
 vis, and dtion (the last two contained in dtset) to make
 molecular dynamics updates.  

PARENTS

CHILDREN

SOURCE

121 subroutine precpred_1geo(ab_mover,ab_xfh,amu_curr,deloc,dt_chkdilatmx,comm_cell,dilatmx,filnam_ds4,hist,hmctt,&
122 & icycle,iexit,itime,mttk_vars,nctime,ncycle,nerr_dilatmx,npsp,ntime,rprimd_orig,skipcycle,usewvl)
123 
124 
125 !This section has been created automatically by the script Abilint (TD).
126 !Do not modify the following lines by hand.
127 #undef ABI_FUNC
128 #define ABI_FUNC 'precpred_1geo'
129 !End of the abilint section
130 
131 implicit none
132 
133 !Arguments ------------------------------------
134 !scalars
135 integer, intent(in) :: comm_cell,dt_chkdilatmx,hmctt,icycle,itime,nctime,npsp,ntime
136 integer, intent(inout) :: iexit,ncycle,nerr_dilatmx
137 integer, intent(in) :: usewvl
138 real(dp), intent(in) :: dilatmx
139 logical, intent(inout) :: skipcycle
140 character(len=fnlen), intent(in) :: filnam_ds4
141 type(ab_xfh_type),intent(inout) :: ab_xfh
142 type(abihist), intent(inout) :: hist
143 type(abimover), intent(in) :: ab_mover
144 type(delocint), intent(inout) :: deloc
145 type(mttk_type), intent(inout) :: mttk_vars
146 !arrays
147 real(dp), intent(in) :: amu_curr(ab_mover%ntypat)
148 real(dp), intent(in) :: rprimd_orig(3,3)
149 
150 !Local variables-------------------------------
151 !scalars
152 integer,parameter :: master=0
153 integer :: ii,me,nloop
154 logical :: DEBUG=.FALSE.
155 character(len=500) :: message
156 character(len=500) :: dilatmx_errmsg
157 character(len=fnlen) :: filename
158 type(abiforstr) :: preconforstr ! Preconditioned forces and stress
159 type(crystal_t) :: crystal
160 !arrays
161 real(dp) :: acell(3),rprimd(3,3)
162 real(dp), allocatable :: xred(:,:)
163 
164 ! ***************************************************************
165 
166  me=xmpi_comm_rank(comm_cell)
167 
168 !Precondition forces, stress and energy
169  call abiforstr_ini(preconforstr,ab_mover%natom)
170 
171  if (ab_mover%goprecon>0 .or. iexit==1)then
172    call prec_simple(ab_mover,preconforstr,hist,icycle,itime,iexit)
173  end if
174 
175 !Call to each predictor
176 !MT->GAF: dirty trick to predict vel(t)
177 !do a double loop: 1- compute vel, 2- exit
178  nloop=1
179 
180  if (nctime>0.and.iexit==1) then
181    iexit=0;nloop=2
182  end if
183 
184  do ii=1,nloop
185    if (ii==2) iexit=1
186 
187    select case (ab_mover%ionmov)
188    case (1)
189      call pred_moldyn(ab_mover,hist,icycle,itime,ncycle,ntime,DEBUG,iexit)
190    case (2,3)
191      call pred_bfgs(ab_mover,ab_xfh,preconforstr,hist,ab_mover%ionmov,itime,DEBUG,iexit)
192    case (4,5)
193      call pred_simple(ab_mover,hist,iexit)
194    case (6,7)
195      call pred_verlet(ab_mover,hist,ab_mover%ionmov,itime,ntime,DEBUG,iexit)
196    case (8)
197      call pred_nose(ab_mover,hist,itime,ntime,DEBUG,iexit)
198    case (9)
199      call pred_langevin(ab_mover,hist,icycle,itime,ncycle,ntime,DEBUG,iexit,skipcycle)
200    case (10,11)
201      call pred_delocint(ab_mover,ab_xfh,deloc,preconforstr,hist,ab_mover%ionmov,itime,DEBUG,iexit)
202    case (12)
203      call pred_isokinetic(ab_mover,hist,itime,ntime,DEBUG,iexit)
204    case (13)
205      call pred_isothermal(ab_mover,hist,itime,mttk_vars,ntime,DEBUG,iexit)
206    case (14)
207      call pred_srkna14(ab_mover,hist,icycle,DEBUG,iexit,skipcycle)
208    case (15)
209      call pred_fire(ab_mover, ab_xfh,preconforstr,hist,ab_mover%ionmov,itime,DEBUG,iexit)
210    case (20)
211      call pred_diisrelax(ab_mover,hist,itime,ntime,DEBUG,iexit)
212    case (21)
213      call pred_steepdesc(ab_mover,preconforstr,hist,itime,DEBUG,iexit)
214    case (22)
215      call pred_lbfgs(ab_mover,ab_xfh,preconforstr,hist,ab_mover%ionmov,itime,DEBUG,iexit)
216 #if defined HAVE_LOTF
217    case (23)
218      call pred_lotf(ab_mover,hist,itime,icycle,DEBUG,iexit)
219 #endif
220    case (24)
221      call pred_velverlet(ab_mover,hist,itime,ntime,DEBUG,iexit)
222    case (25)
223      call pred_hmc(ab_mover,hist,itime,icycle,ntime,hmctt,DEBUG,iexit)
224    case (27)
225      !In case of ionmov 27, all the atomic configurations have been computed at the
226      !begining of the routine in generate_training_set, thus we just need to increase the indexes
227      !in the hist
228      hist%ihist = abihist_findIndex(hist,+1)
229 
230    case default
231      write(message,"(a,i0)") "Wrong value of ionmov: ",ab_mover%ionmov
232      MSG_ERROR(message)
233    end select
234 
235  end do
236 
237  ABI_ALLOCATE(xred,(3,ab_mover%natom))
238  call hist2var(acell,hist,ab_mover%natom,rprimd,xred,DEBUG)
239 
240  ! check dilatmx here and correct if necessary
241  if (usewvl == 0) then
242    call chkdilatmx(dt_chkdilatmx,dilatmx,rprimd,rprimd_orig,dilatmx_errmsg)
243    _IBM6("dilatxm_errmsg: "//TRIM(dilatmx_errmsg))
244    if (LEN_TRIM(dilatmx_errmsg) /= 0) then
245      MSG_WARNING(dilatmx_errmsg)
246      nerr_dilatmx = nerr_dilatmx+1
247      if (nerr_dilatmx > 3) then
248        ! Write last structure before aborting, so that we can restart from it.
249        ! zion is not available, but it's not useful here.
250        if (me == master) then
251          ! Init crystal
252          hist%ihist = abihist_findIndex(hist,-1)
253          call hist2var(acell,hist,ab_mover%natom,rprimd,xred,DEBUG)
254          call crystal_init(amu_curr,crystal,0,ab_mover%natom,&
255 &         npsp,ab_mover%ntypat,ab_mover%nsym,rprimd,ab_mover%typat,xred,&
256 &         [(-one, ii=1,ab_mover%ntypat)],ab_mover%znucl,2,.False.,.False.,"dilatmx_structure",&
257 &         symrel=ab_mover%symrel,tnons=ab_mover%tnons,symafm=ab_mover%symafm)
258 
259 #ifdef HAVE_NETCDF
260          ! Write netcdf file
261          filename = strcat(filnam_ds4, "_DILATMX_STRUCT.nc")
262          NCF_CHECK(crystal_ncwrite_path(crystal, filename))
263 #endif
264          call crystal_free(crystal)
265        end if
266        call xmpi_barrier(comm_cell)
267        write (dilatmx_errmsg, '(a,i0,3a)') &
268         'Dilatmx has been exceeded too many times (', nerr_dilatmx, ')',ch10, &
269         'Restart your calculation from larger lattice vectors and/or a larger dilatmx'
270        MSG_ERROR_CLASS(dilatmx_errmsg, "DilatmxError")
271      end if
272    else
273      nerr_dilatmx=0
274    end if
275  end if
276 
277  call abiforstr_fin(preconforstr)
278  ABI_FREE(xred)
279 
280 end subroutine precpred_1geo