Added FIRE minimization algorithm and plugin flags for ARTn plugin

This commit is contained in:
matic 2021-11-15 16:05:21 +01:00
parent 993f7297f1
commit 95a5a19935
12 changed files with 1135 additions and 303 deletions

View File

@ -1051,11 +1051,11 @@ MODULE input_parameters
CHARACTER(len=80) :: ion_dynamics = 'none'
!! set how ions should be moved
CHARACTER(len=80) :: ion_dynamics_allowed(10)
CHARACTER(len=80) :: ion_dynamics_allowed(11)
!! allowed options for ion\_dynamics.
DATA ion_dynamics_allowed / 'none', 'sd', 'cg', 'langevin', &
'damp', 'verlet', 'bfgs', 'beeman',&
'langevin-smc', 'ipi' /
'langevin-smc', 'ipi', 'fire' /
REAL(DP) :: ion_radius(nsx) = 0.5_DP
!! pseudo-atomic radius of the i-th atomic species (CP only).
@ -1193,6 +1193,17 @@ MODULE input_parameters
REAL(DP) :: w_1 = 0.5E-1_DP
REAL(DP) :: w_2 = 0.5_DP
!
! Parameters for minimization with the FIRE algorithm
!
INTEGER :: fire_nmin = 5 ! minimum number of steps for time step increase
REAL(DP) :: fire_f_inc = 1.1_DP ! factor for time step increase
REAL(DP) :: fire_f_dec = 0.5_DP ! factor for time step decrease
REAL(DP) :: fire_alpha_init = 0.2_DP ! initial value of mixing factor
REAL(DP) :: fire_falpha = 0.99_DP ! modify the mixing factor
REAL(DP) :: fire_dtmax = 10.0_DP ! maximum time step; calculated as dtmax = fire_dtmax*dt
!
!
NAMELIST / ions / ion_dynamics, iesr, ion_radius, ion_damping, &
ion_positions, ion_velocities, ion_temperature, &
@ -1201,7 +1212,10 @@ MODULE input_parameters
refold_pos, upscale, delta_t, pot_extrapolation, &
wfc_extrapolation, nraise, remove_rigid_rot, &
trust_radius_max, trust_radius_min, &
trust_radius_ini, w_1, w_2, bfgs_ndim
trust_radius_ini, w_1, w_2, bfgs_ndim, &
fire_nmin, fire_f_inc, fire_f_dec, fire_alpha_init, &
fire_falpha, fire_dtmax
!=----------------------------------------------------------------------------=!

View File

