TABLE OF CONTENTS


ABINIT/m_rttddft_propagate [ Modules ]

[ Top ] [ Modules ]

NAME

  m_rttddft_propagate

FUNCTION

  Driver to perform electronic or nuclear step

COPYRIGHT

  Copyright (C) 2021-2022 ABINIT group (FB)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module m_rttddft_propagate
23 
24  use defs_basis
25  use defs_abitypes,         only: MPI_type
26  use defs_datatypes,        only: pseudopotential_type
27  
28  use m_dtset,               only: dataset_type
29  use m_errors,              only: msg_hndl
30  use m_hamiltonian,         only: gs_hamiltonian_type
31  use m_rttddft,             only: rttddft_setup_ele_step
32  use m_rttddft_propagators, only: rttddft_propagator_er, &
33                                 & rttddft_propagator_emr
34  use m_rttddft_tdks,        only: tdks_type
35  use m_specialmsg,          only: wrtout
36 
37  implicit none
38 
39  private

m_rttddft/rttddft_propagate_ele [ Functions ]

[ Top ] [ m_rttddft ] [ Functions ]

NAME

  rttddft_propagate_ele

FUNCTION

  Main subroutine to propagate time-dependent KS orbitals

INPUTS

  dtset <type(dataset_type)>=all input variables for this dataset
  istep <integer> = step number
  mpi_enreg <MPI_type> = MPI-parallelisation information
  psps <type(pseudopotential_type)> = variables related to pseudopotentials
  tdks <type(tdks_type)> = the tdks object to initialize

OUTPUT

SOURCE

 67 subroutine rttddft_propagate_ele(dtset, istep, mpi_enreg, psps, tdks)
 68 
 69  implicit none
 70 
 71  !Arguments ------------------------------------
 72  !scalars
 73  integer,                    intent(in)    :: istep
 74  type(dataset_type),         intent(inout) :: dtset
 75  type(MPI_type),             intent(inout) :: mpi_enreg
 76  type(pseudopotential_type), intent(inout) :: psps
 77  type(tdks_type),            intent(inout) :: tdks
 78  
 79  !Local variables-------------------------------
 80  !scalars
 81  character(len=500)        :: msg
 82  type(gs_hamiltonian_type) :: gs_hamk
 83  !arrays
 84  
 85 ! ***********************************************************************
 86 
 87  write(msg,'(a,a,i0)') ch10, '--- Iteration', istep
 88  call wrtout(ab_out,msg)
 89  if (do_write_log) call wrtout(std_out,msg)
 90 
 91  ! Update various quantities after a nuclear step
 92  if (dtset%ionmov /= 0) call rttddft_setup_ele_step(dtset,mpi_enreg,psps,tdks)
 93 
 94  ! Propagate cg
 95  select case (dtset%td_propagator) 
 96    case(0)
 97       call rttddft_propagator_er(dtset,gs_hamk,istep,mpi_enreg,psps,tdks,calc_properties=.true.)
 98    case(1)
 99       call rttddft_propagator_emr(dtset,gs_hamk,istep,mpi_enreg,psps,tdks)  
100    case default
101       write(msg,"(a)") "Unknown Propagator - check the value of td_propagator"
102       ABI_ERROR(msg)
103  end select
104 
105 end subroutine rttddft_propagate_ele

m_rttddft/rttddft_propagate_nuc [ Functions ]

[ Top ] [ m_rttddft ] [ Functions ]

NAME

  rttddft_propagate_nuc

FUNCTION

  Main subroutine to propagate nuclei using Ehrenfest dynamics

INPUTS

  dtset <type(dataset_type)>=all input variables for this dataset
  istep <integer> = step number
  mpi_enreg <MPI_type> = MPI-parallelisation information
  psps <type(pseudopotential_type)> = variables related to pseudopotentials
  tdks <type(tdks_type)> = the tdks object to initialize

OUTPUT

SIDE EFFECTS

SOURCE

128 !subroutine rttddft_propagate_nuc(dtset, istep, mpi_enreg, psps, tdks)
129 
130 !implicit none
131 
132 !!Arguments ------------------------------------
133 !!scalars
134 !type(tdks_type),            intent(inout) :: tdks
135 !integer,                    intent(in)    :: istep
136 !type(dataset_type),         intent(in)    :: dtset
137 !type(MPI_type),             intent(inout) :: mpi_enreg
138 !type(pseudopotential_type), intent(inout) :: psps
139 !
140 !!Local variables-------------------------------
141 !!scalars
142 !character(len=500)   :: msg
143 !!arrays
144 !
145 ! ***********************************************************************
146 
147 !write(msg,'(2a,i5,a)') ch10,'--- Iteration',istep,ch10
148 !call wrtout(ab_out,msg)
149 !if (do_write_log) call wrtout(std_out,msg)
150 
151 !end subroutine rttddft_propagate_nuc