49 use,
intrinsic :: iso_fortran_env
137 type(space_t),
private :: space
138 type(states_elec_dim_t) :: d
139 type(hamiltonian_elec_base_t) :: hm_base
140 type(phase_t) :: phase
141 type(energy_t),
allocatable :: energy
142 type(absorbing_boundaries_t) :: abs_boundaries
143 type(ks_potential_t) :: ks_pot
144 real(real64),
allocatable :: vberry(:,:)
146 type(derivatives_t),
pointer,
private :: der
148 type(nonlocal_pseudopotential_t) :: vnl
150 type(ions_t),
pointer :: ions
151 real(real64) :: exx_coef
153 type(poisson_t) :: psolver
156 logical :: self_induced_magnetic
157 real(real64),
allocatable :: a_ind(:, :)
158 real(real64),
allocatable :: b_ind(:, :)
160 integer :: theory_level
161 type(xc_t),
pointer :: xc
162 type(xc_photons_t),
pointer :: xc_photons
168 logical,
private :: adjoint
171 real(real64),
private :: mass
174 logical,
private :: inh_term
175 type(states_elec_t) :: inh_st
179 type(oct_exchange_t) :: oct_exchange
181 type(scissor_t) :: scissor
183 real(real64) :: current_time
184 logical,
private :: is_applied_packed
187 type(lda_u_t) :: lda_u
188 integer :: lda_u_level
190 logical,
public :: time_zero
192 type(exchange_operator_t),
public :: exxop
194 type(kpoints_t),
pointer,
public :: kpoints => null()
196 type(partner_list_t) :: external_potentials
197 real(real64),
allocatable,
public :: v_ext_pot(:)
198 real(real64),
allocatable,
public :: v_static(:)
200 type(ion_electron_local_potential_t) :: v_ie_loc
203 type(magnetic_constrain_t) :: magnetic_constrain
209 type(mxll_coupling_t) :: mxll
210 type(zora_t),
pointer :: zora => null()
225 integer,
public,
parameter :: &
234 mc, kpoints, need_exchange, xc_photons)
238 type(
grid_t),
target,
intent(inout) :: gr
239 type(
ions_t),
target,
intent(inout) :: ions
242 integer,
intent(in) :: theory_level
243 type(
xc_t),
target,
intent(in) :: xc
246 logical,
optional,
intent(in) :: need_exchange
247 type(
xc_photons_t),
optional,
target,
intent(in) :: xc_photons
250 logical :: need_exchange_
251 real(real64) :: rashba_coupling
258 hm%theory_level = theory_level
261 hm%kpoints => kpoints
287 if (space%dim /= 2)
then
288 write(
message(1),
'(a)')
'Rashba spin-orbit coupling can only be used for two-dimensional systems.'
294 call hm%hm_base%init(hm%d%nspin, hm%mass, rashba_coupling)
297 assert(
associated(gr%der%lapl))
298 hm%hm_base%kinetic => gr%der%lapl
300 safe_allocate(hm%energy)
307 if(
present(xc_photons))
then
308 hm%xc_photons => xc_photons
310 hm%xc_photons => null()
319 call poisson_init(hm%psolver, namespace, space, gr%der, mc, gr%stencil, st%qtot)
322 call epot_init(hm%ep, namespace, gr, hm%ions, hm%psolver, hm%d%ispin, hm%xc%family, hm%kpoints)
323 call kick_init(hm%kick, namespace, space, hm%kpoints, hm%d%ispin)
325 hm%zora =>
zora_t(namespace, hm%der, hm%d, hm%ep, hm%mass)
328 call hm%v_ie_loc%init(gr, hm%psolver, hm%ions, namespace)
330 call hm%nlcc%init(gr, hm%ions)
331 safe_allocate(st%rho_core(1:gr%np))
364 call parse_variable(namespace,
'CalculateSelfInducedMagneticField', .false., hm%self_induced_magnetic)
365 if (hm%self_induced_magnetic)
then
366 safe_allocate(hm%a_ind(1:gr%np_part, 1:space%dim))
367 safe_allocate(hm%b_ind(1:gr%np_part, 1:space%dim))
375 hm%inh_term = .false.
380 call hm%phase%init(gr, hm%d%kpt, hm%kpoints, st%d, space)
400 call lda_u_init(hm%lda_u, namespace, space, hm%lda_u_level, gr, ions, st, mc, &
401 hm%kpoints, hm%phase%is_allocated())
410 if(hm%lda_u_level /=
dft_u_none .and. hm%phase%is_allocated())
then
425 call parse_variable(namespace,
'HamiltonianApplyPacked', .
true., hm%is_applied_packed)
427 if (hm%theory_level ==
hartree_fock .and. st%parallel_in_states)
then
432 .and. st%parallel_in_states)
then
451 if (hm%xc%compute_exchange(hm%theory_level) .or. need_exchange_)
then
456 else if (
bitand(hm%xc%family, xc_family_oep) /= 0 .or. hm%theory_level ==
rdmft)
then
458 hm%kpoints, hm%xc%cam)
459 if (hm%theory_level ==
rdmft) hm%exxop%useACE = .false.
468 if (gr%use_curvilinear)
then
470 hm%is_applied_packed = .false.
471 call messages_write(
'Cannot use CUDA or OpenCL as curvilinear coordinates are used.')
474 call messages_write(
'Cannot use CUDA or OpenCL as curvilinear coordinates are used.', new_line = .
true.)
475 call messages_write(
'Calculation will not be continued. To force execution, set AllowCPUonly = yes.')
505 if (
associated(hm%xc_photons))
then
506 if (hm%xc_photons%wants_to_renormalize_mass())
then
508 call hm%set_mass(namespace, hm%xc_photons%get_renormalized_mass())
512 if (hm%xc%compute_exchange(hm%theory_level) .or. need_exchange_)
call hm%exxop%write_info(namespace)
522 class(*),
pointer :: potential
527 safe_allocate(hm%v_ext_pot(1:gr%np))
528 hm%v_ext_pot(1:gr%np) =
m_zero
530 call iter%start(hm%external_potentials)
531 do while (iter%has_next())
532 potential => iter%get_next()
533 select type (potential)
536 call potential%allocate_memory(gr)
537 call potential%calculate(namespace, gr, hm%psolver)
540 select case (potential%type)
546 message(1) =
"Cannot use static magnetic field with real wavefunctions"
550 if (.not.
allocated(hm%ep%b_field))
then
551 safe_allocate(hm%ep%b_field(1:3))
552 hm%ep%b_field(1:3) =
m_zero
554 hm%ep%b_field(1:3) = hm%ep%b_field(1:3) + potential%b_field(1:3)
556 if (.not.
allocated(hm%ep%a_static))
then
557 safe_allocate(hm%ep%a_static(1:gr%np, 1:space%dim))
558 hm%ep%a_static(1:gr%np, 1:space%dim) =
m_zero
560 call lalg_axpy(gr%np, space%dim,
m_one, potential%a_static, hm%ep%a_static)
563 if (.not.
allocated(hm%ep%e_field))
then
564 safe_allocate(hm%ep%e_field(1:space%dim))
565 hm%ep%e_field(1:space%dim) =
m_zero
567 hm%ep%e_field(1:space%dim) = hm%ep%e_field(1:space%dim) + potential%e_field(1:space%dim)
570 if (space%periodic_dim < space%dim)
then
571 if (.not.
allocated(hm%v_static))
then
572 safe_allocate(hm%v_static(1:gr%np))
573 hm%v_static(1:gr%np) =
m_zero
575 if (.not.
allocated(hm%ep%v_ext))
then
576 safe_allocate(hm%ep%v_ext(1:gr%np_part))
577 hm%ep%v_ext(1:gr%np_part) =
m_zero
583 if (hm%kpoints%use_symmetries)
then
587 message(1) =
"The StaticElectricField breaks (at least) one of the symmetries used to reduce the k-points."
588 message(2) =
"Set SymmetryBreakDir equal to StaticElectricField."
595 call potential%deallocate_memory()
608 class(*),
pointer :: potential
612 call iter%start(hm%external_potentials)
613 do while (iter%has_next())
614 potential => iter%get_next()
615 select type (potential)
619 message(1) =
"Applying StaticElectricField in a periodic direction is only accurate for large supercells."
620 message(2) =
"Single-point Berry phase is not appropriate when k-point sampling is needed."
636 logical :: external_potentials_present
637 logical :: kick_present
641 if (
allocated(hm%ep%e_field) .and. space%is_periodic() .and. .not.
list_has_gauge_field(ext_partners))
then
644 if (any(abs(hm%ep%e_field(1:space%periodic_dim)) >
m_epsilon))
then
645 safe_allocate(hm%vberry(1:gr%np, 1:hm%d%nspin))
656 call pcm_init(hm%pcm, namespace, space, ions, gr, st%qtot, st%val_charge, external_potentials_present, kick_present)
657 if (hm%pcm%run_pcm)
then
680 call hm%hm_base%end()
686 safe_deallocate_a(hm%vberry)
687 safe_deallocate_a(hm%a_ind)
688 safe_deallocate_a(hm%b_ind)
689 safe_deallocate_a(hm%v_ext_pot)
691 safe_deallocate_p(hm%zora)
710 safe_deallocate_a(hm%energy)
712 if (hm%pcm%run_pcm)
call pcm_end(hm%pcm)
714 call hm%v_ie_loc%end()
717 call iter%start(hm%external_potentials)
718 do while (iter%has_next())
719 potential => iter%get_next()
720 safe_deallocate_p(potential)
722 call hm%external_potentials%empty()
723 safe_deallocate_a(hm%v_static)
749 integer,
intent(in) :: terms
770 real(real64),
intent(in) ::
delta(:)
771 real(real64),
intent(in) :: emin
782 hm%spectral_middle_point = (emax + emin) /
m_two
783 hm%spectral_half_span = (emax - emin) /
m_two
818 if (hm%inh_term)
then
820 hm%inh_term = .false.
832 if (.not. hm%adjoint)
then
835 hm%abs_boundaries%mf = -hm%abs_boundaries%mf
852 hm%abs_boundaries%mf = -hm%abs_boundaries%mf
864 class(
mesh_t),
intent(in) :: mesh
866 class(
space_t),
intent(in) :: space
868 real(real64),
optional,
intent(in) :: time
870 integer :: ispin, ip, idir, iatom, ilaser
871 real(real64) :: aa(space%dim), time_
872 real(real64),
allocatable :: vp(:,:)
875 real(real64) :: am(space%dim)
880 this%current_time =
m_zero
881 if (
present(time)) this%current_time = time
886 call this%hm_base%clear(mesh%np)
893 if (
present(time) .or. this%time_zero)
then
896 if(
associated(lasers))
then
897 do ilaser = 1, lasers%no_lasers
898 select case (
laser_kind(lasers%lasers(ilaser)))
900 do ispin = 1, this%d%spin_channels
902 this%hm_base%potential(:, ispin), time_)
908 safe_allocate(vp(1:mesh%np, 1:space%dim))
909 vp(1:mesh%np, 1:space%dim) =
m_zero
913 do idir = 1, space%dim
914 this%hm_base%vector_potential(idir, ip) = this%hm_base%vector_potential(idir, ip) + vp(ip, idir)/
p_c
918 call laser_field(lasers%lasers(ilaser), this%hm_base%uniform_magnetic_field(1:space%dim), time_)
919 safe_deallocate_a(vp)
925 this%hm_base%uniform_vector_potential(1:space%dim) = this%hm_base%uniform_vector_potential(1:space%dim) - aa/
p_c
930 assert(
allocated(this%hm_base%uniform_vector_potential))
932 this%hm_base%uniform_vector_potential(1:space%dim) = this%hm_base%uniform_vector_potential(1:space%dim) - am/
p_c
938 if (
associated(gfield))
then
941 this%hm_base%uniform_vector_potential(1:space%dim) = this%hm_base%uniform_vector_potential(1:space%dim) - aa/
p_c
945 if (
allocated(this%ep%e_field) .and.
associated(gfield))
then
946 this%hm_base%uniform_vector_potential(1:space%periodic_dim) = &
947 this%hm_base%uniform_vector_potential(1:space%periodic_dim) - time_*this%ep%e_field(1:space%periodic_dim)
951 if (
associated(this%xc_photons))
then
952 if(this%xc_photons%lpfmf .and.
allocated(this%xc_photons%mf_vector_potential))
then
955 this%hm_base%uniform_vector_potential(1:space%dim) = &
956 this%hm_base%uniform_vector_potential(1:space%dim) - this%xc_photons%mf_vector_potential(1:space%dim)/
p_c
963 if (
allocated(this%ep%a_static))
then
968 do idir = 1, space%dim
969 this%hm_base%vector_potential(idir, ip) = this%hm_base%vector_potential(idir, ip) + this%ep%a_static(ip, idir)
978 call this%hm_base%accel_copy_pot(mesh)
981 if (
allocated(this%ep%b_field))
then
984 this%hm_base%uniform_magnetic_field(idir) = this%hm_base%uniform_magnetic_field(idir) + this%ep%b_field(idir)
989 call this%hm_base%update_magnetic_terms(mesh, this%ep%gyromagnetic_ratio, this%d%ispin)
994 if (this%mxll%test_equad)
then
1006 integer :: ik, imat, nmat, max_npoints, offset
1008 integer :: iphase, nphase
1012 if ((.not. this%kpoints%gamma_only()) .or.
allocated(this%hm_base%uniform_vector_potential))
then
1016 do iatom = 1, this%ep%natoms
1017 call projector_init_phases(this%ep%proj(iatom), space%dim, this%d, this%der%boundaries, this%kpoints, &
1018 vec_pot = this%hm_base%uniform_vector_potential, vec_pot_var = this%hm_base%vector_potential)
1024 if (
allocated(this%hm_base%uniform_vector_potential))
then
1026 call this%phase%update(mesh, this%d%kpt, this%kpoints, this%d, space, this%hm_base%uniform_vector_potential)
1031 vec_pot = this%hm_base%uniform_vector_potential, vec_pot_var = this%hm_base%vector_potential)
1035 max_npoints = this%vnl%max_npoints
1036 nmat = this%vnl%nprojector_matrices
1039 if (this%phase%is_allocated() .and.
allocated(this%vnl%projector_matrices))
then
1042 if (this%der%boundaries%spiralBC) nphase = 3
1044 if (.not.
allocated(this%vnl%projector_phases))
then
1045 safe_allocate(this%vnl%projector_phases(1:max_npoints, 1:nphase, nmat, this%d%kpt%start:this%d%kpt%end))
1048 type_cmplx, this%vnl%total_points*nphase*this%d%kpt%nlocal)
1051 this%vnl%nphase = nphase
1056 do ik = this%d%kpt%start, this%d%kpt%end
1057 do imat = 1, this%vnl%nprojector_matrices
1058 iatom = this%vnl%projector_to_atom(imat)
1059 do iphase = 1, nphase
1061 do ip = 1, this%vnl%projector_matrices(imat)%npoints
1062 this%vnl%projector_phases(ip, iphase, imat, ik) = this%ep%proj(iatom)%phase(ip, iphase, ik)
1065 if (
accel_is_enabled() .and. this%vnl%projector_matrices(imat)%npoints > 0)
then
1067 this%vnl%projector_matrices(imat)%npoints, this%vnl%projector_phases(1:, iphase, imat, ik), &
1068 offset = offset, async=.
true.)
1070 offset = offset + this%vnl%projector_matrices(imat)%npoints
1091 class(
mesh_t),
intent(in) :: mesh
1092 logical,
optional,
intent(in) :: accumulate
1094 integer :: ispin, ip
1101 do ispin = 1, this%d%nspin
1104 this%hm_base%potential(ip, ispin) =
m_zero
1111 do ispin = 1, this%d%nspin
1112 if (ispin <= 2)
then
1116 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + this%ep%vpsl(ip) + this%v_ext_pot(ip)
1120 if (this%pcm%run_pcm)
then
1121 if (this%pcm%solute)
then
1124 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + &
1125 this%pcm%v_e_rs(ip) + this%pcm%v_n_rs(ip)
1129 if (this%pcm%localf)
then
1132 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + &
1133 this%pcm%v_ext_rs(ip)
1143 this%hm_base%Impotential(ip, ispin) = this%hm_base%Impotential(ip, ispin) + this%abs_boundaries%mf(ip)
1154 call this%zora%update(this%der, this%hm_base%potential)
1158 do ispin = 1, this%d%nspin
1160 if (
allocated(this%hm_base%zeeman_pot))
then
1163 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + this%hm_base%zeeman_pot(ip, ispin)
1169 if (this%mxll%test_equad)
then
1172 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + this%mxll%e_quadrupole_pot(ip)
1181 call this%ks_pot%add_vhxc(this%hm_base%potential)
1183 call this%hm_base%accel_copy_pot(mesh)
1193 type(
grid_t),
intent(in) :: gr
1194 type(
ions_t),
target,
intent(inout) :: ions
1197 real(real64),
optional,
intent(in) :: time
1202 call epot_generate(this%ep, namespace, gr, this%ions, this%d)
1207 if (
allocated(this%ep%e_field) .and. space%periodic_dim < space%dim)
then
1213 this%v_ie_loc%atoms_dist => ions%atoms_dist
1214 this%v_ie_loc%atom => ions%atom
1215 call this%v_ie_loc%calculate()
1218 call lalg_axpy(gr%np,
m_one, this%v_ie_loc%potential(:,1), this%ep%vpsl)
1222 if (this%ep%nlcc)
then
1223 this%nlcc%atoms_dist => ions%atoms_dist
1224 this%nlcc%atom => ions%atom
1225 call this%nlcc%calculate()
1226 call lalg_copy(gr%np, this%nlcc%density(:,1), st%rho_core)
1229 call this%vnl%build(space, gr, this%ep)
1234 if (this%ep%non_local .and. .not. this%vnl%apply_projector_matrices)
then
1236 call messages_write(
'Relativistic pseudopotentials have not been fully implemented for CUDA or OpenCL.')
1239 call messages_write(
'Relativistic pseudopotentials have not been fully implemented for CUDA or OpenCL.',&
1241 call messages_write(
'Calculation will not be continued. To force execution, set AllowCPUonly = yes.')
1248 if (this%pcm%run_pcm)
then
1251 if (this%pcm%solute)
then
1258 if (this%pcm%localf .and.
allocated(this%v_static))
then
1259 call pcm_calc_pot_rs(this%pcm, gr, this%psolver, v_ext = this%ep%v_ext(1:gr%np_part))
1264 call lda_u_update_basis(this%lda_u, space, gr, ions, st, this%psolver, namespace, this%kpoints, &
1265 this%phase%is_allocated())
1275 time = this%current_time
1283 apply = this%is_applied_packed
1291 type(namespace_t),
intent(in) :: namespace
1292 class(space_t),
intent(in) :: space
1293 type(lattice_vectors_t),
intent(in) :: latt
1294 class(species_t),
intent(in) :: species
1295 real(real64),
intent(in) :: pos(1:space%dim)
1296 integer,
intent(in) :: ia
1297 class(mesh_t),
intent(in) :: mesh
1298 complex(real64),
intent(in) :: psi(:,:)
1299 complex(real64),
intent(out) :: vpsi(:,:)
1302 real(real64),
allocatable :: vlocal(:)
1305 safe_allocate(vlocal(1:mesh%np_part))
1307 call epot_local_potential(hm%ep, namespace, space, latt, mesh, species, pos, ia, vlocal)
1309 do idim = 1, hm%d%dim
1310 vpsi(1:mesh%np, idim) = vlocal(1:mesh%np) * psi(1:mesh%np, idim)
1313 safe_deallocate_a(vlocal)
1324 class(space_t),
intent(in) :: space
1325 class(mesh_t),
intent(in) :: mesh
1326 type(partner_list_t),
intent(in) :: ext_partners
1327 real(real64),
intent(in) :: time(1:2)
1328 real(real64),
intent(in) :: mu(1:2)
1330 integer :: ispin, ip, idir, iatom, ilaser, itime
1331 real(real64) :: aa(space%dim), time_
1332 real(real64),
allocatable :: vp(:,:)
1333 real(real64),
allocatable :: velectric(:)
1334 type(lasers_t),
pointer :: lasers
1335 type(gauge_field_t),
pointer :: gfield
1338 call profiling_in(
"HAMILTONIAN_ELEC_UPDATE_EXT_POT")
1340 this%current_time = m_zero
1341 this%current_time = time(1)
1344 call this%hm_base%clear(mesh%np)
1347 call this%hm_base%allocate_field(mesh, field_potential, &
1348 complex_potential = this%abs_boundaries%abtype == imaginary_absorbing)
1353 lasers => list_get_lasers(ext_partners)
1354 if(
associated(lasers))
then
1355 do ilaser = 1, lasers%no_lasers
1356 select case (laser_kind(lasers%lasers(ilaser)))
1357 case (e_field_scalar_potential, e_field_electric)
1358 safe_allocate(velectric(1:mesh%np))
1359 do ispin = 1, this%d%spin_channels
1361 call laser_potential(lasers%lasers(ilaser), mesh, velectric, time_)
1364 this%hm_base%potential(ip, ispin) = this%hm_base%potential(ip, ispin) + mu(itime) * velectric(ip)
1367 safe_deallocate_a(velectric)
1368 case (e_field_magnetic)
1369 call this%hm_base%allocate_field(mesh, field_vector_potential + field_uniform_magnetic_field, .false.)
1371 safe_allocate(vp(1:mesh%np, 1:space%dim))
1372 vp(1:mesh%np, 1:space%dim) = m_zero
1373 call laser_vector_potential(lasers%lasers(ilaser), mesh, vp, time_)
1374 do idir = 1, space%dim
1377 this%hm_base%vector_potential(idir, ip) = this%hm_base%vector_potential(idir, ip) &
1378 - mu(itime) * vp(ip, idir)/p_c
1382 call laser_field(lasers%lasers(ilaser), this%hm_base%uniform_magnetic_field(1:space%dim), time_)
1383 safe_deallocate_a(vp)
1384 case (e_field_vector_potential)
1385 call this%hm_base%allocate_field(mesh, field_uniform_vector_potential, .false.)
1388 call laser_field(lasers%lasers(ilaser), aa, time_)
1389 this%hm_base%uniform_vector_potential(1:space%dim) = this%hm_base%uniform_vector_potential(1:space%dim) &
1390 - mu(itime) * aa/p_c
1396 gfield => list_get_gauge_field(ext_partners)
1397 if (
associated(gfield))
then
1398 call this%hm_base%allocate_field(mesh, field_uniform_vector_potential, .false.)
1399 call gauge_field_get_vec_pot(gfield, aa)
1400 this%hm_base%uniform_vector_potential(1:space%dim) = this%hm_base%uniform_vector_potential(1:space%dim) - aa/p_c
1407 if (
allocated(this%ep%e_field) .and.
associated(gfield))
then
1408 this%hm_base%uniform_vector_potential(1:space%periodic_dim) = &
1409 this%hm_base%uniform_vector_potential(1:space%periodic_dim) - time_*this%ep%e_field(1:space%periodic_dim)
1415 if (
allocated(this%ep%a_static))
then
1416 call this%hm_base%allocate_field(mesh, field_vector_potential, .false.)
1420 do idir = 1, space%dim
1421 this%hm_base%vector_potential(idir, ip) = this%hm_base%vector_potential(idir, ip) + this%ep%a_static(ip, idir)
1427 call this%hm_base%accel_copy_pot(mesh)
1430 if (
allocated(this%ep%b_field))
then
1431 call this%hm_base%allocate_field(mesh, field_uniform_magnetic_field, .false.)
1433 this%hm_base%uniform_magnetic_field(idir) = this%hm_base%uniform_magnetic_field(idir) + this%ep%b_field(idir)
1437 call this%hm_base%update_magnetic_terms(mesh, this%ep%gyromagnetic_ratio, this%d%ispin)
1443 call profiling_out(
"HAMILTONIAN_ELEC_UPDATE_EXT_POT")
1449 integer :: ik, imat, nmat, max_npoints, offset, iphase, nphase
1453 if ((.not. this%kpoints%gamma_only()) .or.
allocated(this%hm_base%uniform_vector_potential))
then
1455 call profiling_in(
'UPDATE_PHASES')
1457 do iatom = 1, this%ep%natoms
1458 call projector_init_phases(this%ep%proj(iatom), space%dim, this%d, this%der%boundaries, this%kpoints, &
1459 vec_pot = this%hm_base%uniform_vector_potential, vec_pot_var = this%hm_base%vector_potential)
1462 call profiling_out(
'UPDATE_PHASES')
1465 if (
allocated(this%hm_base%uniform_vector_potential))
then
1466 call this%phase%update(mesh, this%d%kpt, this%kpoints, this%d, space, this%hm_base%uniform_vector_potential)
1469 max_npoints = this%vnl%max_npoints
1470 nmat = this%vnl%nprojector_matrices
1473 if (this%phase%is_allocated() .and.
allocated(this%vnl%projector_matrices))
then
1476 if (this%der%boundaries%spiralBC) nphase = 3
1478 if (.not.
allocated(this%vnl%projector_phases))
then
1479 safe_allocate(this%vnl%projector_phases(1:max_npoints, nphase, nmat, this%d%kpt%start:this%d%kpt%end))
1480 if (accel_is_enabled())
then
1481 call accel_create_buffer(this%vnl%buff_projector_phases, accel_mem_read_only, &
1482 type_cmplx, this%vnl%total_points*nphase*this%d%kpt%nlocal)
1487 do ik = this%d%kpt%start, this%d%kpt%end
1488 do imat = 1, this%vnl%nprojector_matrices
1489 iatom = this%vnl%projector_to_atom(imat)
1490 do iphase = 1, nphase
1492 do ip = 1, this%vnl%projector_matrices(imat)%npoints
1493 this%vnl%projector_phases(ip, imat, iphase, ik) = this%ep%proj(iatom)%phase(ip, iphase, ik)
1496 if (accel_is_enabled() .and. this%vnl%projector_matrices(imat)%npoints > 0)
then
1497 call accel_write_buffer(this%vnl%buff_projector_phases, &
1498 this%vnl%projector_matrices(imat)%npoints, this%vnl%projector_phases(1:, iphase, imat, ik), &
1501 offset = offset + this%vnl%projector_matrices(imat)%npoints
1515 logical,
intent(in) :: states_are_real
1519 if (hm%self_induced_magnetic)
then
1520 if (.not. states_are_real)
then
1523 message(1) =
'No current density for real states since it is identically zero.'
1524 call messages_warning(1)
1533 type(namespace_t),
intent(in) :: namespace
1534 type(grid_t),
intent(in) :: gr
1535 type(states_elec_t),
intent(inout) :: st
1536 type(states_elec_t),
intent(inout) :: hst
1538 integer :: ik, ib, ist
1539 complex(real64),
allocatable :: psi(:, :)
1540 complex(real64),
allocatable :: psiall(:, :, :, :)
1544 do ik = st%d%kpt%start, st%d%kpt%end
1545 do ib = st%group%block_start, st%group%block_end
1550 if (oct_exchange_enabled(hm%oct_exchange))
then
1552 safe_allocate(psiall(gr%np_part, 1:hst%d%dim, st%st_start:st%st_end, st%d%kpt%start:st%d%kpt%end))
1554 call states_elec_get_state(st, gr, psiall)
1556 call oct_exchange_prepare(hm%oct_exchange, gr, psiall, hm%xc, hm%psolver, namespace)
1558 safe_deallocate_a(psiall)
1560 safe_allocate(psi(gr%np_part, 1:hst%d%dim))
1564 call states_elec_get_state(hst, gr, ist, ik, psi)
1565 call oct_exchange_operator(hm%oct_exchange, namespace, gr, psi, ist, ik)
1566 call states_elec_set_state(hst, gr, ist, ik, psi)
1570 safe_deallocate_a(psi)
1580 subroutine magnus(hm, namespace, mesh, psi, hpsi, ik, vmagnus, set_phase)
1582 type(namespace_t),
intent(in) :: namespace
1583 class(mesh_t),
intent(in) :: mesh
1584 complex(real64),
contiguous,
intent(inout) :: psi(:,:)
1585 complex(real64),
contiguous,
intent(out) :: hpsi(:,:)
1586 integer,
intent(in) :: ik
1587 real(real64),
intent(in) :: vmagnus(:, :, :)
1588 logical,
optional,
intent(in) :: set_phase
1590 complex(real64),
allocatable :: auxpsi(:, :), aux2psi(:, :)
1591 integer :: idim, ispin
1596 if (hm%d%dim /= 1)
then
1597 call messages_not_implemented(
"Magnus with spinors", namespace=namespace)
1600 safe_allocate( auxpsi(1:mesh%np_part, 1:hm%d%dim))
1601 safe_allocate(aux2psi(1:mesh%np, 1:hm%d%dim))
1603 ispin = hm%d%get_spin_index(ik)
1607 terms = term_kinetic + term_non_local_potential, set_phase = set_phase)
1610 do idim = 1, hm%d%dim
1611 call lalg_copy(mesh%np, auxpsi(:, idim), hpsi(:, idim))
1612 hpsi(1:mesh%np, idim) = hpsi(1:mesh%np, idim) + hm%ep%Vpsl(1:mesh%np)*psi(1:mesh%np,idim)
1613 call vborders(mesh, hm, psi(:, idim), hpsi(:, idim))
1615 hpsi(1:mesh%np, 1) = hpsi(1:mesh%np, 1) + vmagnus(1:mesh%np, ispin, 2)*psi(1:mesh%np, 1)
1618 hpsi(1:mesh%np, 1) = hpsi(1:mesh%np, 1) - m_zi*vmagnus(1:mesh%np, ispin, 1)*auxpsi(1:mesh%np, 1)
1621 auxpsi(1:mesh%np, 1) = vmagnus(1:mesh%np, ispin, 1)*psi(1:mesh%np, 1)
1623 terms = term_kinetic + term_non_local_potential, set_phase = set_phase)
1624 hpsi(1:mesh%np, 1) = hpsi(1:mesh%np, 1) + m_zi*aux2psi(1:mesh%np, 1)
1626 safe_deallocate_a(auxpsi)
1627 safe_deallocate_a(aux2psi)
1632 subroutine vborders (mesh, hm, psi, hpsi)
1633 class(mesh_t),
intent(in) :: mesh
1635 complex(real64),
intent(in) :: psi(:)
1636 complex(real64),
intent(inout) :: hpsi(:)
1642 if (hm%abs_boundaries%abtype == imaginary_absorbing)
then
1644 hpsi(ip) = hpsi(ip) + m_zi*hm%abs_boundaries%mf(ip)*psi(ip)
1666 type(namespace_t),
intent(in) :: namespace
1667 real(real64),
intent(in) :: mass
1671 if (parse_is_defined(namespace,
'ParticleMass'))
then
1672 message(1) =
'Attempting to redefine a non-unit electron mass'
1673 call messages_fatal(1)
1690 type(grid_t),
intent(in) :: gr
1691 type(distributed_t),
intent(in) :: kpt
1692 type(wfs_elec_t),
intent(in) :: psib
1693 type(wfs_elec_t),
intent(out) :: psib_with_phase
1695 integer :: k_offset, n_boundary_points
1699 call psib%copy_to(psib_with_phase)
1700 if (hm%phase%is_allocated())
then
1701 call hm%phase%apply_to(gr, gr%np, conjugate = .false., psib = psib_with_phase, src = psib, async=.
true.)
1704 k_offset = psib%ik - kpt%start
1705 n_boundary_points = int(gr%np_part - gr%np)
1706 call boundaries_set(gr%der%boundaries, gr, psib_with_phase, phase_correction = hm%phase%phase_corr(:, psib%ik), &
1707 buff_phase_corr = hm%phase%buff_phase_corr, offset=k_offset*n_boundary_points, async=.
true.)
1709 call psib%copy_data_to(gr%np, psib_with_phase)
1710 call boundaries_set(gr%der%boundaries, gr, psib_with_phase)
1713 call psib_with_phase%do_pack(copy = .
true.)
1721 class(mesh_t),
intent(in) :: mesh
1722 real(real64),
contiguous,
intent(out) :: diag(:,:)
1723 integer,
intent(in) :: ik
1725 integer :: idim, ip, ispin
1727 real(real64),
allocatable :: ldiag(:)
1731 safe_allocate(ldiag(1:mesh%np))
1735 call derivatives_lapl_diag(hm%der, ldiag)
1737 do idim = 1, hm%d%dim
1738 diag(1:mesh%np, idim) = -m_half/hm%mass*ldiag(1:mesh%np)
1741 select case (hm%d%ispin)
1743 case (unpolarized, spin_polarized)
1744 ispin = hm%d%get_spin_index(ik)
1745 diag(1:mesh%np, 1) = diag(1:mesh%np, 1) + hm%ep%vpsl(1:mesh%np)
1749 diag(ip, 1) = diag(ip, 1) + hm%ep%vpsl(ip)
1750 diag(ip, 2) = diag(ip, 2) + hm%ep%vpsl(ip)
1755 call hm%ks_pot%add_vhxc(diag)
1764#include "hamiltonian_elec_inc.F90"
1767#include "complex.F90"
1768#include "hamiltonian_elec_inc.F90"
subroutine build_external_potentials()
subroutine external_potentials_checks()
subroutine build_interactions()
constant times a vector plus a vector
Copies a vector x, to a vector y.
integer, parameter, public imaginary_absorbing
subroutine, public absorbing_boundaries_end(this)
subroutine, public absorbing_boundaries_init(this, namespace, space, gr)
pure logical function, public accel_allow_cpu_only()
subroutine, public accel_finish()
pure logical function, public accel_is_enabled()
integer, parameter, public accel_mem_read_only
This module implements batches of mesh functions.
This module implements common operations on batches of mesh functions.
Module implementing boundary conditions in Octopus.
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
real(real64) function, public derivatives_lapl_get_max_eigenvalue(this)
Get maximum eigenvalue of discrete Laplacian. For the star and star_general stencils,...
logical function, public epot_have_external_potentials(ep)
integer, parameter, public scalar_relativistic_zora
subroutine, public epot_end(ep)
integer, parameter, public fully_relativistic_zora
subroutine, public epot_init(ep, namespace, gr, ions, psolver, ispin, xc_family, kpoints)
subroutine, public epot_generate(ep, namespace, mesh, ions, st_d)
subroutine, public exchange_operator_init(this, namespace, space, st, der, mc, stencil, kpoints, cam)
subroutine, public exchange_operator_end(this)
logical function, public list_has_gauge_field(partners)
type(gauge_field_t) function, pointer, public list_get_gauge_field(partners)
logical function, public list_has_lasers(partners)
type(lasers_t) function, pointer, public list_get_lasers(partners)
integer, parameter, public external_pot_from_file
potential, defined in a file
subroutine, public load_external_potentials(external_potentials, namespace)
integer, parameter, public external_pot_charge_density
user-defined function for charge density
integer, parameter, public external_pot_usdef
user-defined function for local potential
integer, parameter, public external_pot_static_efield
Static electric field.
integer, parameter, public external_pot_static_bfield
Static magnetic field.
subroutine, public gauge_field_get_vec_pot(this, vec_pot)
real(real64), parameter, public m_two
real(real64), parameter, public m_zero
integer, parameter, public rdmft
integer, parameter, public hartree_fock
integer, parameter, public independent_particles
Theory level.
integer, parameter, public generalized_kohn_sham_dft
integer, parameter, public kohn_sham_dft
real(real64), parameter, public m_epsilon
real(real64), parameter, public p_c
Electron gyromagnetic ratio, see Phys. Rev. Lett. 130, 071801 (2023)
real(real64), parameter, public m_one
This module implements the underlying real-space grid.
This module defines an abstract class for Hamiltonians.
integer, parameter, public field_uniform_magnetic_field
integer, parameter, public field_uniform_vector_potential
integer, parameter, public field_vector_potential
integer, parameter, public term_mgga
integer, parameter, public field_potential
subroutine vborders(mesh, hm, psi, hpsi)
pure logical function, public hamiltonian_elec_apply_packed(this)
subroutine, public hamiltonian_elec_set_inh(hm, st)
subroutine, public zvmask(mesh, hm, st)
subroutine, public zhamiltonian_elec_apply_batch(hm, namespace, mesh, psib, hpsib, terms, set_bc)
subroutine, public hamiltonian_elec_adjoint(hm)
subroutine, public hamiltonian_elec_end(hm)
subroutine, public hamiltonian_elec_copy_and_set_phase(hm, gr, kpt, psib, psib_with_phase)
Copy a batch to another batch and apply the Bloch phase to it.
subroutine, public zhamiltonian_elec_apply_single(hm, namespace, mesh, psi, hpsi, ist, ik, terms, set_bc, set_phase)
pure logical function hamiltonian_elec_needs_mgga_term(hm, terms)
logical function, public hamiltonian_elec_has_kick(hm)
subroutine, public magnus(hm, namespace, mesh, psi, hpsi, ik, vmagnus, set_phase)
logical function hamiltonian_elec_hermitian(hm)
subroutine, public dhamiltonian_elec_apply_single(hm, namespace, mesh, psi, hpsi, ist, ik, terms, set_bc, set_phase)
subroutine, public hamiltonian_elec_epot_generate(this, namespace, space, gr, ions, ext_partners, st, time)
real(real64) function, public hamiltonian_elec_get_time(this)
subroutine, public dvmask(mesh, hm, st)
logical function, public hamiltonian_elec_needs_current(hm, states_are_real)
subroutine, public hamiltonian_elec_remove_inh(hm)
subroutine zhamiltonian_elec_magnus_apply(hm, namespace, mesh, psib, hpsib, vmagnus)
subroutine, public zhamiltonian_elec_apply_atom(hm, namespace, space, latt, species, pos, ia, mesh, psi, vpsi)
subroutine, public zhamiltonian_elec_apply_all(hm, namespace, gr, st, hst)
subroutine, public hamiltonian_elec_diagonal(hm, mesh, diag, ik)
subroutine, public hamiltonian_elec_update_pot(this, mesh, accumulate)
Update the KS potential of the electronic Hamiltonian.
integer, parameter, public velocity
subroutine dhamiltonian_elec_magnus_apply(hm, namespace, mesh, psib, hpsib, vmagnus)
subroutine hamiltonian_elec_update(this, mesh, namespace, space, ext_partners, time)
(re-)build the Hamiltonian for the next application:
subroutine hamiltonian_elec_span(hm, delta, emin, namespace)
subroutine, public hamiltonian_elec_init(hm, namespace, space, gr, ions, ext_partners, st, theory_level, xc, mc, kpoints, need_exchange, xc_photons)
subroutine dhamiltonian_elec_apply(hm, namespace, mesh, psib, hpsib, terms, set_bc)
pure logical function, public hamiltonian_elec_inh_term(hm)
subroutine hamiltonian_elec_set_mass(this, namespace, mass)
set the effective electron mass, checking whether it was previously redefined.
subroutine, public dhamiltonian_elec_apply_batch(hm, namespace, mesh, psib, hpsib, terms, set_bc)
subroutine, public hamiltonian_elec_update_with_ext_pot(this, mesh, space, ext_partners, time, mu)
This is an extension of "hamiltonian_elec_update_pot" to be used by the CFM4 propagator....
subroutine, public hamiltonian_elec_not_adjoint(hm)
subroutine zhamiltonian_elec_apply(hm, namespace, mesh, psib, hpsib, terms, set_bc)
This module defines classes and functions for interaction partners.
integer, parameter, public kick_magnon_mode
subroutine, public kick_end(kick)
subroutine, public kick_init(kick, namespace, space, kpoints, nspin)
pure integer function, public kick_get_type(kick)
A module to handle KS potential, without the external potential.
subroutine, public laser_vector_potential(laser, mesh, aa, time)
subroutine, public lasers_nondipole_laser_field_step(this, field, time)
Retrieves the NDSFA vector_potential correction. The nondipole field is obtained for consecutive time...
logical function, public lasers_with_nondipole_field(lasers)
Check if a nondipole SFA correction should be computed for the given laser.
integer, parameter, public e_field_electric
integer, parameter, public e_field_vector_potential
subroutine, public laser_potential(laser, mesh, pot, time)
integer, parameter, public e_field_scalar_potential
integer pure elemental function, public laser_kind(laser)
subroutine, public laser_field(laser, field, time)
Retrieves the value of either the electric or the magnetic field. If the laser is given by a scalar p...
integer, parameter, public e_field_magnetic
integer, parameter, public dft_u_none
subroutine, public lda_u_update_basis(this, space, gr, ions, st, psolver, namespace, kpoints, has_phase)
subroutine, public lda_u_build_phase_correction(this, space, std, boundaries, namespace, kpoints, vec_pot, vec_pot_var)
Build the phase correction to the global phase for all orbitals.
subroutine, public lda_u_end(this)
subroutine, public lda_u_init(this, namespace, space, level, gr, ions, st, mc, kpoints, has_phase)
This module implements fully polymorphic linked lists, and some specializations thereof.
This modules implements the routines for doing constrain DFT for noncollinear magnetism.
subroutine, public magnetic_constrain_end(this)
Releases memory of the magnetic constrain.
subroutine, public magnetic_constrain_init(this, namespace, mesh, std, natoms, min_dist)
Initilializes the magnetic_constrain_t object.
This module is intended to contain "only mathematical" functions and procedures.
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
subroutine, public messages_not_implemented(feature, namespace)
subroutine, public messages_warning(no_lines, all_nodes, namespace)
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public messages_experimental(name, namespace)
This module handles the communicators for the various parallelization strategies.
subroutine, public mxll_coupling_init(this, d, gr, namespace, mass)
Parse variables and initialize Maxwell coupling.
subroutine, public set_electric_quadrupole_pot(this, mesh)
Computes the electric quadrupole potential where .
subroutine, public mxll_coupling_end(this)
Finalize and deallocate Maxwell coupling arrays.
subroutine, public mxll_coupling_calc(this, hm_base, mesh, d, space)
Add the Maxwell coupling to the electronic Hamiltonian.
logical function, public oct_exchange_enabled(this)
subroutine, public oct_exchange_remove(this)
this module contains the low-level part of the output system
Some general things and nomenclature:
logical function, public parse_is_defined(namespace, name)
subroutine, public pcm_calc_pot_rs(pcm, mesh, psolver, ions, v_h, v_ext, kick, time_present, kick_time)
subroutine, public pcm_end(pcm)
real(real64), dimension(:,:), allocatable delta
D_E matrix in JCP 139, 024105 (2013).
subroutine, public pcm_init(pcm, namespace, space, ions, grid, qtot, val_charge, external_potentials_present, kick_present)
Initializes the PCM calculation: reads the VdW molecular cavity and generates the PCM response matrix...
subroutine, public poisson_init(this, namespace, space, der, mc, stencil, qtot, label, solver, verbose, force_serial, force_cmplx)
subroutine, public poisson_end(this)
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
subroutine, public projector_init_phases(this, dim, std, bnd, kpoints, vec_pot, vec_pot_var)
subroutine, public scissor_end(this)
subroutine, public states_set_complex(st)
pure logical function, public states_are_real(st)
This module handles spin dimensions of the states and the k-point distribution.
subroutine, public states_elec_dim_copy(dout, din)
subroutine, public states_elec_dim_end(dim)
subroutine, public states_elec_end(st)
finalize the states_elec_t object
subroutine, public states_elec_copy(stout, stin, exclude_wfns, exclude_eigenval, special)
make a (selective) copy of a states_elec_t object
This module provides routines for communicating states when using states parallelization.
integer pure function, public symmetries_identity_index(this)
integer pure function, public symmetries_number(this)
type(type_t), public type_cmplx
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
This module defines the unit system, used for input and output.
type(unit_system_t), public units_inp
the units systems for reading and writing
type(xc_cam_t), parameter, public cam_exact_exchange
Use only Hartree Fock exact exchange.
integer, parameter, public xc_oep_x_slater
Slater approximation to the exact exchange.
integer, parameter, public func_x
logical pure function, public family_is_mgga_with_exc(xcs)
Is the xc function part of the mGGA family with an energy functional.
logical pure function, public family_is_hybrid(xcs)
Returns true if the functional is an hybrid functional.
This module implements the "photon-free" electron-photon exchange-correlation functional.
This module implements the ZORA terms for the Hamoiltonian.
Extension of space that contains the knowledge of the spin dimension.
Description of the grid, containing information on derivatives, stencil, and symmetries.
The abstract Hamiltonian class defines a skeleton for specific implementations.
abstract class for general interaction partners
iterator for the list of partners
Describes mesh distribution to nodes.
Stores all communicators and groups.
The states_elec_t class contains all electronic wave functions.
This class described the 'photon-exchange' electron-photon xc functionals, based on QEDFT.
This class is responsible for calculating and applying the ZORA.