@ -238,6 +238,7 @@ CONTAINS
CALL delete_if_present( trim( file_path ) // '.update' )
CALL delete_if_present( trim( file_path ) // '.md' )
CALL delete_if_present( trim( file_path ) // '.bfgs' )
CALL delete_if_present( trim( file_path ) // '.fire' )
ENDIF
!
RETURN

View File

@ -30,6 +30,7 @@ SUBROUTINE plugin_arguments()
use_plumed = .false.
use_pw2casino = .false.
use_environ = .false.
use_partn = .false.
!
DO iiarg = 1, nargs
CALL get_command_argument( iiarg, plugin_name)
@ -51,6 +52,9 @@ SUBROUTINE plugin_arguments()
IF ( TRIM(arg)=='environ' ) THEN
use_environ = .true.
ENDIF
IF ( TRIM(arg)=='partn' ) THEN
use_partn = .true.
ENDIF
ENDIF
ENDDO
!
@ -77,6 +81,8 @@ END SUBROUTINE plugin_arguments
!
CALL mp_bcast(use_environ,root,comm)
!
CALL mp_bcast(use_partn,root,comm)
!
! write(0,*) "use_plumed: ", use_plumed
!
RETURN

View File

@ -27,6 +27,7 @@ MODULE plugin_flags
CHARACTER(LEN=256), PUBLIC :: plugin_name
LOGICAL, PUBLIC :: use_plumed
LOGICAL, PUBLIC :: use_pw2casino
LOGICAL, PUBLIC :: use_environ
LOGICAL, PUBLIC :: use_environ
LOGICAL, PUBLIC :: use_partn
!
END MODULE plugin_flags

View File

@ -532,6 +532,14 @@ MODULE read_namelists_module
w_1 = 0.01_DP
w_2 = 0.50_DP
!
! FIRE minimization defaults
!
fire_nmin = 5 ! minimum number of steps P > 0 before dt increas
fire_f_inc = 1.1_DP ! factor for time step increase
fire_f_dec = 0.5_DP ! factor for time step decrease
fire_alpha_init = 0.20_DP ! initial value of mixing factor
fire_falpha = 0.99_DP ! modification of the mixing factor
fire_dtmax = 10.0_DP ! factor for calculating dtmax
RETURN
!
END SUBROUTINE

View File

@ -2409,7 +2409,17 @@ input_description -distribution {Quantum Espresso} -package PWscf -program pw.x
dynamics for structural relaxation
Can be used for constrained
optimisation: see @ref CONSTRAINTS card
}
}
opt -val 'fire' {
use the FIRE minimization algorithm employing the
semi-implicit Euler intergration scheme
see:
Bitzek et al.,PRL, 97, 170201, (2006), doi: 10.1103/PhysRevLett.97.170201
Guenole et al.,CMS, 175, 109584, (2020), doi: 10.1016/j.commatsci.2020.109584
Can be used for constrained
optimisation: see @ref CONSTRAINTS card
}
info {
@b CASE ( @ref calculation == 'md' )
}
@ -2706,6 +2716,50 @@ input_description -distribution {Quantum Espresso} -package PWscf -program pw.x
}
}
}
group {
label {
keywords used only in the FIRE minimization algorithm
}
var fire_alpha_init -type REAL {
default { 0.2D0 }
info {
Initial value of the alpha mixing factor in the FIRE minimization scheme;
recommended values are between 0.1 and 0.3
}
}
var fire_falpha -type REAL {
default { 0.99D0 }
info {
Scaling of the alpha mixing parameter for steps with P > 0;
}
}
var fire_nmin -type INTEGER {
default { 5 }
info {
Minimum number of steps with P > 0 before increase of @ref dt
}
}
var fire_f_inc -type REAL {
default { 1.1D0 }
info {
Factor for increasing @ref dt
}
}
var fire_f_dec -type REAL {
default { 0.5D0 }
info {
Factor for decreasing @ref dt
}
}
var fire_dtmax -type REAL {
default { 10.D0 }
info {
Determines the maximum value of @ref dt in the FIRE minimization;
dtmax = fire_dtmax*@ref dt
}
}
}
}

File diff suppressed because it is too large Load Diff

View File

