Branch data Line data Source code
1 : : !!****m*ABINIT/m_lwf_mover
2 : : !! NAME
3 : : !! m_lwf_mover
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains the lwf mover, which controls how the lattice wannier function move.
7 : : !!
8 : : !!
9 : : !! Datatypes:
10 : : !!
11 : : !! * lwf_mover_t
12 : : !!
13 : : !! Subroutines:
14 : : !!
15 : : !! * lwf_mover_t_initialize
16 : : !! * lwf_mover_t_run_one_step
17 : : !! * lwf_mover_t_run_time
18 : : !! * TODO: update this when F2003 documentation format decided.
19 : : !!
20 : : !!
21 : : !! COPYRIGHT
22 : : !! Copyright (C) 2001-2026 ABINIT group (hexu)
23 : : !! This file is distributed under the terms of the
24 : : !! GNU General Public License, see ~abinit/COPYING
25 : : !! or http://www.gnu.org/copyleft/gpl.txt .
26 : : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
27 : : !!
28 : : !! SOURCE
29 : :
30 : :
31 : : #if defined HAVE_CONFIG_H
32 : : #include "config.h"
33 : : #endif
34 : : #include "abi_common.h"
35 : : module m_lwf_mover
36 : : use defs_basis
37 : : use m_errors
38 : : use m_abicore
39 : : use m_xmpi
40 : : use m_nctk
41 : : use netcdf
42 : : use m_mpi_scheduler, only: mpi_scheduler_t, init_mpi_info
43 : : use m_multibinit_dataset, only: multibinit_dtset_type
44 : : use m_random_xoroshiro128plus, only: set_seed, rand_normal_array, rng_t
45 : : use m_abstract_potential, only: abstract_potential_t
46 : : use m_abstract_mover, only: abstract_mover_t
47 : : use m_hashtable_strval, only: hash_table_t
48 : : use m_multibinit_cell, only: mbcell_t, mbsupercell_t
49 : : use m_lwf_hist, only: lwf_hist_t
50 : : use m_lwf_observables, only: lwf_observables_t
51 : : use m_lwf_ncfile, only: lwf_ncfile_t
52 : :
53 : :
54 : : implicit none
55 : : private
56 : : !!***
57 : :
58 : : type, public, extends(abstract_mover_t):: lwf_mover_t
59 : : type(multibinit_dtset_type), pointer:: params
60 : : real(dp) :: energy
61 : : integer:: nlwf
62 : : real(dp), allocatable:: lwf(:), lwf_force(:), vcart(:)
63 : : type(lwf_ncfile_t):: ncfile
64 : : type(lwf_hist_t):: hist
65 : : real(dp), pointer:: lwf_masses(:) => null()
66 : : real(dp):: Ek = 0.0_dp ! kinetic energy
67 : : real(dp):: T_ob = 0.0_dp ! observed temperature
68 : :
69 : : ! constraints
70 : : integer:: n_fixed_lwf = 0
71 : : integer, allocatable:: fixed_lwf_ids(:)
72 : : real(dp), allocatable:: fixed_lwf_values(:)
73 : :
74 : : contains
75 : : procedure:: initialize
76 : : procedure:: finalize
77 : : procedure:: set_temperature
78 : : procedure:: set_params
79 : : procedure:: set_initial_state
80 : : procedure:: get_T_and_Ek
81 : : procedure:: run_one_step
82 : : procedure:: run_time
83 : : procedure:: run_varT
84 : : procedure:: prepare_ncfile
85 : : procedure:: set_ncfile_name
86 : : procedure:: read_hist_lwf_state
87 : : procedure:: read_lwf_constraints
88 : : procedure:: apply_fixed_lwf
89 : : procedure:: apply_constraints
90 : : end type lwf_mover_t
91 : :
92 : : contains
93 : :
94 : 1 : subroutine initialize(self, params, supercell, rng)
95 : : class(lwf_mover_t), intent(inout):: self
96 : : type(multibinit_dtset_type), target, intent(in):: params
97 : : type(mbsupercell_t), target, intent(in):: supercell
98 : : type(rng_t), target, intent(in):: rng
99 : 1 : self%params=>params
100 : 1 : self%supercell=>supercell
101 : 1 : self%label="LWF Mover"
102 : 1 : call self%set_params(params)
103 : 1 : call self%set_rng(rng)
104 : 1 : self%nlwf = self%supercell%lwf%nlwf
105 [ - + + - : 3 : ABI_MALLOC(self%lwf, (self%nlwf))
+ - + - -
+ - + -
+ ]
106 [ - + + - : 3 : ABI_MALLOC(self%vcart, (self%nlwf))
+ - + - -
+ - + -
+ ]
107 [ - + + - : 3 : ABI_MALLOC(self%lwf_force, (self%nlwf))
+ - + - -
+ - + -
+ ]
108 [ + + ]: 1025 : self%lwf(:) = 0.0_dp
109 [ + + ]: 1025 : self%lwf_force(:) = 0.0_dp
110 [ + + ]: 1025 : self%vcart(:) = 0.0_dp
111 : 1 : self%energy = 0.0_dp
112 : 1 : self%lwf_masses=>self%supercell%lwf%lwf_masses
113 : 1 : call self%hist%initialize(nlwf = self%nlwf, mxhist = 1)
114 : :
115 : : !call self%read_lwf_constraints()
116 : 1 : end subroutine initialize
117 : :
118 : :
119 : 1 : subroutine finalize(self)
120 : : class(lwf_mover_t), intent(inout):: self
121 : 1 : nullify(self%supercell)
122 : 1 : nullify(self%params)
123 [ + - ]: 1 : ABI_SFREE(self%lwf)
124 [ + - ]: 1 : ABI_SFREE(self%vcart)
125 [ + - ]: 1 : ABI_SFREE(self%lwf_force)
126 : 1 : nullify(self%lwf_masses)
127 : 1 : call self%hist%finalize()
128 : :
129 : : ! constraints
130 : 1 : self%n_fixed_lwf = 0
131 [ + - ]: 1 : ABI_SFREE(self%fixed_lwf_ids)
132 [ + - ]: 1 : ABI_SFREE(self%fixed_lwf_values)
133 : :
134 : : !call self%ncfile%finalize()
135 : 1 : end subroutine finalize
136 : :
137 : 1 : subroutine set_params(self, params)
138 : : class(lwf_mover_t), intent(inout):: self
139 : : type(multibinit_dtset_type):: params
140 : 1 : self%dt = params%lwf_dt
141 : 1 : self%total_time = params%lwf_ntime*params%lwf_dt
142 : 1 : self%temperature = params%lwf_temperature
143 : 1 : end subroutine set_params
144 : :
145 : 7 : subroutine set_temperature(self, temperature)
146 : : class(lwf_mover_t), intent(inout):: self
147 : : real(dp), intent(in):: temperature
148 : 7 : self%temperature = temperature
149 : 7 : end subroutine set_temperature
150 : :
151 : : !-------------------------------------------------------------------!
152 : : !get_temperature_and_kinetic_energy
153 : : ! Ek = 1/2 \sum m_i vi^2
154 : : ! T = 2 Ek/nlwf (in a.u.)
155 : : !-------------------------------------------------------------------!
156 : 21000 : subroutine get_T_and_Ek(self)
157 : : class(lwf_mover_t), intent(inout):: self
158 [ + + ]: 21525000 : self%Ek = 0.5*sum(self%lwf_masses * (self%vcart*self%vcart))
159 : 21000 : self%T_ob = 2.0*self%Ek/self%nlwf
160 : 21000 : end subroutine get_T_and_Ek
161 : :
162 : :
163 : :
164 : 0 : subroutine run_one_step(self, effpot, displacement, strain, spin, lwf, energy_table)
165 : : class(lwf_mover_t), intent(inout):: self
166 : : real(dp), optional, intent(inout):: displacement(:,:), strain(:,:), spin(:,:), lwf(:)
167 : : class(abstract_potential_t), intent(inout):: effpot
168 : : type(hash_table_t), optional, intent(inout):: energy_table
169 : 0 : ABI_UNUSED_A(self)
170 : 0 : ABI_UNUSED_A(effpot)
171 : 0 : ABI_UNUSED_A(displacement)
172 : 0 : ABI_UNUSED_A(strain)
173 : 0 : ABI_UNUSED_A(spin)
174 : 0 : ABI_UNUSED_A(lwf)
175 : 0 : ABI_UNUSED_A(energy_table)
176 : 0 : ABI_ERROR("run_one_step not implemented for this mover")
177 : 0 : end subroutine run_one_step
178 : :
179 : : !-------------------------------------------------------------------!
180 : : ! run from begining to end.
181 : : !-------------------------------------------------------------------!
182 [ - + - + : 7 : subroutine run_time(self, effpot, displacement, strain, spin, lwf, energy_table)
- - - - ]
183 : : ! run one step. (For MC also?)
184 : : class(lwf_mover_t), intent(inout):: self
185 : : ! array of effective potentials so that there can be multiple of them.
186 : : class(abstract_potential_t), intent(inout):: effpot
187 : : real(dp), optional, intent(inout):: displacement(:,:), strain(:,:), spin(:,:), lwf(:)
188 : : type(hash_table_t), optional, intent(inout):: energy_table
189 : : integer:: i, nstep
190 : : character(len = 90):: msg
191 [ - + ]: 7 : if(present(lwf)) then
192 : 0 : ABI_ERROR("lwf should not be input for lwf mover")
193 : : end if
194 : 7 : ABI_UNUSED_A(displacement)
195 : 7 : ABI_UNUSED_A(strain)
196 : 7 : ABI_UNUSED_A(spin)
197 : 7 : ABI_UNUSED_A(energy_table)
198 : :
199 : :
200 : 7 : msg = repeat("=", 90)
201 : 7 : call wrtout(std_out, msg, 'COLL')
202 : 7 : call wrtout(ab_out, msg, 'COLL')
203 : 7 : write(msg, '(A22)') "LWF dynamic steps:"
204 : 7 : call wrtout(std_out, msg, 'COLL')
205 : 7 : call wrtout(ab_out, msg, 'COLL')
206 : 7 : msg = repeat("=", 90)
207 : 7 : call wrtout(std_out, msg, 'COLL')
208 : 7 : call wrtout(ab_out, msg, 'COLL')
209 : :
210 : : write(msg, "(A13, 4X, A15, 4X, A15, 4X, A15, 4X, A15)") &
211 : 7 : & "Iteration", "temperature(K)", "Ekin(Ha/uc)", &
212 : 14 : & "Epot(Ha/uc)", "ETOT(Ha/uc)"
213 : 7 : call wrtout(std_out, msg, 'COLL')
214 : 7 : call wrtout(ab_out, msg, 'COLL')
215 : :
216 : :
217 : 7 : call self%apply_constraints(self%lwf)
218 : :
219 [ - + ]: 7 : nstep = floor(self%total_time/self%dt)
220 [ + + ]: 21007 : do i = 1, nstep
221 : : !print *, "Step: ", i, " T: ", self%T_ob*Ha_K, " Ek:", self%Ek, "Ev", self%energy, "Etot", self%energy+self%Ek
222 [ + - ]: 42000 : call self%run_one_step(effpot = effpot, spin = spin, lwf = self%lwf, energy_table = energy_table)
223 : :
224 : :
225 : :
226 : :
227 : 21000 : call self%hist%set_hist(lwf = self%lwf, vcart = self%vcart, energy = self%energy )
228 [ + + ]: 21007 : if(modulo(i, self%params%lwf_nctime)==0) then
229 : 210 : call self%ncfile%write_one_step(self%hist)
230 : : !print *, "Step: ", i, "Ev", self%energy, "Etot"
231 : :
232 : 210 : write(msg, "(I13, 4X, F15.5, 4X, ES15.5, 4X, ES15.5, 4X, ES15.5)") i, self%T_ob*Ha_K, &
233 : 210 : & self%Ek/self%supercell%ncell, self%energy/self%supercell%ncell, &
234 : 420 : & (self%Ek+self%energy)/self%supercell%ncell
235 : 210 : call wrtout(std_out, msg, 'COLL')
236 : 210 : call wrtout(ab_out, msg, 'COLL')
237 : : !write(msg, "(I13, 4X, ES15.5)") i, self%energy/self%supercell%ncell
238 : : ! & (self%Ek+self%energy)/self%supercell%ncell
239 : :
240 : : end if
241 : : !TODO: output, observables
242 : : end do
243 : :
244 : 7 : msg = repeat("=", 90)
245 : 7 : call wrtout(std_out, msg, 'COLL')
246 : 7 : call wrtout(ab_out, msg, 'COLL')
247 : :
248 : 7 : end subroutine run_time
249 : :
250 : :
251 : : !----------------------------------------------------------------------
252 : : !> @brief set initial state.
253 : : !>
254 : : !> @param[in] mode: a integer to define the kind of initial state.
255 : : !----------------------------------------------------------------------
256 : 2 : subroutine set_initial_state(self, mode)
257 : : ! set initial positions, spin, etc
258 : : class(lwf_mover_t), intent(inout):: self
259 : : integer, optional, intent(in):: mode
260 : : integer:: i
261 : : real(dp):: tmp
262 : : real(dp):: kpoint(3)
263 : :
264 [ + + ]: 2050 : self%lwf(:)=0.0
265 : 2 : select case(mode)
266 : : case(0) ! using k-vector
267 : : ! TODO make it more generic
268 : 0 : kpoint(:)=[0.5_dp, 0.0_dp, 0.5_dp]
269 [ # # ]: 0 : do i = 1, self%supercell%ncell
270 : : tmp = 0.2*real(exp(cmplx(0.0, two_pi, kind = dp) * &
271 [ # # ]: 0 : &dot_product(kpoint, self%supercell%supercell_maker%rvecs(:, i))), kind = dp)
272 : 0 : self%lwf(i*2-1)=tmp
273 : 0 : self%lwf(i*2)=tmp
274 : : enddo
275 : : ! random number between-0.1 0.1
276 : : case(1)
277 : 2 : call self%rng%rand_unif_01_array(self%lwf, self%nlwf)
278 [ + + ]: 2050 : self%lwf=(self%lwf-0.5)*0.1
279 : : ! zero
280 : : case(2)
281 [ # # ]: 0 : self%lwf(:)=0.0
282 : : ! read from lwf hist file
283 : : case(4)
284 [ - + - - : 2 : call self%read_hist_lwf_state(self%params%lwf_init_hist_fname)
- ]
285 : : end select
286 : :
287 : 2 : call self%rng%rand_normal_array(self%vcart(:), self%nlwf)
288 [ + + ]: 2050 : do i = 1, self%nlwf
289 : 2050 : self%vcart(i) = self%vcart(i) *sqrt(self%temperature/self%lwf_masses(i))
290 : : end do
291 : :
292 : 2 : call self%hist%set_hist(lwf = self%lwf, vcart = self%vcart, energy = 0.0_dp)
293 : :
294 : 2 : end subroutine set_initial_state
295 : :
296 : : !-------------------------------------------------------------------!
297 : : ! read_hist_lwf_state
298 : : ! read the last step of spin from hist file.
299 : : !-------------------------------------------------------------------!
300 : 0 : subroutine read_hist_lwf_state(self, fname)
301 : : class(lwf_mover_t), intent(inout):: self
302 : : character(len=*), intent(in):: fname
303 : : integer:: ierr, ncid, varid
304 : : integer:: nlwf, ntime
305 : : character(len = 118):: msg
306 : : ! open file
307 : :
308 : 0 : ierr = nf90_open(trim(fname), NF90_NOWRITE, ncid)
309 [ # # # # ]: 0 : NCF_CHECK_MSG(ierr, "The lwf_init_state is set to 4. But opening netcdf file "//trim(fname)//" Failed. ")
310 : :
311 : : ! sanity check. If the hist file is consistent with the current calculation
312 : 0 : ierr = nctk_get_dim(ncid, "nlwf" , nlwf)
313 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading nlwf")
314 : :
315 : : msg="The number of lwfs in histfile is not equal & & to the present calculation." // &
316 : 0 : & " Please check if the file is consistent."
317 [ # # ]: 0 : if (nlwf /= self%nlwf) then
318 : 0 : ABI_ERROR(msg)
319 : : end if
320 : :
321 : 0 : ierr = nctk_get_dim(ncid, "ntime", ntime)
322 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading ntime")
323 : :
324 : : ! read lwf and set as initial state
325 : 0 : ierr = nf90_inq_varid(ncid, "lwf", varid)
326 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading lwf.")
327 : :
328 : : ierr = nf90_get_var(ncid = ncid, varid = varid, values = self%lwf(:), &
329 [ # # # # ]: 0 : & start=(/ 1, ntime/), count=(/nlwf, 1/))
330 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading lwf from lwf hist file")
331 : :
332 : : ! close file
333 : 0 : ierr = nf90_close(ncid)
334 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "Close netcdf file")
335 : :
336 : 0 : end subroutine read_hist_lwf_state
337 : :
338 : :
339 : : !-------------------------------------------------------------------!
340 : : ! read_lwf_constrain
341 : : ! read the constrains from file
342 : : !-------------------------------------------------------------------!
343 : 1 : subroutine read_lwf_constraints(self, fname)
344 : : class(lwf_mover_t), intent(inout):: self
345 : : character(len=*), intent(in):: fname
346 : : integer:: ierr, ncid, id_id, value_id
347 : : ! open file
348 [ + - ]: 1 : if (self%params%lwf_constraint .eq. 0) then
349 : 1 : self%n_fixed_lwf = 0
350 [ - + - + ]: 1 : ABI_MALLOC(self%fixed_lwf_ids, (self%n_fixed_lwf))
351 [ - + - + ]: 1 : ABI_MALLOC(self%fixed_lwf_values, (self%n_fixed_lwf))
352 : : else
353 : :
354 : 0 : ierr = nf90_open(trim(fname), NF90_NOWRITE, ncid)
355 [ # # # # ]: 0 : NCF_CHECK_MSG(ierr, "Trying to read constrain from netcdf file "//trim(fname)//" Failed. ")
356 : :
357 : : ! sanity check. If the hist file is consistent with the current calculation
358 : 0 : ierr = nctk_get_dim(ncid, "n_fixed_lwf" , self%n_fixed_lwf)
359 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading n_fixed_lwf")
360 : :
361 [ # # ]: 0 : if (self%n_fixed_lwf .ne. 0) then
362 [ # # # # : 0 : ABI_MALLOC(self%fixed_lwf_ids, (self%n_fixed_lwf))
# # # # #
# # # #
# ]
363 [ # # # # : 0 : ABI_MALLOC(self%fixed_lwf_values, (self%n_fixed_lwf))
# # # # #
# # # #
# ]
364 : :
365 : :
366 : : ! read fixed lwf ids and values
367 : 0 : ierr = nf90_inq_varid(ncid, "fixed_lwf_ids", id_id)
368 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading fixed_lwf_ids.")
369 : :
370 : : ierr = nf90_get_var(ncid = ncid, varid = id_id, values = self%fixed_lwf_ids, &
371 [ # # ]: 0 : & start=[1], count=[self%n_fixed_lwf])
372 [ # # # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading fixed_lwf_ids from file "//trim(fname)//". " )
373 : :
374 : 0 : ierr = nf90_inq_varid(ncid, "fixed_lwf_values", value_id)
375 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading fixed_lwf_values.")
376 : :
377 : : ierr = nf90_get_var(ncid = ncid, varid = value_id, values = self%fixed_lwf_values, &
378 [ # # ]: 0 : & start=[1], count=[self%n_fixed_lwf])
379 [ # # # # ]: 0 : NCF_CHECK_MSG(ierr, "when reading fixed_lwf_values from file "//trim(fname)//". " )
380 : : endif
381 : :
382 : :
383 : : ! close file
384 : 0 : ierr = nf90_close(ncid)
385 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "Close netcdf file")
386 : : endif
387 : :
388 : 1 : end subroutine read_lwf_constraints
389 : :
390 : :
391 : :
392 : 7 : subroutine prepare_ncfile(self, params, fname)
393 : : class(lwf_mover_t), intent(inout):: self
394 : : type(multibinit_dtset_type):: params
395 : : character(len=*), intent(in):: fname
396 : : integer:: master, my_rank, comm, nproc
397 : : logical:: iam_master
398 : 7 : ABI_UNUSED_A(params)
399 : 7 : call init_mpi_info(master, iam_master, my_rank, comm, nproc)
400 [ + - ]: 7 : if(iam_master) then
401 : 7 : call self%ncfile%initialize(fname, 1)
402 : 7 : call self%ncfile%write_cell(self%supercell)
403 : 7 : call self%ncfile%def_lwf_var(self%hist)
404 : : end if
405 : 7 : end subroutine prepare_ncfile
406 : :
407 : : !-------------------------------------------------------------------!
408 : : !set_ncfile_name :
409 : : !-------------------------------------------------------------------!
410 : 0 : subroutine set_ncfile_name(self, params, fname)
411 : : class(lwf_mover_t), intent(inout):: self
412 : : type(multibinit_dtset_type):: params
413 : : character(len = fnlen), intent(in):: fname
414 : : integer:: master, my_rank, comm, nproc
415 : : logical:: iam_master
416 : 0 : call init_mpi_info(master, iam_master, my_rank, comm, nproc)
417 [ # # ]: 0 : if (iam_master) then
418 [ # # ]: 0 : call self%prepare_ncfile(params, trim(fname)//'_lwfhist.nc')
419 : 0 : call self%ncfile%write_one_step(self%hist)
420 : : endif
421 : 0 : end subroutine set_ncfile_name
422 : :
423 : :
424 : :
425 : : !!****f*m_lwf_mover/run_varT
426 : : !!
427 : : !! NAME
428 : : !! run_varT
429 : : !!
430 : : !! FUNCTION
431 : : !! run M vs Temperature
432 : : !!
433 : : !! INPUTS
434 : : !! pot: potential
435 : : !! T_start, Tend, T_nstep
436 : : !u
437 : : !! OUTPUT
438 : : !!
439 : : !! SOURCE
440 [ - + - - : 1 : subroutine run_varT(self, pot, ncfile_prefix, displacement, strain, spin, lwf, energy_table)
- - - + -
- - - - +
- - - - -
+ - - -
- ]
441 : : class(lwf_mover_t), intent(inout):: self
442 : : class(abstract_potential_t), intent(inout):: pot
443 : : real(dp), optional, intent(inout):: displacement(:,:), strain(:,:), lwf(:), spin(:,:)
444 : : character(fnlen), intent(inout):: ncfile_prefix
445 : : type(hash_table_t), optional, intent(inout):: energy_table
446 : : real(dp):: T_start, T_end
447 : : integer:: T_nstep
448 : : !type(lwf_ncfile_t):: lwf_ncfile
449 : : character(len = 4):: post_fname
450 : : real(dp):: T, T_step
451 : : integer:: i
452 : : !integer:: Tfile, iostat
453 : : character(len = 90):: msg
454 : : !character(len = 4200):: Tmsg ! to write to var T file
455 : : !character(len = 150):: iomsg
456 : : !character(fnlen):: Tfname ! file name for output various T calculation
457 : : !real(dp), allocatable:: Tlist(:), chi_list(:), Cv_list(:), binderU4_list(:)
458 : : !real(dp), allocatable:: Mst_sub_norm_list(:, :)
459 : : !real(dp), allocatable :: Mst_norm_total_list(:)
460 : :
461 : : integer:: master, my_rank, comm, nproc, ierr
462 : : logical:: iam_master
463 : 1 : call init_mpi_info(master, iam_master, my_rank, comm, nproc)
464 : :
465 [ + - ]: 1 : if (iam_master) then
466 : 1 : T_start = self%params%lwf_temperature_start
467 : 1 : T_end = self%params%lwf_temperature_end
468 : 1 : T_nstep = self%params%lwf_temperature_nstep
469 : : !Tfile = get_unit()
470 : : !Tfname = trim(ncfile_prefix)//'.varT'
471 : : !iostat = open_file(file = Tfname, unit = Tfile, iomsg = iomsg )
472 [ + - ]: 1 : if (T_nstep <= 1) then
473 : : T_step = 0.0
474 : : else
475 : 1 : T_step=(T_end-T_start)/(T_nstep-1)
476 : : endif
477 : : write(msg, "(A52, ES13.5, A11, ES13.5, A1)") &
478 : 1 : & "Starting temperature dependent calculations. T from ", &
479 : 2 : & T_start*Ha_K, "K to ", T_end*Ha_K, " K."
480 : 1 : call wrtout(std_out, msg, "COLL")
481 : 1 : call wrtout(ab_out, msg, "COLL")
482 : : end if
483 : :
484 : 1 : call xmpi_bcast(T_nstep, 0, comm, ierr)
485 [ + + ]: 8 : do i = 1, T_nstep
486 [ + - ]: 7 : if(iam_master) then
487 : 7 : T = T_start+(i-1)*T_step
488 : 7 : msg = repeat("=", 79)
489 : 7 : call wrtout(std_out, msg, "COLL")
490 : 7 : call wrtout(ab_out, msg, "COLL")
491 : :
492 : 7 : write(msg, "(A13, 5X, ES13.5, A3)") "Temperature: ", T*Ha_K, " K."
493 : 7 : call wrtout(std_out, msg, "COLL")
494 : 7 : call wrtout(ab_out, msg, "COLL")
495 : :
496 : 7 : call self%hist%reset(array_to_zero=.False.)
497 : : ! set temperature
498 : : ! TODO make this into a subroutine set_params
499 : : endif
500 : 7 : call self%set_temperature(temperature = T)
501 [ + - ]: 7 : if(iam_master) then
502 [ + + ]: 7 : if(i == 1) then
503 : 1 : call self%set_initial_state(mode = self%params%lwf_init_state)
504 : : endif
505 : 7 : call self%apply_constraints(self%lwf)
506 : :
507 : 7 : write(post_fname, "(I4.4)") i
508 : : call self%prepare_ncfile( self%params, &
509 [ + - ]: 7 : & trim(ncfile_prefix)//'_T'//post_fname//'_lwfhist.nc')
510 : 7 : call self%ncfile%write_one_step(self%hist)
511 : : endif
512 : :
513 : : call self%run_time(pot, displacement = displacement, strain = strain, spin = spin, &
514 [ + - + - : 35 : & lwf = lwf, energy_table = energy_table)
+ - + - ]
515 : :
516 [ + - ]: 8 : if(iam_master) then
517 : 7 : call self%ncfile%finalize()
518 : : endif
519 : : end do
520 : :
521 : 1 : end subroutine run_varT
522 : : !!***
523 : :
524 : :
525 [ - + ]: 42014 : subroutine apply_fixed_lwf(self, lwf)
526 : : class(lwf_mover_t), intent(inout):: self
527 : : real(dp), intent(inout):: lwf(:)
528 : : integer:: i
529 [ - + ]: 42014 : do i = 1, self%n_fixed_lwf
530 : 0 : self%vcart(self%fixed_lwf_ids(i))=0.0_dp
531 : 42014 : lwf(self%fixed_lwf_ids(i))=self%fixed_lwf_values(i)
532 : : end do
533 : 42014 : end subroutine apply_fixed_lwf
534 : :
535 [ - + ]: 42014 : subroutine apply_constraints(self, lwf)
536 : : class(lwf_mover_t), intent(inout):: self
537 : : real(dp), intent(inout):: lwf(:)
538 : 42014 : call self%apply_fixed_lwf(lwf)
539 : 42014 : end subroutine apply_constraints
540 : :
541 [ - + - - : 3 : end module m_lwf_mover
+ + - + +
- - + - +
+ - - + -
+ - + - +
- + - + -
+ - + ]
542 : :
|