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 states_elec_group_oct_m 22 use batch_oct_m 23 use global_oct_m 24 use loct_pointer_oct_m 25 use messages_oct_m 26 use profiling_oct_m 27 use states_elec_dim_oct_m 28 use wfs_elec_oct_m 29 30 implicit none 31 32 private 33 34 public :: & 35 states_elec_group_t, & 36 states_elec_group_null, & 37 states_elec_group_end, & 38 states_elec_group_copy 39 40 type states_elec_group_t 41 ! Components are public by default 42 type(wfs_elec_t), pointer :: psib(:, :) !< A set of wave-functions blocks 43 integer :: nblocks !< The number of blocks 44 integer :: block_start !< The lowest index of local blocks 45 integer :: block_end !< The highest index of local blocks 46 integer, pointer :: iblock(:, :) !< A map, that for each state index, returns the index of block containing it 47 integer, pointer :: block_range(:, :) !< Each block contains states from block_range(:, 1) to block_range(:, 2) 48 integer, pointer :: block_size(:) !< The number of states in each block. 49 logical, pointer :: block_is_local(:, :) !< It is true if the block is in this node. 50 integer, allocatable :: block_node(:) !< The node that contains each block 51 integer, allocatable :: rma_win(:, :) !< The MPI window for one side communication 52 logical :: block_initialized !< For keeping track of the blocks to avoid memory leaks 53 end type states_elec_group_t 54 55contains 56 57 ! --------------------------------------------------------- 58 subroutine states_elec_group_null(this) 59 type(states_elec_group_t), intent(out) :: this 60 61 PUSH_SUB(states_elec_group_null) 62 63 nullify(this%psib) 64 nullify(this%iblock) 65 nullify(this%block_range) 66 nullify(this%block_size) 67 nullify(this%block_is_local) 68 69 this%block_initialized = .false. 70 71 POP_SUB(states_elec_group_null) 72 end subroutine states_elec_group_null 73 74 ! --------------------------------------------------------- 75 subroutine states_elec_group_end(this, d) 76 type(states_elec_group_t), intent(inout) :: this 77 type(states_elec_dim_t), intent(in) :: d 78 79 integer :: ib, iq 80 81 PUSH_SUB(states_elec_group_end) 82 83 if (this%block_initialized) then 84 do ib = 1, this%nblocks 85 do iq = d%kpt%start, d%kpt%end 86 if (this%block_is_local(ib, iq)) then 87 call this%psib(ib, iq)%end() 88 end if 89 end do 90 end do 91 92 SAFE_DEALLOCATE_P(this%psib) 93 94 SAFE_DEALLOCATE_P(this%iblock) 95 SAFE_DEALLOCATE_P(this%block_range) 96 SAFE_DEALLOCATE_P(this%block_size) 97 SAFE_DEALLOCATE_P(this%block_is_local) 98 SAFE_DEALLOCATE_A(this%block_node) 99 this%block_initialized = .false. 100 end if 101 102 POP_SUB(states_elec_group_end) 103 end subroutine states_elec_group_end 104 105 !--------------------------------------------------------- 106 107 subroutine states_elec_group_copy(d, group_in, group_out, copy_data) 108 type(states_elec_dim_t), intent(in) :: d 109 type(states_elec_group_t), intent(in) :: group_in 110 type(states_elec_group_t), intent(out) :: group_out 111 logical, optional, intent(in) :: copy_data 112 113 integer :: qn_start, qn_end, ib, iqn 114 115 PUSH_SUB(states_elec_group_copy) 116 117 call states_elec_group_null(group_out) 118 119 120 group_out%nblocks = group_in%nblocks 121 group_out%block_start = group_in%block_start 122 group_out%block_end = group_in%block_end 123 group_out%block_initialized = group_in%block_initialized 124 125 if(group_out%block_initialized) then 126 127 ASSERT(associated(group_in%psib)) 128 129 qn_start = d%kpt%start 130 qn_end = d%kpt%end 131 132 SAFE_ALLOCATE(group_out%psib(1:group_out%nblocks, qn_start:qn_end)) 133 134 do iqn = qn_start, qn_end 135 do ib = group_out%block_start, group_out%block_end 136 call group_in%psib(ib, iqn)%copy_to(group_out%psib(ib, iqn), copy_data = optional_default(copy_data, .true.)) 137 end do 138 end do 139 140 call loct_pointer_copy(group_out%iblock, group_in%iblock) 141 call loct_pointer_copy(group_out%block_range, group_in%block_range) 142 call loct_pointer_copy(group_out%block_size, group_in%block_size) 143 call loct_pointer_copy(group_out%block_is_local, group_in%block_is_local) 144 call loct_allocatable_copy(group_out%block_node, group_in%block_node) 145 call loct_allocatable_copy(group_out%rma_win, group_in%rma_win) 146 147 end if 148 149 POP_SUB(states_elec_group_copy) 150 end subroutine states_elec_group_copy 151 152end module states_elec_group_oct_m 153 154 155!! Local Variables: 156!! mode: f90 157!! coding: utf-8 158!! End: 159