Line data Source code
1 : !!****m* ABINIT/m_lwf_berendsen_mover
2 : !! NAME
3 : !! m_lwf_berendsen_mover
4 : !!
5 : !! FUNCTION
6 : !! This module contains the lwf berensen NVT mover
7 : !!
8 : !!
9 : !! Datatypes:
10 : !!
11 : !! * lwf_berendsen_mover_t
12 : !!
13 : !! Subroutines:
14 : !!
15 : !! * TODO: update this when F2003 documentation format decided.
16 : !!
17 : !!
18 : !! COPYRIGHT
19 : !! Copyright (C) 2001-2026 ABINIT group (hexu)
20 : !! This file is distributed under the terms of the
21 : !! GNU General Public License, see ~abinit/COPYING
22 : !! or http://www.gnu.org/copyleft/gpl.txt .
23 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
24 : !!
25 : !! SOURCE
26 :
27 :
28 : #if defined HAVE_CONFIG_H
29 : #include "config.h"
30 : #endif
31 :
32 : #include "abi_common.h"
33 :
34 : module m_lwf_berendsen_mover
35 : use defs_basis
36 : use m_errors
37 : use m_abicore
38 : use m_xmpi
39 : use m_nctk
40 : use netcdf
41 : use m_mpi_scheduler, only: mpi_scheduler_t, init_mpi_info
42 : use m_multibinit_dataset, only: multibinit_dtset_type
43 : use m_random_xoroshiro128plus, only: set_seed, rand_normal_array, rng_t
44 : use m_abstract_potential, only: abstract_potential_t
45 : use m_abstract_mover, only: abstract_mover_t
46 : use m_hashtable_strval, only: hash_table_t
47 : use m_multibinit_cell, only: mbcell_t, mbsupercell_t
48 : use m_lwf_hist, only: lwf_hist_t
49 : use m_lwf_observables, only: lwf_observables_t
50 : use m_lwf_ncfile, only: lwf_ncfile_t
51 : use m_lwf_mover, only: lwf_mover_t
52 :
53 : implicit none
54 : private
55 : !!***
56 :
57 : type, public, extends(lwf_mover_t) :: lwf_berendsen_mover_t
58 : real(dp) :: taut ! the characteristic time of the relaxation of velocity.
59 : contains
60 : procedure :: set_params
61 : procedure :: scale_velocities
62 : procedure :: run_one_step
63 : end type lwf_berendsen_mover_t
64 :
65 : contains
66 :
67 1 : subroutine set_params(self, params)
68 : class(lwf_berendsen_mover_t), intent(inout) :: self
69 : type(multibinit_dtset_type) :: params
70 1 : call self%lwf_mover_t%set_params(params)
71 1 : self%taut=params%lwf_taut
72 1 : end subroutine set_params
73 :
74 :
75 21000 : subroutine run_one_step(self, effpot, displacement, strain, spin, lwf, energy_table)
76 : class(lwf_berendsen_mover_t), intent(inout) :: self
77 : real(dp), optional, intent(inout) :: displacement(:,:), strain(:,:), spin(:,:), lwf(:)
78 : class(abstract_potential_t), intent(inout) :: effpot
79 : type(hash_table_t),optional, intent(inout) :: energy_table
80 : integer :: i
81 : character(len=40) :: key
82 21000 : ABI_UNUSED_A(lwf)
83 : ! scale the velocity.
84 21000 : self%energy=0.0
85 21525000 : self%lwf_force(:) =0.0
86 : call effpot%calculate( displacement=displacement, strain=strain, &
87 : & spin=spin, lwf=self%lwf, lwf_force=self%lwf_force, &
88 84000 : & energy=self%energy, energy_table=energy_table)
89 :
90 21525000 : do i=1, self%nlwf
91 : self%vcart(i) = self%vcart(i) + &
92 21525000 : & (0.5_dp * self%dt) * self%lwf_force(i)/self%lwf_masses(i)
93 : end do
94 21000 : call self%scale_velocities()
95 21546000 : self%lwf= self%lwf+self%vcart * self%dt
96 21000 : call self%apply_constraints(self%lwf)
97 :
98 :
99 21000 : self%energy=0.0
100 21525000 : self%lwf_force(:)=0.0
101 : call effpot%calculate( displacement=displacement, strain=strain, &
102 : & spin=spin, lwf=self%lwf, lwf_force=self%lwf_force, &
103 84000 : & energy=self%energy, energy_table=energy_table)
104 : !call effpot%calculate( displacement=displacement, strain=strain, &
105 : ! & spin=spin, lwf=self%lwf, lwf_force=self%lwf_force, &
106 : ! & energy=self%energy, energy_table=energy_table)
107 21525000 : do i=1, self%nlwf
108 : self%vcart(i) = self%vcart(i) + &
109 21525000 : & (0.5_dp * self%dt) * self%lwf_force(i)/self%lwf_masses(i)
110 : end do
111 : !call self%force_stationary()
112 21000 : call self%scale_velocities()
113 21546000 : self%lwf= self%lwf+self%vcart * self%dt
114 21000 : call self%apply_constraints(self%lwf)
115 21000 : call self%get_T_and_Ek()
116 :
117 21000 : if (present(energy_table)) then
118 21000 : key = 'Lwf kinetic energy'
119 21000 : call energy_table%put(key, self%Ek)
120 : end if
121 :
122 21000 : end subroutine run_one_step
123 :
124 :
125 : !-------------------------------------------------------------------!
126 : ! scale_velocities:
127 : ! scale the velocities so that they get close to the required temperture
128 : !
129 : !-------------------------------------------------------------------!
130 42000 : subroutine scale_velocities(self)
131 : class(lwf_berendsen_mover_t), intent(inout) :: self
132 : real(dp) :: tautscl, old_temperature, scale_temperature, tmp
133 42000 : tautscl = self%dt / self%taut
134 42000 : old_temperature=self%T_ob
135 : if (old_temperature< 1e-19) then
136 : old_temperature=1e-19
137 : end if
138 42000 : tmp=1.0 +(self%temperature / old_temperature - 1.0) * tautscl
139 42000 : if(tmp< 0.0) then
140 0 : ABI_ERROR("The time scale for the Berendsen algorithm should be at least larger than lwf_dt")
141 : else
142 42000 : scale_temperature=sqrt(tmp)
143 : end if
144 : ! Limit the velocity scaling to reasonable values
145 42000 : if( scale_temperature > 1.1) then
146 : scale_temperature = 1.1
147 : elseif (scale_temperature < 0.9) then
148 : scale_temperature = 0.9
149 : endif
150 43050000 : self%vcart = self%vcart * scale_temperature
151 42000 : end subroutine scale_velocities
152 :
153 :
154 1 : end module m_lwf_berendsen_mover
155 :
156 :
|