@ -3,7 +3,7 @@
------------------------------------------------------------------------
INPUT FILE DESCRIPTION
Program: pw.x / PWscf / Quantum Espresso (version: 6.8)
Program: pw.x / PWscf / Quantum Espresso
------------------------------------------------------------------------
@ -1406,6 +1406,31 @@ NAMELIST: &SYSTEM
can be used only when "lda_plus_u_kind" = 2.
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: dmft
Type: LOGICAL
Default: .FALSE.
Status: Requires compilation with hdf5 support
Description: If true, nscf calculation will exit in restart mode, scf calculation
will restart from there if DMFT updates are provided as hdf5 archive.
Scf calculation should be used only with "electron_maxstep" = 1.
"K_POINTS" have to be identical and given explicitly with "nosym".
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: dmft_prefix
Type: CHARACTER
Default: "prefix"
Description: prepended to hdf5 archive: dmft_prefix.h5
DMFT update should be provided in group/dataset as:
- dft_misc_input/band_window with dimension [1, number of k-points, 2 (real + complex)]
- dft_update/delta_N with dimension [number of k-points, number of correlated orbitals,
number of correlated orbitals, 2 (real + complex)]
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: ensemble_energies
@ -2227,6 +2252,14 @@ NAMELIST: &ELECTRONS
'paro', 'ParO' :
ParO iterative diagonalization
'rmm-davidson', 'rmm-paro' :
RMM-DIIS iterative diagonalization.
To stabilize the SCF loop
RMM-DIIS is alternated with calls to Davidson or
ParO solvers depending on the string used.
Other variables that can be used to tune the behavior of
RMM-DIIS are: "diago_rmm_ndim" and "diago_rmm_conv"
+--------------------------------------------------------------------
+--------------------------------------------------------------------
@ -2278,6 +2311,23 @@ NAMELIST: &ELECTRONS
total energy, forces, and other ground-state properties).
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: diago_rmm_ndim
Type: INTEGER
Default: 4
Description: Max dimension of the iterative subspace for RMM-DIIS diagonalization
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: diago_rmm_conv
Type: LOGICAL
Default: .FALSE.
Description: If .TRUE. during the SCF loop the RMM-DIIS is reiterated until all bands
are converged or up to a max of 8 times.
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: efield
@ -2459,6 +2509,16 @@ NAMELIST: &IONS
Can be used for constrained
optimisation: see "CONSTRAINTS" card
'fire' :
use the FIRE minimization algorithm employing the
semi-implicit Euler intergration scheme
see:
Bitzek et al.,PRL, 97, 170201, (2006), doi: 10.1103/PhysRevLett.97.170201
Guenole et al.,CMS, 175, 109584, (2020), doi: 10.1016/j.commatsci.2020.109584
Can be used for constrained
optimisation: see "CONSTRAINTS" card
CASE ( "calculation" == 'md' )
'verlet' :
@ -2767,6 +2827,61 @@ NAMELIST: &IONS
\\\---
///---
KEYWORDS USED ONLY IN THE FIRE MINIMIZATION ALGORITHM
+--------------------------------------------------------------------
Variable: fire_alpha_init
Type: REAL
Default: 0.2D0
Description: Initial value of the alpha mixing factor in the FIRE minimization scheme;
recommended values are between 0.1 and 0.3
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: fire_falpha
Type: REAL
Default: 0.99D0
Description: Scaling of the alpha mixing parameter for steps with P > 0;
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: fire_nmin
Type: INTEGER
Default: 5
Description: Minimum number of steps with P > 0 before increase of "dt"
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: fire_f_inc
Type: REAL
Default: 1.1D0
Description: Factor for increasing "dt"
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: fire_f_dec
Type: REAL
Default: 0.5D0
Description: Factor for decreasing "dt"
+--------------------------------------------------------------------
+--------------------------------------------------------------------
Variable: fire_dtmax
Type: REAL
Default: 10.D0
Description: Determines the maximum value of "dt" in the FIRE minimization;
dtmax = fire_dtmax*"dt"
+--------------------------------------------------------------------
\\\---
===END OF NAMELIST======================================================
@ -3817,4 +3932,4 @@ CARD: ATOMIC_FORCES
===END OF CARD==========================================================
This file has been created by helpdoc utility on Fri Jul 16 11:28:04 CEST 2021
This file has been created by helpdoc utility on Mon Nov 15 13:00:30 CET 2021

View File

