1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch 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 21module hgh_projector_oct_m 22 use atom_oct_m 23 use blas_oct_m 24 use comm_oct_m 25 use global_oct_m 26 use hardware_oct_m 27 use mesh_oct_m 28 use messages_oct_m 29 use profiling_oct_m 30 use ps_oct_m 31 use species_oct_m 32 use submesh_oct_m 33 34 implicit none 35 36 private 37 public :: & 38 hgh_projector_t, & 39 hgh_projector_null, & 40 hgh_projector_init, & 41 dhgh_project, zhgh_project, & 42 dhgh_project_bra, & 43 zhgh_project_bra, & 44 dhgh_project_ket, & 45 zhgh_project_ket, & 46 hgh_projector_end 47 48 49 50 type hgh_projector_t 51 private 52 integer :: n_s !< number of points inside the sphere 53 FLOAT, pointer, public :: dp(:, :) !< projectors 54 CMPLX, pointer, public :: zp(:, :) 55 FLOAT, public :: h(3, 3) !< parameters 56 FLOAT, public :: k(3, 3) !< spin-orbit parameters 57 end type hgh_projector_t 58 59 60contains 61 62 ! --------------------------------------------------------- 63 subroutine hgh_projector_null(hgh_p) 64 type(hgh_projector_t), intent(out) :: hgh_p 65 66 PUSH_SUB(hgh_projector_null) 67 68 nullify(hgh_p%dp) 69 nullify(hgh_p%zp) 70 hgh_p%h = M_ZERO 71 hgh_p%k = M_ZERO 72 73 POP_SUB(hgh_projector_null) 74 end subroutine hgh_projector_null 75 76 ! --------------------------------------------------------- 77 subroutine hgh_projector_init(hgh_p, sm, reltyp, a, l, lm, so_strength) 78 type(hgh_projector_t), intent(inout) :: hgh_p 79 type(submesh_t), intent(in) :: sm 80 integer, intent(in) :: reltyp 81 type(atom_t), target, intent(in) :: a 82 integer, intent(in) :: l, lm 83 FLOAT, intent(in) :: so_strength 84 85 integer :: is, i 86 FLOAT :: v, x(1:3) 87 type(ps_t), pointer :: ps 88 89 PUSH_SUB(hgh_projector_init) 90 91 hgh_p%n_s = sm%np 92 if(reltyp == 0) then 93 SAFE_ALLOCATE(hgh_p%dp(1:hgh_p%n_s, 1:3)) 94 x = M_ZERO 95 do is = 1, hgh_p%n_s 96 x(1:3) = sm%x(is, 1:3) 97 98 do i = 1, 3 99 call species_real_nl_projector(a%species, x, l, lm, i, v) 100 hgh_p%dp (is, i) = v 101 end do 102 end do 103 else 104 SAFE_ALLOCATE(hgh_p%zp(1:hgh_p%n_s, 1:3)) 105 do i = 1, 3 106 call species_nl_projector(a%species, hgh_p%n_s, sm%x(:, 0:3), l, lm, i, hgh_p%zp(:, i)) 107 end do 108 end if 109 110 ps => species_ps(a%species) 111 hgh_p%h(:, :) = ps%h(l, :, :) 112 hgh_p%k(:, :) = ps%k(l, :, :)*so_strength 113 nullify(ps) 114 115 POP_SUB(hgh_projector_init) 116 end subroutine hgh_projector_init 117 118 ! --------------------------------------------------------- 119 subroutine hgh_projector_end(hgh_p) 120 type(hgh_projector_t), intent(inout) :: hgh_p 121 122 PUSH_SUB(hgh_projector_end) 123 124 SAFE_DEALLOCATE_P(hgh_p%dp) 125 SAFE_DEALLOCATE_P(hgh_p%zp) 126 127 POP_SUB(hgh_projector_end) 128 end subroutine hgh_projector_end 129 130 131#include "undef.F90" 132#include "real.F90" 133#include "hgh_projector_inc.F90" 134 135#include "undef.F90" 136#include "complex.F90" 137#include "hgh_projector_inc.F90" 138 139end module hgh_projector_oct_m 140 141!! Local Variables: 142!! mode: f90 143!! coding: utf-8 144!! End: 145