tdinit.f90

Go to the documentation of this file.

! Copyright (C) 2022 J. K. Dewhurst, S. Sharma and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine tdinit
use modmain
use modtddft
use moddftu
use modmpi
use modomp
use modramdisk
implicit none
! local variables
integer ik
! allocatable arrays
real(8), allocatable :: evalfv(:,:)
complex(8), allocatable :: evecfv(:,:,:),evecsv(:,:)
! flag for starting at t=0 or a restart
tdt0=((task == 460).or.(task == 462))
! determine the static density and charge and write to file
if (tdt0) call rhostatic
! store original parameters
tforce0=tforce
tfav00=tfav0
tjr0=tjr
! currents should be calculated with forces
if (tforce) tjr=.true.
! average force can be non-zero (allow for translation of atomic basis)
tfav0=.false.
! ensure eigenvectors are written to disk during initialisation
wrtdsk0=wrtdsk
wrtdsk=.true.
! initialise global variables
call init0
call init1
! read the charge density and potentials from file
call readstate
call genvsig
call gencore
call energykncr
call linengy
call genapwlofr
call gensocfr
if (tdt0) then
! generate eigenvalues and eigenvectors only at t=0 (not for the restart) for
! the k-point set reduced with the symmetries which leave A(t) invariant for all
! time steps
  call genevfsv
else
! only read in the second-variational eigenvalues for restarts
  call readevalsv
! read in the static density and charge
  call readrhos
end if
! compute the occupation numbers
call occupy
! DFT+U
if (dftu /= 0) then
  call gendmatmt
  call genvmatmt
  call vmatmtsc
  if (tmwrite) call genwkpr0
end if
! generate the kinetic matrix elements in the second-variational basis
call genkmat(.false.,.false.)
! write the momentum matrix elements in the second-variational basis
call genpmat
! copy EVALFV.OUT, EVECFV.OUT, OCCSV.OUT and EVECSV.OUT to _TD.OUT extension
if (tdt0) then
  allocate(evalfv(nstfv,nspnfv),evecfv(nmatmax,nstfv,nspnfv))
  allocate(evecsv(nstsv,nstsv))
  do ik=1,nkpt
! distribute among MPI processes
    if (mod(ik-1,np_mpi) /= lp_mpi) cycle
    call getevalfv('.OUT',ik,vkl(:,ik),evalfv)
    call putevalfv('_TD.OUT',ik,evalfv)
    call getevecfv('.OUT',ik,vkl(:,ik),vgkl(:,:,:,ik),evecfv)
    call putevecfv('_TD.OUT',ik,evecfv)
    call putoccsv('_TD.OUT',ik,occsv(:,ik))
    call getevecsv('.OUT',ik,vkl(:,ik),evecsv)
! randomise eigenvectors at t=0 if required
    call rndevsv(rndevt0,evecsv)
    call putevecsv('_TD.OUT',ik,evecsv)
  end do
  deallocate(evalfv,evecfv,evecsv)
end if
! set global file extension
filext='_TD.OUT'
! output the new k-point set to file
if (mp_mpi) call writekpts
! synchronise MPI processes
call mpi_barrier(mpicom,ierror)
! Ehrenfest dynamics
if ((task == 462).or.(task == 463)) then
! forces should not be calculated
  tforce=.false.
! enable small amplitude displacements
  tatdisp=.true.
! zero the displacements and velocities
  atdvc(:,:,:,:)=0.d0
! generate the gradient of the nucleus and static density Coulomb potential
  call gengvnsmt
end if
if (tdt0) then
! start from t=0
  itimes0=1
else
! restart if required
  call tdrestart
end if
! deallocate the static density if not required
if (.not.tjr) deallocate(rhosmt,rhosir)
! read the forces calculated during the previous TDDFT run
if (tatdisp) call readforcet
! read the atomic displacements and velocities
if (trdatdv) call readatdvc
wrtdsk=wrtdsk0
! make the time-dependent Kohn-Sham orbitals strictly orthogonal after each
! time step for complex time evolution
if (tdphi /= 0.d0) ntsorth=1

end subroutine