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
23 use batch_oct_m
25 use blas_oct_m
26 use comm_oct_m
27 use debug_oct_m
29 use epot_oct_m
30 use global_oct_m
35 use math_oct_m
36 use mesh_oct_m
38 use mpi_oct_m
42 use ps_oct_m
44 use space_oct_m
48 use types_oct_m
50
51 implicit none
52
53 private
54
55 public :: &
58
62 private
63 type(projector_matrix_t), allocatable, public :: projector_matrices(:)
64 integer, public :: nprojector_matrices
65 logical, public :: apply_projector_matrices
66 logical, public :: has_non_local_potential
67 integer :: full_projection_size
68 integer, public :: max_npoints
69 integer, public :: total_points
70 integer :: max_nprojs
71 logical :: projector_mix
72 complex(real64), allocatable, public :: projector_phases(:, :, :, :)
73 integer, allocatable, public :: projector_to_atom(:)
74 integer :: nregions
75 integer, allocatable :: regions(:)
76 integer, public :: nphase
81 type(accel_mem_t) :: buff_offsets
82 type(accel_mem_t) :: buff_matrices
83 type(accel_mem_t) :: buff_maps
84 type(accel_mem_t) :: buff_scals
85 type(accel_mem_t) :: buff_position
86 type(accel_mem_t) :: buff_pos
87 type(accel_mem_t) :: buff_invmap
88 type(accel_mem_t), public :: buff_projector_phases
89 type(accel_mem_t) :: buff_mix
90 logical :: projector_self_overlap
91 real(real64), pointer, public :: spin(:,:,:) => null()
92 contains
93
94 procedure :: init => nonlocal_pseudopotential_init
95
96 procedure :: build => nonlocal_pseudopotential_build_proj
97
99
100 procedure :: has_self_overlap => nonlocal_pseudopotential_self_overlap
101
102 procedure :: dstart => dnonlocal_pseudopotential_start
103
104 procedure :: zstart => znonlocal_pseudopotential_start
105
106 procedure :: dfinish => dnonlocal_pseudopotential_finish
107
108 procedure :: zfinish => znonlocal_pseudopotential_finish
109
110 procedure :: dforce => dnonlocal_pseudopotential_force
111
112 procedure :: zforce => znonlocal_pseudopotential_force
113
114 procedure :: dposition_commutator => dnonlocal_pseudopotential_position_commutator
115
116 procedure :: zposition_commutator => znonlocal_pseudopotential_position_commutator
117
118 procedure :: dr_vn_local => dnonlocal_pseudopotential_r_vnlocal
119
120 procedure :: zr_vn_local => znonlocal_pseudopotential_r_vnlocal
121
123
125 !
126 type projection_t
127 private
128 real(real64), allocatable :: dprojection(:, :)
129 complex(real64), allocatable :: zprojection(:, :)
130 type(accel_mem_t) :: buff_projection
131 type(accel_mem_t) :: buff_spin_to_phase
132 end type projection_t
133
134contains
135
136 ! ---------------------------------------------------------
139 subroutine nonlocal_pseudopotential_init(this)
140 class(nonlocal_pseudopotential_t), intent(inout) :: this
141
143
144 this%apply_projector_matrices = .false.
145 this%has_non_local_potential = .false.
146 this%nprojector_matrices = 0
147
148 this%projector_self_overlap = .false.
149
151 end subroutine nonlocal_pseudopotential_init
152
153
154 !--------------------------------------------------------
157 class(nonlocal_pseudopotential_t), intent(inout) :: this
159 integer :: iproj
163 if (allocated(this%projector_matrices)) then
166 call accel_release_buffer(this%buff_offsets)
167 call accel_release_buffer(this%buff_matrices)
168 call accel_release_buffer(this%buff_maps)
169 call accel_release_buffer(this%buff_scals)
170 call accel_release_buffer(this%buff_position)
171 call accel_release_buffer(this%buff_pos)
172 call accel_release_buffer(this%buff_invmap)
173 if (this%projector_mix) call accel_release_buffer(this%buff_mix)
174 if (allocated(this%projector_phases)) call accel_release_buffer(this%buff_projector_phases)
175 end if
177 do iproj = 1, this%nprojector_matrices
178 call projector_matrix_deallocate(this%projector_matrices(iproj))
179 end do
180 safe_deallocate_a(this%regions)
181 safe_deallocate_a(this%projector_matrices)
182 safe_deallocate_a(this%projector_phases)
183 safe_deallocate_a(this%projector_to_atom)
184 end if
185
188
189 !-----------------------------------------------------------------
196 subroutine nonlocal_pseudopotential_build_proj(this, space, mesh, epot)
197 class(nonlocal_pseudopotential_t), target, intent(inout) :: this
198 class(space_t), intent(in) :: space
199 class(mesh_t), intent(in) :: mesh
200 type(epot_t), target, intent(in) :: epot
202 integer :: iatom, iproj, ll, lmax, lloc, mm, ic, jc
203 integer :: nmat, imat, ip, iorder
204 integer :: nregion, jatom, katom, iregion
205 integer, allocatable :: order(:), head(:), region_count(:)
206 logical, allocatable :: atom_counted(:)
207 logical :: overlap
208 type(projector_matrix_t), pointer :: pmat
209 type(kb_projector_t), pointer :: kb_p
210 type(rkb_projector_t), pointer :: rkb_p
211 type(hgh_projector_t), pointer :: hgh_p
212
214
215 call profiling_in("ATOM_COLORING")
216
217 ! this is most likely a very inefficient algorithm, O(natom**2) or
218 ! O(natom**3), probably it should be replaced by something better.
220 safe_allocate(order(1:epot%natoms)) ! order(iregion) = ?
221 safe_allocate(head(1:epot%natoms + 1)) ! head(iregion) points to the first atom in region iregion
222 safe_allocate(region_count(1:epot%natoms)) ! region_count(iregion): number of atoms in that region
223 safe_allocate(atom_counted(1:epot%natoms))
225 this%projector_self_overlap = .false.
226 atom_counted = .false.
227 order = -1
228
229 head(1) = 1
230 nregion = 0
231 do
232 nregion = nregion + 1
233 assert(nregion <= epot%natoms)
234
235 region_count(nregion) = 0
236
237 do iatom = 1, epot%natoms
238 if (atom_counted(iatom)) cycle
239
240 overlap = .false.
241
242 if (.not. projector_is(epot%proj(iatom), proj_none)) then
243 assert(associated(epot%proj(iatom)%sphere%mesh))
244 do jatom = 1, region_count(nregion)
245 katom = order(head(nregion) + jatom - 1)
246 if (projector_is(epot%proj(katom), proj_none)) cycle
247 overlap = submesh_overlap(epot%proj(iatom)%sphere, epot%proj(katom)%sphere, space)
248 if (overlap) exit
249 end do
250 end if
251
252 if (.not. overlap) then
253 ! iatom did not overlap with any previously counted atoms:
254 ! iatom will be added to the current region
255 region_count(nregion) = region_count(nregion) + 1
256 order(head(nregion) - 1 + region_count(nregion)) = iatom
257 atom_counted(iatom) = .true.
258 end if
259
260 end do
261
262 head(nregion + 1) = head(nregion) + region_count(nregion)
263
264 if (all(atom_counted)) exit
265 end do
266
267 safe_deallocate_a(atom_counted)
268 safe_deallocate_a(region_count)
269
270 call messages_write('The atoms can be separated in ')
271 call messages_write(nregion)
272 call messages_write(' non-overlapping groups.')
273 call messages_info(debug_only=.true.)
274
275 do iregion = 1, nregion
276 do iatom = head(iregion), head(iregion + 1) - 1
277 if (.not. projector_is(epot%proj(order(iatom)), proj_kb)) cycle
278 do jatom = head(iregion), iatom - 1
279 if (.not. projector_is(epot%proj(order(jatom)), proj_kb)) cycle
280 assert(.not. submesh_overlap(epot%proj(order(iatom))%sphere, epot%proj(order(jatom))%sphere, space))
281 end do
282 end do
283 end do
284
285 call profiling_out("ATOM_COLORING")
286
287 ! deallocate previous projectors
288 call this%end()
290 ! count projectors
291 this%nprojector_matrices = 0
292 this%apply_projector_matrices = .false.
293 this%has_non_local_potential = .false.
294 this%nregions = nregion
295
296 !We determine if we have only local potential or not.
297 do iorder = 1, epot%natoms
298 iatom = order(iorder)
299
300 if (.not. projector_is_null(epot%proj(iatom))) then
301 this%has_non_local_potential = .true.
302 exit
303 end if
304 end do
305
306 do iorder = 1, epot%natoms
307 iatom = order(iorder)
308
309 if (.not. projector_is_null(epot%proj(iatom))) then
310 this%nprojector_matrices = this%nprojector_matrices + 1
311 this%apply_projector_matrices = .true.
312 end if
313 end do
314
315 ! This is currently the only not supported case
316 if (mesh%use_curvilinear) this%apply_projector_matrices = .false.
317
318 if (.not. this%apply_projector_matrices) then
319 safe_deallocate_a(order)
320 safe_deallocate_a(head)
321
323 return
324 end if
325
326
327 safe_allocate(this%projector_matrices(1:this%nprojector_matrices))
328 safe_allocate(this%regions(1:this%nprojector_matrices + 1))
329 safe_allocate(this%projector_to_atom(1:epot%natoms))
330
331 this%full_projection_size = 0
332 this%regions(this%nregions + 1) = this%nprojector_matrices + 1
333
334 this%projector_mix = .false.
335
336 iproj = 0
337 do iregion = 1, this%nregions
338 this%regions(iregion) = iproj + 1
339 do iorder = head(iregion), head(iregion + 1) - 1
340
341 iatom = order(iorder)
342
343 if (projector_is(epot%proj(iatom), proj_none)) cycle
344
345 iproj = iproj + 1
346
347 pmat => this%projector_matrices(iproj)
348
349 this%projector_to_atom(iproj) = iatom
350
351 lmax = epot%proj(iatom)%lmax
352 lloc = epot%proj(iatom)%lloc
353
354 if (projector_is(epot%proj(iatom), proj_kb)) then
355
356 ! count the number of projectors for this matrix
357 nmat = 0
358 do ll = 0, lmax
359 if (ll == lloc) cycle
360 do mm = -ll, ll
361 nmat = nmat + epot%proj(iatom)%kb_p(ll, mm)%n_c
362 end do
363 end do
364
365 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, has_mix_matrix = .false.)
366
367 ! generate the matrix
368 pmat%dprojectors = m_zero
369 imat = 1
370 do ll = 0, lmax
371 if (ll == lloc) cycle
372 do mm = -ll, ll
373 kb_p => epot%proj(iatom)%kb_p(ll, mm)
374 do ic = 1, kb_p%n_c
375 call lalg_copy(pmat%npoints, kb_p%p(:, ic), pmat%dprojectors(:, imat))
376 pmat%scal(imat) = kb_p%e(ic)*mesh%vol_pp(1)
377 imat = imat + 1
378 end do
379 end do
380 end do
381
382 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
383
384 else if (projector_is(epot%proj(iatom), proj_hgh)) then
385
386 this%projector_mix = .true.
387
388 ! count the number of projectors for this matrix
389 nmat = 0
390 do ll = 0, lmax
391 if (ll == lloc) cycle
392 do mm = -ll, ll
393 nmat = nmat + 3
394 end do
395 end do
396
397 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, &
398 has_mix_matrix = .true., is_cmplx = (epot%reltype == spin_orbit))
399
400 ! generate the matrix
401 if (epot%reltype == spin_orbit) then
402 pmat%zprojectors = m_zero
403 pmat%zmix = m_zero
404 else
405 pmat%dprojectors = m_zero
406 pmat%dmix = m_zero
407 end if
408
409 imat = 1
410 do ll = 0, lmax
411 if (ll == lloc) cycle
412 do mm = -ll, ll
413 hgh_p => epot%proj(iatom)%hgh_p(ll, mm)
414
415 ! HGH pseudos mix different components, so we need to
416 ! generate a matrix that mixes the projections
417 if (epot%reltype == spin_orbit .or. epot%reltype == fully_relativistic_zora) then
418 do ic = 1, 3
419 do jc = 1, 3
420 pmat%zmix(imat - 1 + ic, imat - 1 + jc, 1) = hgh_p%h(ic, jc) + m_half*mm*hgh_p%k(ic, jc)
421 pmat%zmix(imat - 1 + ic, imat - 1 + jc, 2) = hgh_p%h(ic, jc) - m_half*mm*hgh_p%k(ic, jc)
422
423 if (mm < ll) then
424 pmat%zmix(imat - 1 + ic, imat + 3 - 1 + jc, 3) = m_half*hgh_p%k(ic, jc) * &
425 sqrt(real(ll*(ll+1)-mm*(mm+1), real64))
426 end if
427
428 if (-mm < ll) then
429 pmat%zmix(imat - 1 + ic, imat - 3 - 1 + jc, 4) = m_half*hgh_p%k(ic, jc) * &
430 sqrt(real(ll*(ll+1)-mm*(mm-1), real64))
431 end if
432 end do
433 end do
434 else
435 do ic = 1, 3
436 do jc = 1, 3
437 pmat%dmix(imat - 1 + ic, imat - 1 + jc) = hgh_p%h(ic, jc)
438 end do
439 end do
440 end if
441
442 do ic = 1, 3
443 if (epot%reltype == spin_orbit .or. epot%reltype == fully_relativistic_zora) then
444 call lalg_copy(pmat%npoints, hgh_p%zp(:, ic), pmat%zprojectors(:, imat))
445 else
446 call lalg_copy(pmat%npoints, hgh_p%dp(:, ic), pmat%dprojectors(:, imat))
447 end if
448 pmat%scal(imat) = mesh%volume_element
449 imat = imat + 1
450 end do
451
452 end do
453 end do
454
455 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
456
457 else if (projector_is(epot%proj(iatom), proj_rkb)) then
458 assert(epot%reltype == spin_orbit)
459
460 this%projector_mix = .true.
461
462 ! count the number of projectors for this matrix
463 nmat = 0
464 if (lloc /= 0) nmat = nmat + epot%proj(iatom)%kb_p(1, 1)%n_c
465
466 do ll = 1, lmax
467 if (ll == lloc) cycle
468 do mm = -ll, ll
469 nmat = nmat + epot%proj(iatom)%rkb_p(ll, mm)%n_c
470 end do
471 end do
472
473 call projector_matrix_allocate(pmat, nmat, epot%proj(iatom)%sphere, &
474 has_mix_matrix = .true., is_cmplx = .true.)
475
476 pmat%zprojectors = m_zero
477 pmat%zmix = m_zero
478
479 imat = 1
480 if (lloc /= 0) then
481 kb_p => epot%proj(iatom)%kb_p(1, 1)
482
483 do ic = 1, kb_p%n_c
484 pmat%zmix(ic, ic, 1:2) = kb_p%e(ic)
485 do ip = 1, pmat%npoints
486 pmat%zprojectors(ip, ic) = kb_p%p(ip, ic)
487 end do
488 pmat%scal(ic) = mesh%volume_element
489 end do
490 imat = kb_p%n_c + 1
491 nullify(kb_p)
492 end if
493
494 do ll = 1, lmax
495 if (ll == lloc) cycle
496 do mm = -ll, ll
497 rkb_p => epot%proj(iatom)%rkb_p(ll, mm)
498
499 ! See rkb_projector.F90 for understanding the indices
500 do ic = 0, rkb_p%n_c/2-1
501 pmat%zmix(imat + ic*2, imat + ic*2, 1) = rkb_p%f(ic*2+1, 1, 1)
502 pmat%zmix(imat + ic*2, imat + ic*2, 2) = rkb_p%f(ic*2+1, 2, 2)
503
504 pmat%zmix(imat + ic*2+1, imat + ic*2+1, 1) = rkb_p%f(ic*2+2, 1, 1)
505 pmat%zmix(imat + ic*2+1, imat + ic*2+1, 2) = rkb_p%f(ic*2+2, 2, 2)
506
507 if (mm < ll) then
508 pmat%zmix(imat + ic*2+rkb_p%n_c, imat + ic*2, 4) = rkb_p%f(ic*2+1, 2, 1)
509 pmat%zmix(imat + ic*2+1+rkb_p%n_c, imat + ic*2+1, 4) = rkb_p%f(ic*2+2, 2, 1)
510 end if
511
512 if (-mm < ll) then
513 pmat%zmix(imat + ic*2-rkb_p%n_c, imat + ic*2, 3) = rkb_p%f(ic*2+1, 1, 2)
514 pmat%zmix(imat + ic*2+1-rkb_p%n_c, imat + ic*2+1, 3) = rkb_p%f(ic*2+2, 1, 2)
515 end if
516 end do
517
518 do ic = 1, rkb_p%n_c
519 call lalg_copy(pmat%npoints, rkb_p%ket(:, ic, 1, 1), pmat%zprojectors(:, imat))
520 pmat%scal(imat) = mesh%volume_element
521 imat = imat + 1
522 end do
523 end do
524
525 nullify(rkb_p)
526 end do
527
528 this%projector_self_overlap = this%projector_self_overlap .or. epot%proj(iatom)%sphere%overlap
529
530 else
531 cycle
532 end if
533
534 pmat%map => epot%proj(iatom)%sphere%map
535 pmat%position => epot%proj(iatom)%sphere%rel_x
536
537 pmat%regions = epot%proj(iatom)%sphere%regions
538
539 this%full_projection_size = this%full_projection_size + pmat%nprojs
540
541 end do
542 end do
543
544 if (mesh%parallel_in_domains) then
545 call mesh%mpi_grp%allreduce_inplace(this%projector_self_overlap, 1, mpi_logical, mpi_lor)
546 end if
547
548 safe_deallocate_a(order)
549 safe_deallocate_a(head)
550
551 this%total_points = 0
552 this%max_npoints = 0
553 this%max_nprojs = 0
554 do imat = 1, this%nprojector_matrices
555 pmat => this%projector_matrices(imat)
556
557 this%max_npoints = max(this%max_npoints, pmat%npoints)
558 this%max_nprojs = max(this%max_nprojs, pmat%nprojs)
559 this%total_points = this%total_points + pmat%npoints
560 end do
561
562 if (accel_is_enabled()) call build_accel()
563
565
566 contains
567
568 subroutine build_accel()
569
570 integer :: matrix_size, scal_size
571 integer, allocatable :: cnt(:), invmap(:, :), invmap2(:), pos(:)
572 integer, allocatable :: offsets(:, :)
573 integer, parameter :: OFFSET_SIZE = 6 ! also defined in share/opencl/projectors.cl
574 integer, parameter :: POINTS = 1, projs = 2, matrix = 3, map = 4, scal = 5, mix = 6 ! update OFFSET_SIZE
575 integer :: ip, is, ii, ipos, mix_offset
576
578
579 safe_allocate(offsets(1:offset_size, 1:this%nprojector_matrices))
580 safe_allocate(cnt(1:mesh%np))
581
582 cnt = 0
583
584 ! Here we construct the offsets for accessing various arrays within the GPU kernels.
585 ! The offset(:,:) array contains a number of sizes and offsets, describing how to address the arrays.
586 ! This allows to transfer all these numbers to the GPU in one memory transfer.
587 !
588 ! For each projection matrix (addressed by imap), we have:
589 !
590 ! offset(POINTS, imap) : number of points of the sphere imap
591 ! offset(PROJS, imap) : number of projectors for imap
592 ! offset(MATRIX, imap) : address offset: cumulative of pmat%npoints * pmat%nprojs
593 ! offset(MAP, imap) : address offset: cumulative of pmat%npoints for each imap
594 ! offset(SCAL, imap) : address_offset: cumulative of pmat%nprojs
595 ! offset(MIX, imap) : address_offset: cumulative of pmat%nprojs**2 or 4*pmat%nprojs**2 for complex mixing
596
597 ! first we count
598 matrix_size = 0
599 this%total_points = 0
600 scal_size = 0
601 this%max_npoints = 0
602 this%max_nprojs = 0
603 mix_offset = 0
604 do imat = 1, this%nprojector_matrices
605 pmat => this%projector_matrices(imat)
606
607 this%max_npoints = max(this%max_npoints, pmat%npoints)
608 this%max_nprojs = max(this%max_nprojs, pmat%nprojs)
609
610 offsets(points, imat) = pmat%npoints
611 offsets(projs, imat) = pmat%nprojs
612
613 offsets(matrix, imat) = matrix_size
614 matrix_size = matrix_size + pmat%npoints*pmat%nprojs
615
616 offsets(map, imat) = this%total_points
617 this%total_points = this%total_points + pmat%npoints
618
619 offsets(scal, imat) = scal_size
620 scal_size = scal_size + pmat%nprojs
621
622 offsets(mix, imat) = mix_offset
623 if (allocated(pmat%dmix)) then
624 mix_offset = mix_offset + pmat%nprojs**2
625 else if (allocated(pmat%zmix)) then
626 mix_offset = mix_offset + 4 * pmat%nprojs**2
627 else
628 offsets(mix, imat) = -1
629 end if
630
631 do is = 1, pmat%npoints
632 ip = pmat%map(is)
633 cnt(ip) = cnt(ip) + 1
634 end do
635 end do
636
637 safe_allocate(invmap(1:max(maxval(cnt), 1), 1:mesh%np))
638 safe_allocate(invmap2(1:max(maxval(cnt)*mesh%np, 1)))
639 safe_allocate(pos(1:mesh%np + 1))
640
641 cnt = 0
642 ii = 0
643 do imat = 1, this%nprojector_matrices
644 pmat => this%projector_matrices(imat)
645 do is = 1, pmat%npoints
646 ip = pmat%map(is)
647 cnt(ip) = cnt(ip) + 1
648 invmap(cnt(ip), ip) = ii
649 ii = ii + 1
650 end do
651 end do
652
653 ipos = 0
654 pos(1) = 0
655 do ip = 1, mesh%np
656 do ii = 1, cnt(ip)
657 ipos = ipos + 1
658 invmap2(ipos) = invmap(ii, ip)
659 end do
660 pos(ip + 1) = ipos
661 end do
662
663 ! allocate
664 if (this%projector_matrices(1)%is_cmplx) then
665 call accel_create_buffer(this%buff_matrices, accel_mem_read_only, type_cmplx, matrix_size)
666 else
667 call accel_create_buffer(this%buff_matrices, accel_mem_read_only, type_float, matrix_size)
668 end if
669 call accel_create_buffer(this%buff_maps, accel_mem_read_only, type_integer, this%total_points)
670 call accel_create_buffer(this%buff_position, accel_mem_read_only, type_float, 3*this%total_points)
671 call accel_create_buffer(this%buff_scals, accel_mem_read_only, type_float, scal_size)
672
673 if (mix_offset > 0) then
674 if (allocated(this%projector_matrices(1)%zmix)) then
675 call accel_create_buffer(this%buff_mix, accel_mem_read_only, type_cmplx, mix_offset)
676 else
677 call accel_create_buffer(this%buff_mix, accel_mem_read_only, type_float, mix_offset)
678 end if
679 end if
680
681 ! now copy
682 do imat = 1, this%nprojector_matrices
683 pmat => this%projector_matrices(imat)
684 ! in parallel some spheres might not have points
685 if (pmat%npoints > 0) then
686 if (pmat%is_cmplx) then
687 call accel_write_buffer(this%buff_matrices, pmat%nprojs*pmat%npoints, pmat%zprojectors, offset = offsets(matrix, imat))
688 else
689 call accel_write_buffer(this%buff_matrices, pmat%nprojs*pmat%npoints, pmat%dprojectors, offset = offsets(matrix, imat))
690 end if
691 call accel_write_buffer(this%buff_maps, pmat%npoints, pmat%map, offset = offsets(map, imat))
692 call accel_write_buffer(this%buff_position, 3*pmat%npoints, pmat%position, offset = 3*offsets(map, imat))
693 end if
694 call accel_write_buffer(this%buff_scals, pmat%nprojs, pmat%scal, offset = offsets(scal, imat))
695 if (offsets(mix, imat) /= -1) then
696 if (allocated(pmat%zmix)) then
697 call accel_write_buffer(this%buff_mix, 4*pmat%nprojs**2, pmat%zmix, offset = offsets(mix, imat))
698 else
699 call accel_write_buffer(this%buff_mix, pmat%nprojs**2, pmat%dmix, offset = offsets(mix, imat))
700 end if
701 end if
702 end do
703
704 ! write the offsets
705 call accel_create_buffer(this%buff_offsets, accel_mem_read_only, type_integer, offset_size*this%nprojector_matrices)
706 call accel_write_buffer(this%buff_offsets, offset_size*this%nprojector_matrices, offsets)
707
708 ! the inverse map
709 call accel_create_buffer(this%buff_pos, accel_mem_read_only, type_integer, mesh%np + 1)
710 call accel_write_buffer(this%buff_pos, mesh%np + 1, pos)
711
712 call accel_create_buffer(this%buff_invmap, accel_mem_read_only, type_integer, ipos)
713 call accel_write_buffer(this%buff_invmap, ipos, invmap2)
714
715 safe_deallocate_a(offsets)
716 safe_deallocate_a(cnt)
717 safe_deallocate_a(invmap)
718 safe_deallocate_a(invmap2)
719 safe_deallocate_a(pos)
720
722 end subroutine build_accel
723
725
726 ! ----------------------------------------------------------------------------------
728 !
729 logical pure function nonlocal_pseudopotential_self_overlap(this) result(projector_self_overlap)
730 class(nonlocal_pseudopotential_t), intent(in) :: this
731
732 projector_self_overlap = this%projector_self_overlap
734
735#include "undef.F90"
736#include "real.F90"
737#include "nonlocal_pseudopotential_inc.F90"
738
739#include "undef.F90"
740#include "complex.F90"
741#include "nonlocal_pseudopotential_inc.F90"
742
744
745!! Local Variables:
746!! mode: f90
747!! coding: utf-8
748!! 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)
Definition: accel.F90:1246
pure logical function, public accel_is_enabled()
Definition: accel.F90:400
integer, parameter, public accel_mem_read_only
Definition: accel.F90:183
type(accel_kernel_t), pointer head
Definition: accel.F90:394
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:776
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)