@ -15,7 +15,7 @@ $ECHO "of a simple molecule, CO, and of an Al (001) slab."
$ECHO "In the latter case the relaxation is performed in two ways:"
$ECHO "1) using the quasi-Newton BFGS algorithm"
$ECHO "2) using a damped dynamics algorithm."
$ECHO "3) using the FIRE algorithm."
# set the needed environment variables
. ../../../environment_variables
@ -121,6 +121,50 @@ $ECHO " cleaning $TMP_DIR...\c"
rm -rf $TMP_DIR/CO*
$ECHO " done"
cat > co.fire.in << EOF
&CONTROL
calculation = "relax",
prefix = "CO",
pseudo_dir = "$PSEUDO_DIR",
outdir = "$TMP_DIR",
/
&SYSTEM
ibrav = 0,
nat = 2,
ntyp = 2,
ecutwfc = 24.D0,
ecutrho = 144.D0,
/
&ELECTRONS
conv_thr = 1.D-7,
mixing_beta = 0.7D0,
/
&IONS
ion_dynamics = 'fire'
/
CELL_PARAMETERS bohr
12.0 0.0 0.0
0.0 12.0 0.0
0.0 0.0 12.0
ATOMIC_SPECIES
O 1.00 O.pz-rrkjus.UPF
C 1.00 C.pz-rrkjus.UPF
ATOMIC_POSITIONS {bohr}
C 2.256 0.0 0.0
O 0.000 0.0 0.0 0 0 0
K_POINTS {Gamma}
EOF
$ECHO " running the geometry relaxation for CO using FIRE...\c"
$PW_COMMAND < co.fire.in > co.fire.out
check_failure $?
$ECHO " done"
# clean TMP_DIR
$ECHO " cleaning $TMP_DIR...\c"
rm -rf $TMP_DIR/CO*
$ECHO " done"
# self-consistent calculation
cat > al001.rx.in << EOF
&CONTROL
@ -224,6 +268,59 @@ $PW_COMMAND < al001.mm.in > al001.mm.out
check_failure $?
$ECHO " done"
# clean TMP_DIR
$ECHO " cleaning $TMP_DIR...\c"
rm -rf $TMP_DIR/Al*
$ECHO " done"
cat > al001.fire.in << EOF
&CONTROL
calculation = "relax",
dt = 30.D0,
pseudo_dir = "$PSEUDO_DIR",
outdir = "$TMP_DIR",
prefix = "Al"
/
&SYSTEM
ibrav = 6,
celldm(1) = 5.3033D0,
celldm(3) = 8.D0,
nat = 7,
ntyp = 1,
ecutwfc = 12.D0,
occupations = "smearing",
smearing = "marzari-vanderbilt",
degauss = 0.05D0,
/
&ELECTRONS
conv_thr = 1.D-7,
mixing_beta = 0.3D0,
/
&IONS
ion_dynamics = "fire",
pot_extrapolation = "second_order",
wfc_extrapolation = "second_order",
/
ATOMIC_SPECIES
Al 1.D0 Al.pz-vbc.UPF
ATOMIC_POSITIONS alat
Al 0.5000000 0.5000000 -2.121320
Al 0.0000000 0.0000000 -1.414213
Al 0.5000000 0.5000000 -0.707107
Al 0.0000000 0.0000000 0.000000
Al 0.5000000 0.5000000 0.707107
Al 0.0000000 0.0000000 1.414213
Al 0.5000000 0.5000000 2.121320
K_POINTS
3
0.125 0.125 0.0 1.0
0.125 0.375 0.0 2.0
0.375 0.375 0.0 1.0
EOF
$ECHO " running the geometry relaxation for Al (001) using the FIRE algorithm...\c"
$PW_COMMAND < al001.fire.in > al001.fire.out
check_failure $?
$ECHO " done"
# clean TMP_DIR
$ECHO " cleaning $TMP_DIR...\c"
rm -rf $TMP_DIR/Al*

View File

