Octopus
nonlocal_pseudopotential.F90
Go to the documentation of this file.
1!! Copyright (C) 2009 X. Andrade
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19#include "global.h"
20
22 use accel_oct_m
24 use batch_oct_m
26 use blas_oct_m
27 use comm_oct_m
28 use debug_oct_m
30 use epot_oct_m
31 use global_oct_m
36 use math_oct_m
37 use mesh_oct_m
39 use mpi_oct_m
43 use ps_oct_m
45 use space_oct_m
49 use types_oct_m
51
52 implicit none
53
54 private
55
56 public :: &
59
63 private
64 type(projector_matrix_t), allocatable, public :: projector_matrices(:)
65 integer, public :: nprojector_matrices
66 logical, public :: apply_projector_matrices
67 logical, public :: has_non_local_potential
68 integer :: full_projection_size
69 integer, public :: max_npoints
70 integer, public :: total_points
71 integer :: max_nprojs
72 logical :: projector_mix
73 complex(real64), allocatable, public :: projector_phases(:, :, :, :)
74 integer, allocatable, public :: projector_to_atom(:)
75 integer :: nregions
76 integer, allocatable :: regions(:)
77 integer, public :: nphase
82 type(accel_mem_t) :: buff_offsets
83 type(accel_mem_t) :: buff_matrices
84 type(accel_mem_t) :: buff_maps
85 type(accel_mem_t) :: buff_scals
86 type(accel_mem_t) :: buff_position
87 type(accel_mem_t) :: buff_pos
88 type(accel_mem_t) :: buff_invmap
89 type(accel_mem_t), public :: buff_projector_phases
90 type(accel_mem_t) :: buff_mix
91 logical :: projector_self_overlap
92 real(real64), pointer, public :: spin(:,:,:) => null()
93 contains
94
95 procedure :: init => nonlocal_pseudopotential_init
96
97 procedure :: build => nonlocal_pseudopotential_build_proj
98
100
101 procedure :: has_self_overlap => nonlocal_pseudopotential_self_overlap
102
103 procedure :: dstart => dnonlocal_pseudopotential_start
104
105 procedure :: zstart => znonlocal_pseudopotential_start
106
107 procedure :: dfinish => dnonlocal_pseudopotential_finish
108
109 procedure :: zfinish => znonlocal_pseudopotential_finish
110
111 procedure :: dforce => dnonlocal_pseudopotential_force
112
113 procedure :: zforce => znonlocal_pseudopotential_force
115 procedure :: dposition_commutator => dnonlocal_pseudopotential_position_commutator
116
117 procedure :: zposition_commutator => znonlocal_pseudopotential_position_commutator
118
119 procedure :: dr_vn_local => dnonlocal_pseudopotential_r_vnlocal
120
121 procedure :: zr_vn_local => znonlocal_pseudopotential_r_vnlocal
122
124
126 !
127 type projection_t
128 private
129 real(real64), allocatable :: dprojection(:, :)
130 complex(real64), allocatable :: zprojection(:, :)
131 type(accel_mem_t) :: buff_projection
132 type(accel_mem_t) :: buff_spin_to_phase
133 type(accel_mem_t), allocatable :: buff_phasepsi(:)
134 type(accel_mem_t), allocatable :: buff_projection_temp(:)
135 end type projection_t
136
137contains
138
139 ! ---------------------------------------------------------
142 subroutine nonlocal_pseudopotential_init(this)
143 class(nonlocal_pseudopotential_t), intent(inout) :: this
144
146
147 this%apply_projector_matrices = .false.
148 this%has_non_local_potential = .false.
149 this%nprojector_matrices = 0
150
151 this%projector_self_overlap = .false.
152
154 end subroutine nonlocal_pseudopotential_init
156
157 !--------------------------------------------------------
160 class(nonlocal_pseudopotential_t), intent(inout) :: this
162 integer :: iproj
166 if (allocated(this%projector_matrices)) then
169 call accel_release_buffer(this%buff_offsets)
170 call accel_release_buffer(this%buff_matrices)
171 call accel_release_buffer(this%buff_maps)
172 call accel_release_buffer(this%buff_scals)
173 call accel_release_buffer(this%buff_position)
174 call accel_release_buffer(this%buff_pos)
175 call accel_release_buffer(this%buff_invmap)
176 if (this%projector_mix) call accel_release_buffer(this%buff_mix)
177 if (allocated(this%projector_phases)) call accel_release_buffer(this%buff_projector_phases)
178 end if
180 do iproj = 1, this%nprojector_matrices
181 call projector_matrix_deallocate(this%projector_matrices(iproj))
182 end do
183 safe_deallocate_a(this%regions)
184 safe_deallocate_a(this%projector_matrices)
185 safe_deallocate_a(this%projector_phases)
186 safe_deallocate_a(this%projector_to_atom)
187 end if
191
192 !-----------------------------------------------------------------
199 subroutine nonlocal_pseudopotential_build_proj(this, space, mesh, epot)
200 class(nonlocal_pseudopotential_t), target, intent(inout) :: this
201 class(space_t), intent(in) :: space
202 class(mesh_t), intent(in) :: mesh
203 type(epot_t), target, intent(in) :: epot
205 integer :: iatom, iproj, ll, lmax, lloc, mm, ic, jc
206 integer :: nmat, imat, ip, iorder
207 integer :: nregion, jatom, katom, iregion
208 integer, allocatable :: order(:), head(:), region_count(:)
209 logical, allocatable :: atom_counted(:)
210 logical :: overlap
211 type(projector_matrix_t), pointer :: pmat
212 type(kb_projector_t), pointer :: kb_p
213 type(rkb_projector_t), pointer :: rkb_p
214 type(hgh_projector_t), pointer :: hgh_p
215
217
218 call profiling_in("ATOM_COLORING")
219
220 ! this is most likely a very inefficient algorithm, O(natom**2) or
221 ! O(natom**3), probably it should be replaced by something better.
223 safe_allocate(order(1:epot%natoms)) ! order(iregion) = ?
224 safe_allocate(head(1:epot%natoms + 1)) ! head(iregion) points to the first atom in region iregion
225 safe_allocate(region_count(1:epot%natoms)) ! region_count(iregion): number of atoms in that region
226 safe_allocate(atom_counted(1:epot%natoms))
228 this%projector_self_overlap = .false.
229 atom_counted = .false.
230 order = -1
231
232 head(1) = 1
233 nregion = 0
234 do
235 nregion = nregion + 1
236 assert(nregion <= epot%natoms)
237
238 region_count(nregion) = 0
239
240 do iatom = 1, epot%natoms
241 if (atom_counted(iatom)) cycle
242
243 overlap = .false.
244
245 if (.not. projector_is(epot%proj(iatom), proj_none)) then
246 assert(associated(epot%proj(iatom)%sphere%mesh))
247 do jatom = 1, region_count(nregion)
248 katom = order(head(nregion) + jatom - 1)
249 if (projector_is(epot%proj(katom), proj_none)) cycle
250 overlap = submesh_overlap(epot%proj(iatom)%sphere, epot%proj(katom)%sphere, space)
251 if (overlap) exit
252 end do
253 end if
254
255 if (.not. overlap) then
256 ! iatom did not overlap with any previously counted atoms:
257 ! iatom will be added to the current region
258 region_count(nregion) = region_count(nregion) + 1
259 order(head(nregion) - 1 + region_count(nregion)) = iatom
260 atom_counted(iatom) = .true.
261 end if
262
263 end do
264
265 head(nregion + 1) = head(nregion) + region_count(nregion)
266
267 if (all(atom_counted)) exit
268 end do
269
270 safe_deallocate_a(atom_counted)
271 safe_deallocate_a(region_count)
272
273 call messages_write('The atoms can be separated in ')
274 call messages_write(nregion)
275 call messages_write(' non-overlapping groups.')
276 call messages_info(debug_only=.true.)
277
278 do iregion = 1, nregion
279 do iatom = head(iregion), head(iregion + 1) - 1
280 if (.not. projector_is(epot%proj(order(iatom)), proj_kb)) cycle
281 do jatom = head(iregion), iatom - 1
282 if (.not. projector_is(epot%proj(order(jatom)), proj_kb)) cycle
283 assert(.not. submesh_overlap(epot%proj(order(iatom))%sphere, epot%proj(order(jatom))%sphere, space))
284 end do
285 end do
286 end do
287
288 call profiling_out("ATOM_COLORING")
289
290 ! deallocate previous projectors
291 call this%end()
293 ! count projectors
294 this%nprojector_matrices = 0
295 this%apply_projector_matrices = .false.
296 this%has_non_local_potential = .false.
297 this%nregions = nregion
298
299 !We determine if we have only local potential or not.
300 do iorder = 1, epot%natoms
301 iatom = order(iorder)
302
303 if (.not. projector_is_null(epot%proj(iatom))) then
304 this%has_non_local_potential = .true.
305 exit
306 end if
307 end do
308
309 do iorder = 1, epot%natoms
310 iatom = order(iorder)
311
312 if (.not. projector_is_null(epot%proj(iatom))) then
313 this%nprojector_matrices = this%nprojector_matrices + 1
314 this%apply_projector_matrices = .true.
315 end if
316 end do
317
318 ! This is currently the only not supported case
319 if (mesh%use_curvilinear) this%apply_projector_matrices = .false.
320
321 if (.not. this%apply_projector_matrices) then
322 safe_deallocate_a(order)
323 safe_deallocate_a(head)
324
326 return
327 end if
328
329
330 safe_allocate(this%projector_matrices(1:this%nprojector_matrices))
331 safe_allocate(this%regions(1:this%nprojector_matrices + 1))
332 safe_allocate(this%projector_to_atom(1:epot%natoms))
333
334 this%full_projection_size = 0
335 this%regions(this%nregions + 1) = this%nprojector_matrices + 1
336
337 this%projector_mix = .false.
338
339 iproj = 0
340 do iregion = 1, this%nregions
341 this%regions(iregion) = iproj + 1
342 do iorder = head(iregion), head(iregion + 1) - 1
343
344 iatom = order(iorder)
345
346 if (projector_is(epot%proj(iatom), proj_none)) cycle
347
348 iproj = iproj + 1
349
350 pmat => this%projector_matrices(iproj)
351
352 this%projector_to_atom(iproj) = iatom
353
354 lmax = epot%proj(iatom)%lmax
355 lloc = epot%proj(iatom)%lloc
356
357 if (projector_is(epot%proj(iatom), proj_kb)) then
358
359 ! count the number of projectors for this matrix
360 nmat = 0
361 do ll = 0, lmax
362 if (ll == lloc) cycle
363 do mm = -ll, ll
364 nmat = nmat + epot%proj(iatom)%kb_p(ll, mm)%n_c
365 end do
366 end do
367
368 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, has_mix_matrix = .false.)
369
370 ! generate the matrix
371 pmat%dprojectors = m_zero
372 imat = 1
373 do ll = 0, lmax
374 if (ll == lloc) cycle
375 do mm = -ll, ll
376 kb_p => epot%proj(iatom)%kb_p(ll, mm)
377 do ic = 1, kb_p%n_c
378 call lalg_copy(pmat%npoints, kb_p%p(:, ic), pmat%dprojectors(:, imat))
379 pmat%scal(imat) = kb_p%e(ic)*mesh%vol_pp(1)
380 imat = imat + 1
381 end do
382 end do
383 end do
384
385 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
386
387 else if (projector_is(epot%proj(iatom), proj_hgh)) then
388
389 this%projector_mix = .true.
390
391 ! count the number of projectors for this matrix
392 nmat = 0
393 do ll = 0, lmax
394 if (ll == lloc) cycle
395 do mm = -ll, ll
396 nmat = nmat + 3
397 end do
398 end do
399
400 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, &
401 has_mix_matrix = .true., is_cmplx = (epot%reltype == spin_orbit))
402
403 ! generate the matrix
404 if (epot%reltype == spin_orbit) then
405 pmat%zprojectors = m_zero
406 pmat%zmix = m_zero
407 else
408 pmat%dprojectors = m_zero
409 pmat%dmix = m_zero
410 end if
411
412 imat = 1
413 do ll = 0, lmax
414 if (ll == lloc) cycle
415 do mm = -ll, ll
416 hgh_p => epot%proj(iatom)%hgh_p(ll, mm)
417
418 ! HGH pseudos mix different components, so we need to
419 ! generate a matrix that mixes the projections
420 if (epot%reltype == spin_orbit .or. epot%reltype == fully_relativistic_zora) then
421 do ic = 1, 3
422 do jc = 1, 3
423 pmat%zmix(imat - 1 + ic, imat - 1 + jc, 1) = hgh_p%h(ic, jc) + m_half*mm*hgh_p%k(ic, jc)
424 pmat%zmix(imat - 1 + ic, imat - 1 + jc, 2) = hgh_p%h(ic, jc) - m_half*mm*hgh_p%k(ic, jc)
425
426 if (mm < ll) then
427 pmat%zmix(imat - 1 + ic, imat + 3 - 1 + jc, 3) = m_half*hgh_p%k(ic, jc) * &
428 sqrt(real(ll*(ll+1)-mm*(mm+1), real64))
429 end if
430
431 if (-mm < ll) then
432 pmat%zmix(imat - 1 + ic, imat - 3 - 1 + jc, 4) = m_half*hgh_p%k(ic, jc) * &
433 sqrt(real(ll*(ll+1)-mm*(mm-1), real64))
434 end if
435 end do
436 end do
437 else
438 do ic = 1, 3
439 do jc = 1, 3
440 pmat%dmix(imat - 1 + ic, imat - 1 + jc) = hgh_p%h(ic, jc)
441 end do
442 end do
443 end if
444
445 do ic = 1, 3
446 if (epot%reltype == spin_orbit .or. epot%reltype == fully_relativistic_zora) then
447 call lalg_copy(pmat%npoints, hgh_p%zp(:, ic), pmat%zprojectors(:, imat))
448 else
449 call lalg_copy(pmat%npoints, hgh_p%dp(:, ic), pmat%dprojectors(:, imat))
450 end if
451 pmat%scal(imat) = mesh%volume_element
452 imat = imat + 1
453 end do
454
455 end do
456 end do
457
458 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
459
460 else if (projector_is(epot%proj(iatom), proj_rkb)) then
461 assert(epot%reltype == spin_orbit)
462
463 this%projector_mix = .true.
464
465 ! count the number of projectors for this matrix
466 nmat = 0
467 if (lloc /= 0) nmat = nmat + epot%proj(iatom)%kb_p(1, 1)%n_c
468
469 do ll = 1, lmax
470 if (ll == lloc) cycle
471 do mm = -ll, ll
472 nmat = nmat + epot%proj(iatom)%rkb_p(ll, mm)%n_c
473 end do
474 end do
475
476 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, &
477 has_mix_matrix = .true., is_cmplx = .true.)
478
479 pmat%zprojectors = m_zero
480 pmat%zmix = m_zero
481
482 imat = 1
483 if (lloc /= 0) then
484 kb_p => epot%proj(iatom)%kb_p(1, 1)
485
486 do ic = 1, kb_p%n_c
487 pmat%zmix(ic, ic, 1:2) = kb_p%e(ic)
488 do ip = 1, pmat%npoints
489 pmat%zprojectors(ip, ic) = kb_p%p(ip, ic)
490 end do
491 pmat%scal(ic) = mesh%volume_element
492 end do
493 imat = kb_p%n_c + 1
494 nullify(kb_p)
495 end if
496
497 do ll = 1, lmax
498 if (ll == lloc) cycle
499 do mm = -ll, ll
500 rkb_p => epot%proj(iatom)%rkb_p(ll, mm)
501
502 ! See rkb_projector.F90 for understanding the indices
503 do ic = 0, rkb_p%n_c/2-1
504 pmat%zmix(imat + ic*2, imat + ic*2, 1) = rkb_p%f(ic*2+1, 1, 1)
505 pmat%zmix(imat + ic*2, imat + ic*2, 2) = rkb_p%f(ic*2+1, 2, 2)
506
507 pmat%zmix(imat + ic*2+1, imat + ic*2+1, 1) = rkb_p%f(ic*2+2, 1, 1)
508 pmat%zmix(imat + ic*2+1, imat + ic*2+1, 2) = rkb_p%f(ic*2+2, 2, 2)
509
510 if (mm < ll) then
511 pmat%zmix(imat + ic*2+rkb_p%n_c, imat + ic*2, 4) = rkb_p%f(ic*2+1, 2, 1)
512 pmat%zmix(imat + ic*2+1+rkb_p%n_c, imat + ic*2+1, 4) = rkb_p%f(ic*2+2, 2, 1)
513 end if
514
515 if (-mm < ll) then
516 pmat%zmix(imat + ic*2-rkb_p%n_c, imat + ic*2, 3) = rkb_p%f(ic*2+1, 1, 2)
517 pmat%zmix(imat + ic*2+1-rkb_p%n_c, imat + ic*2+1, 3) = rkb_p%f(ic*2+2, 1, 2)
518 end if
519 end do
520
521 do ic = 1, rkb_p%n_c
522 call lalg_copy(pmat%npoints, rkb_p%ket(:, ic, 1, 1), pmat%zprojectors(:, imat))
523 pmat%scal(imat) = mesh%volume_element
524 imat = imat + 1
525 end do
526 end do
527
528 nullify(rkb_p)
529 end do
530
531 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
532
533 else
534 cycle
535 end if
536
537 pmat%map => epot%proj(iatom)%sphere%map
538 pmat%position => epot%proj(iatom)%sphere%rel_x
539
540 pmat%regions = epot%proj(iatom)%sphere%regions
541
542 this%full_projection_size = this%full_projection_size + pmat%nprojs
543
544 end do
545 end do
546
547 if (mesh%parallel_in_domains) then
548 call mesh%mpi_grp%allreduce_inplace(this%projector_self_overlap, 1, mpi_logical, mpi_lor)
549 end if
550
551 safe_deallocate_a(order)
552 safe_deallocate_a(head)
553
554 this%total_points = 0
555 this%max_npoints = 0
556 this%max_nprojs = 0
557 do imat = 1, this%nprojector_matrices
558 pmat => this%projector_matrices(imat)
559
560 this%max_npoints = max(this%max_npoints, pmat%npoints)
561 this%max_nprojs = max(this%max_nprojs, pmat%nprojs)
562 this%total_points = this%total_points + pmat%npoints
563 end do
564
565 if (accel_is_enabled()) call build_accel()
566
568
569 contains
570
571 subroutine build_accel()
572
573 integer :: matrix_size, scal_size
574 integer, allocatable :: cnt(:), invmap(:, :), invmap2(:), pos(:)
575 integer, allocatable :: offsets(:, :)
576 integer, parameter :: OFFSET_SIZE = 6 ! also defined in share/opencl/projectors.cl
577 integer, parameter :: POINTS = 1, projs = 2, matrix = 3, map = 4, scal = 5, mix = 6 ! update OFFSET_SIZE
578 integer :: ip, is, ii, ipos, mix_offset
579
581
582 safe_allocate(offsets(1:offset_size, 1:this%nprojector_matrices))
583 safe_allocate(cnt(1:mesh%np))
584
585 cnt = 0
586
587 ! Here we construct the offsets for accessing various arrays within the GPU kernels.
588 ! The offset(:,:) array contains a number of sizes and offsets, describing how to address the arrays.
589 ! This allows to transfer all these numbers to the GPU in one memory transfer.
590 !
591 ! For each projection matrix (addressed by imap), we have:
592 !
593 ! offset(POINTS, imap) : number of points of the sphere imap
594 ! offset(PROJS, imap) : number of projectors for imap
595 ! offset(MATRIX, imap) : address offset: cumulative of pmat%npoints * pmat%nprojs
596 ! offset(MAP, imap) : address offset: cumulative of pmat%npoints for each imap
597 ! offset(SCAL, imap) : address_offset: cumulative of pmat%nprojs
598 ! offset(MIX, imap) : address_offset: cumulative of pmat%nprojs**2 or 4*pmat%nprojs**2 for complex mixing
599
600 ! first we count
601 matrix_size = 0
602 this%total_points = 0
603 scal_size = 0
604 this%max_npoints = 0
605 this%max_nprojs = 0
606 mix_offset = 0
607 do imat = 1, this%nprojector_matrices
608 pmat => this%projector_matrices(imat)
609
610 this%max_npoints = max(this%max_npoints, pmat%npoints)
611 this%max_nprojs = max(this%max_nprojs, pmat%nprojs)
612
613 offsets(points, imat) = pmat%npoints
614 offsets(projs, imat) = pmat%nprojs
615
616 offsets(matrix, imat) = matrix_size
617 matrix_size = matrix_size + pmat%npoints*pmat%nprojs
618
619 offsets(map, imat) = this%total_points
620 this%total_points = this%total_points + pmat%npoints
621
622 offsets(scal, imat) = scal_size
623 scal_size = scal_size + pmat%nprojs
624
625 offsets(mix, imat) = mix_offset
626 if (allocated(pmat%dmix)) then
627 mix_offset = mix_offset + pmat%nprojs**2
628 else if (allocated(pmat%zmix)) then
629 mix_offset = mix_offset + 4 * pmat%nprojs**2
630 else
631 offsets(mix, imat) = -1
632 end if
633
634 do is = 1, pmat%npoints
635 ip = pmat%map(is)
636 cnt(ip) = cnt(ip) + 1
637 end do
638 end do
639
640 safe_allocate(invmap(1:max(maxval(cnt), 1), 1:mesh%np))
641 safe_allocate(invmap2(1:max(maxval(cnt)*mesh%np, 1)))
642 safe_allocate(pos(1:mesh%np + 1))
643
644 cnt = 0
645 ii = 0
646 do imat = 1, this%nprojector_matrices
647 pmat => this%projector_matrices(imat)
648 do is = 1, pmat%npoints
649 ip = pmat%map(is)
650 cnt(ip) = cnt(ip) + 1
651 invmap(cnt(ip), ip) = ii
652 ii = ii + 1
653 end do
654 end do
655
656 ipos = 0
657 pos(1) = 0
658 do ip = 1, mesh%np
659 do ii = 1, cnt(ip)
660 ipos = ipos + 1
661 invmap2(ipos) = invmap(ii, ip)
662 end do
663 pos(ip + 1) = ipos
664 end do
665
666 ! allocate
667 if (this%projector_matrices(1)%is_cmplx) then
668 call accel_create_buffer(this%buff_matrices, accel_mem_read_only, type_cmplx, matrix_size)
669 else
670 call accel_create_buffer(this%buff_matrices, accel_mem_read_only, type_float, matrix_size)
671 end if
672 call accel_create_buffer(this%buff_maps, accel_mem_read_only, type_integer, this%total_points)
673 call accel_create_buffer(this%buff_position, accel_mem_read_only, type_float, 3*this%total_points)
674 call accel_create_buffer(this%buff_scals, accel_mem_read_only, type_float, scal_size)
675
676 if (mix_offset > 0) then
677 if (allocated(this%projector_matrices(1)%zmix)) then
678 call accel_create_buffer(this%buff_mix, accel_mem_read_only, type_cmplx, mix_offset)
679 else
680 call accel_create_buffer(this%buff_mix, accel_mem_read_only, type_float, mix_offset)
681 end if
682 end if
683
684 ! now copy
685 do imat = 1, this%nprojector_matrices
686 pmat => this%projector_matrices(imat)
687 ! in parallel some spheres might not have points
688 if (pmat%npoints > 0) then
689 if (pmat%is_cmplx) then
690 call accel_write_buffer(this%buff_matrices, pmat%nprojs*pmat%npoints, pmat%zprojectors, offset = offsets(matrix, imat))
691 else
692 call accel_write_buffer(this%buff_matrices, pmat%nprojs*pmat%npoints, pmat%dprojectors, offset = offsets(matrix, imat))
693 end if
694 call accel_write_buffer(this%buff_maps, pmat%npoints, pmat%map, offset = offsets(map, imat))
695 call accel_write_buffer(this%buff_position, 3*pmat%npoints, pmat%position, offset = 3*offsets(map, imat))
696 end if
697 call accel_write_buffer(this%buff_scals, pmat%nprojs, pmat%scal, offset = offsets(scal, imat))
698 if (offsets(mix, imat) /= -1) then
699 if (allocated(pmat%zmix)) then
700 call accel_write_buffer(this%buff_mix, 4*pmat%nprojs**2, pmat%zmix, offset = offsets(mix, imat))
701 else
702 call accel_write_buffer(this%buff_mix, pmat%nprojs**2, pmat%dmix, offset = offsets(mix, imat))
703 end if
704 end if
705 end do
706
707 ! write the offsets
708 call accel_create_buffer(this%buff_offsets, accel_mem_read_only, type_integer, offset_size*this%nprojector_matrices)
709 call accel_write_buffer(this%buff_offsets, offset_size*this%nprojector_matrices, offsets)
710
711 ! the inverse map
712 call accel_create_buffer(this%buff_pos, accel_mem_read_only, type_integer, mesh%np + 1)
713 call accel_write_buffer(this%buff_pos, mesh%np + 1, pos)
714
715 call accel_create_buffer(this%buff_invmap, accel_mem_read_only, type_integer, ipos)
716 call accel_write_buffer(this%buff_invmap, ipos, invmap2)
717
718 safe_deallocate_a(offsets)
719 safe_deallocate_a(cnt)
720 safe_deallocate_a(invmap)
721 safe_deallocate_a(invmap2)
722 safe_deallocate_a(pos)
723
725 end subroutine build_accel
726
728
729 ! ----------------------------------------------------------------------------------
731 !
732 logical pure function nonlocal_pseudopotential_self_overlap(this) result(projector_self_overlap)
733 class(nonlocal_pseudopotential_t), intent(in) :: this
734
735 projector_self_overlap = this%projector_self_overlap
737
738#include "undef.F90"
739#include "real.F90"
740#include "nonlocal_pseudopotential_inc.F90"
741
742#include "undef.F90"
743#include "complex.F90"
744#include "nonlocal_pseudopotential_inc.F90"
745
747
748!! Local Variables:
749!! mode: f90
750!! coding: utf-8
751!! End:
Copies a vector x, to a vector y.
Definition: lalg_basic.F90:186
double sqrt(double __x) __attribute__((__nothrow__
subroutine, public accel_release_buffer(this, async)
Definition: accel.F90:1331
pure logical function, public accel_is_enabled()
Definition: accel.F90:427
integer, parameter, public accel_mem_read_only
Definition: accel.F90:187
type(accel_kernel_t), pointer head
Definition: accel.F90:421
This module implements batches of mesh functions.
Definition: batch.F90:133
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:116
This module contains interfaces for BLAS routines You should not use these routines directly....
Definition: blas.F90:118
integer, parameter, public spin_orbit
Definition: epot.F90:166
integer, parameter, public fully_relativistic_zora
Definition: epot.F90:166
real(real64), parameter, public m_zero
Definition: global.F90:188
real(real64), parameter, public m_half
Definition: global.F90:194
This module is intended to contain "only mathematical" functions and procedures.
Definition: math.F90:115
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:616
subroutine nonlocal_pseudopotential_destroy_proj(this)
Destroy the data of nonlocal_pseudopotential_t.
subroutine dnonlocal_pseudopotential_force(this, mesh, st, spiral_bnd, iqn, ndim, psi1b, psi2b, force)
calculate contribution to forces, from non-local potentials
subroutine znonlocal_pseudopotential_position_commutator(this, mesh, std, spiral_bnd, psib, commpsib, async)
apply the commutator between the non-local potential and the position to the wave functions.
subroutine znonlocal_pseudopotential_force(this, mesh, st, spiral_bnd, iqn, ndim, psi1b, psi2b, force)
calculate contribution to forces, from non-local potentials
subroutine dnonlocal_pseudopotential_start(this, mesh, std, spiral_bnd, psib, projection, async)
Start application of non-local potentials (stored in the Hamiltonian) to the wave functions.
subroutine dnonlocal_pseudopotential_finish(this, mesh, spiral_bnd, std, projection, vpsib)
finish the application of non-local potentials.
subroutine dnonlocal_pseudopotential_position_commutator(this, mesh, std, spiral_bnd, psib, commpsib, async)
apply the commutator between the non-local potential and the position to the wave functions.
subroutine znonlocal_pseudopotential_r_vnlocal(this, mesh, std, spiral_bnd, psib, commpsib)
Accumulates to commpsib the result of x V_{nl} | psib >
logical pure function nonlocal_pseudopotential_self_overlap(this)
Returns .true. if the Hamiltonian contains projectors, which overlap with themself.
subroutine nonlocal_pseudopotential_init(this)
initialize the nonlocal_pseudopotential_t object
subroutine dnonlocal_pseudopotential_r_vnlocal(this, mesh, std, spiral_bnd, psib, commpsib)
Accumulates to commpsib the result of x V_{nl} | psib >
subroutine znonlocal_pseudopotential_start(this, mesh, std, spiral_bnd, psib, projection, async)
Start application of non-local potentials (stored in the Hamiltonian) to the wave functions.
subroutine nonlocal_pseudopotential_build_proj(this, space, mesh, epot)
build the projectors for the application of pseudo-potentials
subroutine znonlocal_pseudopotential_finish(this, mesh, spiral_bnd, std, projection, vpsib)
finish the application of non-local potentials.
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
Definition: profiling.F90:623
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
Definition: profiling.F90:552
subroutine, public projector_matrix_deallocate(this)
subroutine, public projector_matrix_allocate(this, nprojs, sphere, has_mix_matrix, is_cmplx)
logical elemental function, public projector_is(p, type)
Definition: projector.F90:208
logical elemental function, public projector_is_null(p)
Definition: projector.F90:201
Definition: ps.F90:114
integer, parameter, public proj_hgh
Definition: ps.F90:167
integer, parameter, public proj_rkb
Definition: ps.F90:167
integer, parameter, public proj_none
Definition: ps.F90:167
integer, parameter, public proj_kb
Definition: ps.F90:167
This module handles spin dimensions of the states and the k-point distribution.
logical function, public submesh_overlap(sm1, sm2, space)
Definition: submesh.F90:777
type(type_t), public type_float
Definition: types.F90:133
type(type_t), public type_cmplx
Definition: types.F90:134
type(type_t), public type_integer
Definition: types.F90:135
subroutine build_accel()
Describes mesh distribution to nodes.
Definition: mesh.F90:186
Class for projections of wave functions.
A set of projectors defined on a submesh.
The rkb_projector data type holds the KB projectors build with total angular momentum eigenfunctions....
int true(void)