28 use,
intrinsic :: iso_fortran_env
57 integer,
public,
parameter :: &
58 SM_POISSON_DIRECT = 0, &
72 class(species_t),
intent(in) :: species
73 integer,
optional,
intent(in) :: iselect
75 integer :: iorb, ii, ll, mm
78 do iorb = 1, species%get_niwfs()
79 call species%get_iwf_ilm(iorb, 1, ii, ll, mm)
80 if (
present(iselect))
then
81 if (ii == iselect) norb = norb + 1
88 subroutine orbitalset_init_intersite(this, namespace, space, ind, ions, der, psolver, os, nos, maxnorbs, &
89 rcut, kpt, has_phase, sm_poisson, basis_from_states, combine_j_orbitals)
90 type(orbitalset_t),
intent(inout) :: this
91 type(namespace_t),
intent(in) :: namespace
92 class(space_t),
intent(in) :: space
93 integer,
intent(in) :: ind
94 type(ions_t),
intent(in) :: ions
95 type(derivatives_t),
intent(in) :: der
96 type(poisson_t),
intent(in) :: psolver
97 type(orbitalset_t),
intent(inout) :: os(:)
98 integer,
intent(in) :: nos, maxnorbs
99 real(real64),
intent(in) :: rcut
100 type(distributed_t),
intent(in) :: kpt
101 logical,
intent(in) :: has_phase
102 integer,
intent(in) :: sm_poisson
103 logical,
intent(in) :: basis_from_states
104 logical,
intent(in) :: combine_j_orbitals
107 type(lattice_iterator_t) :: latt_iter
108 real(real64) :: xat(space%dim), xi(space%dim)
110 integer :: inn, ist, jst
111 integer :: ip, ios, ip2, is1, is2
112 type(submesh_t) :: sm
113 real(real64),
allocatable :: tmp(:), vv(:), nn(:)
114 complex(real64),
allocatable :: ztmp(:), zvv(:,:), znn(:)
115 real(real64),
allocatable :: dorb(:,:,:,:)
116 complex(real64),
allocatable :: zorb(:,:,:,:)
117 real(real64),
parameter :: TOL_INTERSITE = 1.e-5_real64
118 type(distributed_t) :: dist
119 integer :: sm_poisson_
125 sm_poisson_ = sm_poisson
128 if (this%iatom /= -1)
then
129 xat = ions%pos(:, this%iatom)
131 xat = this%sphere%center
138 do inn = 1, latt_iter%n_cells
140 xi = os(ios)%sphere%center(1:space%dim) + latt_iter%get(inn)
144 if( rr >rcut + tol_intersite ) cycle
146 if( ios == ind .and. rr < tol_intersite) cycle
148 this%nneighbors = this%nneighbors +1
154 safe_allocate(this%V_ij(1:this%nneighbors, 0:space%dim+1))
155 this%V_ij(1:this%nneighbors, 0:space%dim+1) =
m_zero
156 safe_allocate(this%map_os(1:this%nneighbors))
157 this%map_os(1:this%nneighbors) = 0
159 safe_allocate(this%phase_shift(1:this%nneighbors, kpt%start:kpt%end))
164 do inn = 1, latt_iter%n_cells
165 xi = os(ios)%sphere%center(1:space%dim) + latt_iter%get(inn)
168 if( rr > rcut + tol_intersite ) cycle
169 if( ios == ind .and. rr < tol_intersite) cycle
171 this%nneighbors = this%nneighbors +1
173 this%V_ij(this%nneighbors, 1:space%dim) = xi(1:space%dim) -os(ios)%sphere%center(1:space%dim)
174 this%V_ij(this%nneighbors, space%dim+1) = rr
176 this%map_os(this%nneighbors) = ios
180 write(
message(1),
'(a, i3, a)')
'Intersite interaction will be computed for ', this%nneighbors,
' neighboring atoms.'
184 if (this%ndim == 1)
then
185 safe_allocate(this%coulomb_IIJJ(1:this%norbs,1:this%norbs,1:maxnorbs,1:maxnorbs,1:this%nneighbors))
186 this%coulomb_IIJJ =
m_zero
188 safe_allocate(this%zcoulomb_IIJJ(1:this%norbs,1:this%norbs,1:maxnorbs,1:maxnorbs, 1:this%ndim, 1:this%ndim, 1:this%nneighbors))
189 this%zcoulomb_IIJJ =
m_zero
192 if(this%nneighbors == 0)
then
200 if(.not. der%mesh%parallel_in_domains)
then
205 do inn = dist%start, dist%end
207 ios = this%map_os(inn)
209 if(.not. basis_from_states)
then
211 call submesh_merge(sm, space, der%mesh, this%sphere, os(ios)%sphere, &
212 shift = this%V_ij(inn, 1:space%dim))
214 write(
message(1),
'(a, i3, a, f6.3, a, i7, a)')
'Neighbor ', inn,
' is located at ', &
215 this%V_ij(inn, space%dim+1),
' Bohr and has ', sm%np,
' grid points.'
218 if (this%ndim == 1)
then
219 safe_allocate(dorb(1:sm%np, 1:1, 1:max(this%norbs, os(ios)%norbs), 1:2))
221 safe_allocate(zorb(1:sm%np, 1:2, 1:max(this%norbs, os(ios)%norbs), 1:2))
225 if (this%ndim == 1)
then
226 do ist = 1, this%norbs
227 call dget_orbital(this%spec, sm, this%ii, this%ll, this%jj, ist, 1, dorb(:, :, ist, 1), combine_j_orbitals)
230 do ist = 1, this%norbs
231 call zget_orbital(this%spec, sm, this%ii, this%ll, this%jj, ist, 2, zorb(:, :, ist, 1), combine_j_orbitals)
235 call submesh_shift_center(sm, space, this%V_ij(inn, 1:space%dim)+os(ios)%sphere%center(1:space%dim))
238 if (this%ndim == 1)
then
239 do ist = 1, os(ios)%norbs
240 call dget_orbital(os(ios)%spec, sm, os(ios)%ii, os(ios)%ll, os(ios)%jj, ist, 1, dorb(:, :, ist, 2), combine_j_orbitals)
243 do ist = 1, os(ios)%norbs
244 call zget_orbital(os(ios)%spec, sm, os(ios)%ii, os(ios)%ll, os(ios)%jj, ist, 2, zorb(:, :, ist, 2), combine_j_orbitals)
250 assert(this%ndim == 1)
252 call submesh_merge(sm, ions%space, der%mesh, this%sphere, os(ios)%sphere, &
253 shift = this%V_ij(inn, 1:ions%space%dim))
255 write(
message(1),
'(a, i3, a, f6.3, a, i5, a)')
'Neighbor ', inn,
' is located at ', &
256 this%V_ij(inn, ions%space%dim+1),
' Bohr and has ', sm%np,
' grid points.'
259 safe_allocate(dorb(1:sm%np, 1:1, 1:max(this%norbs,os(ios)%norbs), 1:2))
263 do ist = 1, this%norbs
264 if(
allocated(this%dorb))
then
265 dorb(1:this%sphere%np, 1, ist, 1) = this%dorb(:,1,ist)
267 dorb(1:this%sphere%np, 1, ist, 1) = real(this%zorb(:,1,ist), real64)
271 call submesh_shift_center(sm, ions%space, this%V_ij(inn, 1:ions%space%dim)+os(ios)%sphere%center(1:ions%space%dim))
275 do ist = 1, os(ios)%norbs
276 if(
allocated(this%dorb))
then
279 do ip = 1, os(ios)%sphere%np
280 if(all(abs(sm%rel_x(1:ions%space%dim, ip2)-os(ios)%sphere%rel_x(1:ions%space%dim, ip)) < 1e-6_real64))
then
281 dorb(ip2, 1, ist, 2) = os(ios)%dorb(ip, 1, ist)
288 do ip = 1, os(ios)%sphere%np
289 if(all(abs(sm%rel_x(1:ions%space%dim, ip2)-os(ios)%sphere%rel_x(1:ions%space%dim, ip)) < 1e-6_real64))
then
290 dorb(ip2, 1, ist, 2) = real(os(ios)%zorb(ip, 1, ist), real64)
299 select case (sm_poisson_)
300 case(sm_poisson_direct)
310 force_cmplx=(this%ndim==2))
313 if (this%ndim == 1)
then
314 safe_allocate(tmp(1:sm%np))
315 safe_allocate(nn(1:sm%np))
316 safe_allocate(vv(1:sm%np))
318 do ist = 1, this%norbs
321 nn(ip) = dorb(ip, 1, ist, 1)*dorb(ip, 1, ist, 1)
329 do jst = 1, os(ios)%norbs
333 tmp(ip) = vv(ip)*dorb(ip, 1, jst, 2)*dorb(ip, 1, jst, 2)
337 this%coulomb_IIJJ(ist, ist, jst, jst, inn) =
dsm_integrate(der%mesh, sm, tmp, reduce = .false.)
341 safe_deallocate_a(nn)
342 safe_deallocate_a(vv)
343 safe_deallocate_a(tmp)
346 safe_allocate(ztmp(1:sm%np))
347 safe_allocate(znn(1:sm%np))
348 safe_allocate(zvv(1:sm%np, 1:this%ndim))
350 do ist = 1, this%norbs
351 do is1 = 1, this%ndim
354 znn(ip) = conjg(zorb(ip, is1, ist, 1))*zorb(ip, is1, ist, 1)
363 do jst = 1, os(ios)%norbs
365 do is1 = 1, this%ndim
366 do is2 = 1, this%ndim
369 ztmp(ip) = zvv(ip, is1)*conjg(zorb(ip, is2, jst, 2))*zorb(ip, is2, jst, 2)
373 this%zcoulomb_IIJJ(ist, ist, jst, jst, is1, is2, inn) =
zsm_integrate(der%mesh, sm, ztmp)
379 safe_deallocate_a(znn)
380 safe_deallocate_a(zvv)
381 safe_deallocate_a(ztmp)
390 safe_deallocate_a(dorb)
391 safe_deallocate_a(zorb)
394 if(der%mesh%parallel_in_domains .and. this%ndim == 1)
then
395 call der%mesh%allreduce(this%coulomb_IIJJ)
398 if(dist%parallel)
then
399 if (this%ndim == 1)
then
402 do inn = 1, this%nneighbors
403 do is2 = 1, this%ndim
404 do is1 = 1, this%ndim
405 call comm_allreduce(dist%mpi_grp, this%zcoulomb_IIJJ(:,:,:,:, is1, is2, inn))
421#include "orbitalset_utils_inc.F90"
424#include "complex.F90"
425#include "orbitalset_utils_inc.F90"
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
subroutine, public distributed_end(this)
subroutine, public distributed_nullify(this, total)
subroutine, public distributed_init(this, total, comm, tag, scalapack_compat)
Distribute N instances across M processes of communicator comm
real(real64), parameter, public m_zero
This module is intended to contain "only mathematical" functions and procedures.
This module defines the meshes, which are used in Octopus.
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
character(len=512), private msg
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
integer, parameter, public sm_poisson_psolver
subroutine zget_orbital(species, sm, i, l, j_in, iorb, ndim, orb, combine_j_orbitals)
Returns an orbital on the intersite submesh.
subroutine, public dorbitalset_utils_getorbitals(namespace, os, ions, mesh, use_mesh, normalize, index_shift)
integer, parameter, public sm_poisson_isf
subroutine, public zorbitalset_get_center_of_mass(os, space, mesh, latt)
subroutine, public zorbitalset_utils_getorbitals(namespace, os, ions, mesh, use_mesh, normalize, index_shift)
subroutine, public orbitalset_init_intersite(this, namespace, space, ind, ions, der, psolver, os, nos, maxnorbs, rcut, kpt, has_phase, sm_poisson, basis_from_states, combine_j_orbitals)
subroutine dget_orbital(species, sm, i, l, j_in, iorb, ndim, orb, combine_j_orbitals)
Returns an orbital on the intersite submesh.
subroutine, public dorbitalset_get_center_of_mass(os, space, mesh, latt)
integer, parameter, public sm_poisson_fft
integer function, public orbitalset_utils_count(species, iselect)
Count the number of orbital sets we have for a given atom.
subroutine, public zpoisson_solve_sm(this, namespace, sm, pot, rho, all_nodes)
Calculates the Poisson equation. Given the density returns the corresponding potential.
integer, parameter, public poisson_psolver
integer, parameter, public poisson_fft
subroutine, public dpoisson_solve_sm(this, namespace, sm, pot, rho, all_nodes)
Calculates the Poisson equation. Given the density returns the corresponding potential.
integer, parameter, public poisson_direct_sum
subroutine, public poisson_init_sm(this, namespace, space, main, der, sm, method, force_cmplx)
integer, parameter, public poisson_isf
subroutine, public poisson_end(this)
subroutine, public submesh_end_global(this)
complex(real64) function, public zsm_integrate(mesh, sm, ff, reduce)
subroutine, public submesh_shift_center(this, space, newcenter)
real(real64) function, public dsm_integrate(mesh, sm, ff, reduce)
subroutine, public submesh_merge(this, space, mesh, sm1, sm2, shift)
subroutine, public submesh_end_cube_map(sm)
subroutine, public submesh_end(this)
subroutine, public submesh_build_global(this, space)
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
The following class implements a lattice iterator. It allows one to loop over all cells that are with...