@ -39,10 +39,12 @@ MODULE dynamics_module
!
SAVE
PRIVATE
PUBLIC :: verlet, proj_verlet, terminate_verlet, &
PUBLIC :: verlet, proj_verlet, terminate_verlet, fire, &
langevin_md, smart_MC, allocate_dyn_vars, deallocate_dyn_vars
PUBLIC :: temperature, refold_pos, vel
PUBLIC :: dt, delta_t, nraise, control_temp, thermostat
! FIRE parameters
PUBLIC :: fire_nmin, fire_f_inc, fire_f_dec, fire_alpha_init, fire_falpha, fire_dtmax
!
!
REAL(DP) :: dt
@ -69,6 +71,18 @@ MODULE dynamics_module
!! true if this is the first ionic iteration
CHARACTER(len=10) :: thermostat
!! the thermostat used to control the temperature
INTEGER :: fire_nmin
!! FIRE: minimum number of steps for time step increase
REAL(DP) :: fire_f_inc
!! FIRE: factor for time step increase
REAL(DP) :: fire_f_dec
!! FIRE: factor for time step decrease
REAL(DP) :: fire_alpha_init
!! FIRE: initial value of mixing factor
REAL(DP) :: fire_falpha
!! FIRE: modify the mixing factor
REAL(DP) :: fire_dtmax
!! FIRE: factor for max time step
REAL(DP), ALLOCATABLE :: tau_smart(:,:)
!! used for smart Monte Carlo to store the atomic position of the
!! previous step.
@ -1018,6 +1032,299 @@ CONTAINS
DEALLOCATE( step )
!
END SUBROUTINE proj_verlet
SUBROUTINE fire( conv_ions )
!------------------------------------------------------------------------
!! This routine performs one step of structural relaxation using
!! the FIRE (Fast Inertial Relaxation Engine) algorithm using the
!! semi-implicit Euler integration scheme with an energy monitor;
!!
!! References: (1) Bitzek et al., Phys. Rev. Lett., 97, 170201, (2006),
!! doi: 10.1103/PhysRevLett.97.170201
!! (2) Shuang et al., Comp. Mater. Sci., 156, 135-141 (2019),
!! doi: 10.1016/j.commatsci.2018.09.049
!! (3) (FIRE 2.0) Guénolé et al., Comp. Mater. Sci., 175, 109584, (2020),
!! doi: 10.1016/j.commatsci.2020.109584
!!
!! There are seven global parameters that can be modified by the user:
!!
!! dt .... initial time step of calculation
!! fire_nmin ... number of steps with P > 0 before dt is increased
!! fire_f_inc ... factor for the increase of the time step
!! fire_f_dec ... factor for the decrease of the time step
!! fire_alpha_init ... initial value of the velocity mixing factor
!! fire_falpha ... modifies the velocity mixing factor
!! fire_dtmax ... maximum time step; calculated as dtmax = fire_dtmax*dt
!!
!! Defaults are (taken from ref (2)):
!! dt = 20.0 (the default time step in PW == 20.0 a.u. or 0.9674 fs )
!! fire_f_inc = 1.1
!! fire_f_dec = 0.5
!! fire_f_alpha_init = 0.2
!! fire_falpha = 0.99
!!-----------------------------------------------------------------------
USE ions_base, ONLY : nat, ityp, tau, if_pos
USE cell_base, ONLY : alat
USE ener, ONLY : etot
USE force_mod, ONLY : force
USE relax, ONLY : epse, epsf
USE control_flags, only : istep, lconstrain
!
USE constraints_module, ONLY : remove_constr_force, remove_constr_vec, check_constraint
! TB
USE extfield, ONLY : relaxz
!
IMPLICIT NONE
!
LOGICAL, INTENT(OUT) :: conv_ions
!
REAL(DP), ALLOCATABLE :: step(:,:)
REAL(DP) :: norm_step, etotold, delta(3)
INTEGER :: na,i
LOGICAL :: file_exists
!
REAL(DP), PARAMETER :: step_max = 0.6D0 ! bohr
!
REAL(DP), EXTERNAL :: dnrm2,ddot
!
! FIRE local variables
!
REAL(DP) :: P, alpha, fmax, dt_max, dt_curr
INTEGER :: nsteppos
!
! FIRE parameters; read from input ...
!
INTEGER :: Nmin ! minimum number of steps for time step increase
REAL(DP) :: f_inc ! factor for time step increase
REAL(DP) :: f_dec ! factor for time step decrease
REAL(DP) :: alpha_init ! initial value of mixing factor
REAL(DP) :: falpha ! modify the mixing factor
!
ALLOCATE( step( 3, nat ) )
!
! set up local variables for global input parameters (better readability) ...
!
Nmin = fire_nmin
f_inc = fire_f_inc
f_dec = fire_f_dec
alpha_init = fire_alpha_init
falpha = fire_falpha
!
! initialize alpha
!
tau_new(:,:) = 0.D0
alpha = alpha_init
! maximum time step
dt_curr = dt
dt_max = fire_dtmax*dt
! acc_old and vel_old are here to keep track of acceleration/velocity in the previous time step
nsteppos = 0
conv_ions = .FALSE.
!
CALL seqopn( 4, 'fire', 'FORMATTED', file_exists )
!
!
IF ( file_exists ) THEN
!
! ... the file is read ...
!
READ( UNIT = 4, FMT = * ) etotold, nsteppos, dt_curr, alpha
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
!
ELSE
!
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
! ... atoms are refolded in the central box
!
IF ( refold_pos ) CALL refold_tau()
!
vel(:,:) = 0.D0
acc(:,:) = 0.D0
etotold = etot
istep = 0
WRITE( UNIT = stdout, &
FMT = '(/,5X,"Minimization using the FIRE algorithm")' )
!
! write out the input parameters
WRITE (UNIT = stdout, &
FMT = '(/,5X,"FIRE input parameters:")')
WRITE (UNIT = stdout, &
FMT = '(/,5X," fire_nmin = ",I2," fire_f_inc = ", F4.2, &
" fire_f_dec = ",F4.2," fire_alpha = ",F4.2, &
" fire_falpha = ",F4.2, " dtmax = ",F5.1 )') &
fire_nmin, fire_f_inc, fire_f_dec, &
fire_alpha_init, fire_falpha, dt_max
!
ENDIF
!
IF ( lconstrain ) THEN
!
! ... we first remove the component of the force along the
! ... constraint gradient (this constitutes the initial guess
! ... for the calculation of the lagrange multipliers)
!
write (*,*) "Called remove constr"
CALL remove_constr_force( nat, tau, if_pos, ityp, alat, force )
!
#if ! defined (__REDUCE_OUTPUT)
!
WRITE( stdout, '(/,5X,"Constrained forces (Ry/au):",/)')
!
DO na = 1, nat
!
WRITE( stdout, &
'(5X,"atom ",I3," type ",I2,3X,"force = ",3F14.8)' ) &
na, ityp(na), force(:,na)
!
ENDDO
!
WRITE( stdout, &
'(/5X,"Total force = ",F12.6)') dnrm2( 3*nat, force, 1 )
!
#endif
!
ENDIF
!
! ... check if convergence for structural minimization is achieved
!
conv_ions = ( etotold - etot ) < epse
conv_ions = conv_ions .and. ( MAXVAL( ABS( force ) ) < epsf )
!
IF ( conv_ions ) THEN
!
WRITE( UNIT = stdout, &
FMT = '(/,5X,"FIRE: convergence achieved in " &
& ,I3," steps")' ) istep
WRITE( UNIT = stdout, &
FMT = '(/,5X,"End of FIRE minimization")' )
WRITE( UNIT = stdout, &
FMT = '(/,5X,"Final energy = ",F18.10," Ry"/)' ) etot
!
CALL output_tau( .TRUE., .TRUE. )
!
RETURN
!
ENDIF
!
istep = istep + 1
WRITE( stdout, '(/,5X,"Entering FIRE :",&
& T28,"iteration",T37," = ",I5)' ) istep
!
! calculate acceleration
!
acc(:,:) = force(:,:) / alat / amu_ry
!
! calculate the projection of the velocity on the force
P = ddot(3*nat,force, 1, vel, 1)
!
step(:,:) = 0.0_DP
!
IF ( P < 0.0_DP ) THEN
! FIRE 2.0 algorithm: if P < 0 go back by half a step
! for details see reference (2), doi: 10.1016/j.commatsci.2020.109584
step(:,:) = step(:,:) - 0.5_DP*dt_curr*vel(:,:)
ENDIF
!
! ... manipulate the time step ...
!
! NOTEs:
! in original FIRE the condition is P > 0,
! however to prevent the time step decrease in the first step where v=0 (P=0 and etot=etotold)
! the equality was changed to P >= 0
! The energy difference criterion is also added to prevent the minimization from going uphill
!
IF ( P >= 0.0_DP .AND. (etot - etotold) <= 0.D0 ) THEN
!
!
nsteppos = nsteppos + 1
! increase time step and modify mixing factor only after Nmin steps in positive direction
IF ( nsteppos > Nmin ) THEN
dt_curr = MIN(dt_curr*f_inc, dt_max )
alpha = alpha*falpha
END IF
ELSE
!
! set velocity to 0; return alpha to the initial value; reduce time step
!
vel(:,:) = 0.D0
alpha = alpha_init
nsteppos = 0
dt_curr = dt_curr*f_dec
END IF
! report current parameters
WRITE (stdout, '(/,5X, "FIRE Parameters: P = ", F10.8 ", dt = " F5.2", &
alpha = " F5.3, " nsteppos = ", I3, " at step", I3, /)' ) P, dt_curr, alpha, nsteppos, istep
!
! calculate v(t+dt) = v(t) + a(t)*dt
!
vel(:,:) = vel(:,:) + dt_curr*acc(:,:)
!
! velocity mixing
!
vel(:,:) = (1.D0 - alpha)*vel(:,:) + alpha*force(:,:)*dnrm2(3*nat,vel,1)/dnrm2(3*nat,force,1)
!
!
! ... the velocity of fixed ions must be zero
!
vel = vel * DBLE( if_pos )
!
IF ( lconstrain ) THEN
!
! apply constraints to the velocity as well
!
CALL remove_constr_vec( nat, tau, if_pos, ityp, alat, vel )
ENDIF
!
!
! calculate the displacement x(t+dt) = x(t) + v(t+dt)*dt
!
step(:,:) = step(:,:) + vel(:,:)*dt_curr
!
norm_step = dnrm2( 3*nat, step, 1 )
!
step(:,:) = step(:,:) / norm_step
!
! keep the step within a threshold (taken from damped dynamics)
!
tau_new(:,:) = tau(:,:) + step(:,:)*MIN(norm_step, step_max/alat)
!
IF ( lconstrain ) THEN
!
! ... check if the new positions satisfy the constrain equation
!
CALL check_constraint( nat, tau_new, tau, &
force, if_pos, ityp, alat, dt, amu_ry )
!
ENDIF
!
! ... save the needed quantities to a file
!
CALL seqopn( 4, 'fire', 'FORMATTED', file_exists )
!
WRITE( UNIT = 4, FMT = * ) etot, nsteppos, dt_curr, alpha
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
!
! ... displace tau, and output new positions
!
tau(:,:) = tau_new(:,:)
!
!
#if ! defined (__REDUCE_OUTPUT)
!
CALL output_tau( .FALSE., .FALSE. )
!
#endif
!
DEALLOCATE( step )
!
END SUBROUTINE fire
!
!
!------------------------------------------------------------------------

View File

@ -56,7 +56,13 @@ SUBROUTINE iosys()
dt_ => dt, &
delta_t_ => delta_t, &
nraise_ => nraise, &
refold_pos_ => refold_pos
refold_pos_ => refold_pos, &
fire_nmin_ => fire_nmin, &
fire_f_inc_ => fire_f_inc, &
fire_f_dec_ => fire_f_dec, &
fire_alpha_init_ => fire_alpha_init, &
fire_falpha_ => fire_falpha, &
fire_dtmax_ => fire_dtmax
!
USE extfield, ONLY : tefield_ => tefield, &
dipfield_ => dipfield, &
@ -287,7 +293,9 @@ SUBROUTINE iosys()
refold_pos, remove_rigid_rot, upscale, &
pot_extrapolation, wfc_extrapolation, &
w_1, w_2, trust_radius_max, trust_radius_min, &
trust_radius_ini, bfgs_ndim
trust_radius_ini, bfgs_ndim, &
fire_nmin, fire_f_inc, fire_f_dec, &
fire_alpha_init, fire_falpha, fire_dtmax
!
! ... CELL namelist
!
@ -392,6 +400,20 @@ SUBROUTINE iosys()
!
ntcheck = nstep + 1
!
CASE ( 'fire' )
!
lmd = .true.
calc = 'fi'
! set fire variables
fire_nmin_ = fire_nmin
fire_f_inc_ = fire_f_inc
fire_f_dec_ = fire_f_dec
fire_alpha_init_ = fire_alpha_init
fire_falpha_ = fire_falpha
fire_dtmax_ = fire_dtmax
!
ntcheck = nstep + 1
!
CASE ( 'ipi' )
!
CONTINUE

View File

@ -45,8 +45,10 @@ SUBROUTINE move_ions( idone, ions_status, optimizer_failed )
USE mp, ONLY : mp_bcast
USE bfgs_module, ONLY : bfgs, terminate_bfgs
USE basic_algebra_routines, ONLY : norm
USE dynamics_module, ONLY : verlet, terminate_verlet, proj_verlet
USE dynamics_module, ONLY : smart_MC, langevin_md
USE dynamics_module, ONLY : verlet, terminate_verlet, proj_verlet, fire
USE dynamics_module, ONLY : smart_MC, langevin_md, dt
USE dynamics_module, ONLY : fire_nmin, fire_f_inc, fire_f_dec, &
fire_alpha_init, fire_falpha, fire_dtmax
USE klist, ONLY : nelec, tot_charge
USE fcp_module, ONLY : lfcp, fcp_eps, fcp_mu, fcp_relax, &
fcp_verlet, fcp_terminate, output_fcp
@ -299,6 +301,16 @@ SUBROUTINE move_ions( idone, ions_status, optimizer_failed )
!
ENDIF
!
ELSEIF ( calc == 'fi' ) THEN
!
CALL fire( conv_ions)
!
IF ( .NOT. conv_ions .AND. idone >= nstep ) THEN
WRITE( UNIT = stdout, FMT = &
'(/,5X,"The maximum number of steps has been reached.")' )
WRITE( UNIT = stdout, &
FMT = '(/,5X,"End of FIRE minimization")' )
ENDIF
ELSEIF ( calc(1:1) == 'l' ) THEN
!
! ... for smart monte carlo method