LCOV - code coverage report
Current view: top level - src/78_effpot - m_lwf_berendsen_mover.F90 (source / functions) Coverage Total Hit
Test: coverage.info Lines: 97.4 % 38 37
Test Date: 2026-09-19 17:42:43 Functions: 100.0 % 4 4

            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              : 
        

Generated by: LCOV version 2.3-1