1! 2 3 module class_Dist 4#ifdef MPI 5 use mpi 6#endif 7 8 use class_BlockCyclicDist 9 use class_PEXSIDist 10 11 implicit none 12 13 character(len=*), parameter :: mod_name=__FILE__ 14 15!--------------------------------------------- 16 type, private :: distKind 17 private 18 integer :: i 19 end type distKind 20 21 type(distKind), parameter, public :: TYPE_NULL = distKind(-1) 22 type(distKind), parameter, public :: TYPE_BLOCK_CYCLIC = distKind(1) 23 type(distKind), parameter, public :: TYPE_PEXSI = distKind(2) 24 25 public operator(==) 26 interface operator(==) 27 module procedure isequal_ 28 end interface 29!--------------------------------------------- 30 31 type Dist_ 32 integer :: refCount = 0 33 character(len=36) :: id 34 !------------------------ 35 character(len=256) :: name = "null Dist" 36 !------------------------ 37 integer :: group = -1 ! MPI group 38 integer :: node = -1 ! MPI rank in group (my_proc) 39 integer :: nodes = 0 ! MPI size of group (nprocs) 40 integer :: node_io = -1 ! Node capable of IO 41 !------------------------ 42 type(distKind) :: dist_type 43 type(BlockCyclicDist) :: bdist 44 type(PEXSIDist) :: pdist 45 !------------------------------------------ 46 end type Dist_ 47 48 type Dist 49 type(Dist_), pointer :: data => null() 50 end type Dist 51 52 public :: newDistribution, distType, group 53 public :: num_local_elements, node_handling_element 54 public :: index_local_to_global, index_global_to_local 55 56 interface newDistribution 57 module procedure newDistribution_ 58 end interface 59 interface distType 60 module procedure dist_type 61 end interface 62 interface group 63 module procedure group_ 64 end interface 65 interface num_local_elements 66 module procedure num_local_elements_ 67 end interface 68 interface index_local_to_global 69 module procedure index_local_to_global_ 70 end interface 71 interface index_global_to_local 72 module procedure index_global_to_local_ 73 end interface 74 interface node_handling_element 75 module procedure node_handling_element_ 76 end interface 77 78!==================================== 79#define TYPE_NAME Dist 80#include "basic_type.inc" 81!==================================== 82 83 subroutine delete_Data(spdata) 84 type(Dist_) :: spdata 85 if (spdata%dist_type == TYPE_BLOCK_CYCLIC) then 86 call delete(spdata%bdist) 87 else if (spdata%dist_type == TYPE_PEXSI) then 88 call delete(spdata%pdist) 89 else if (spdata%dist_type == TYPE_NULL) then 90 ! do nothing 91 else 92 ! do nothing 93 end if 94 95 end subroutine delete_Data 96 97 function dist_type(this) result(dtype) 98 type(Dist) :: this 99 type(distKind) :: dtype 100 101 dtype = this%data%dist_type 102 103 end function dist_type 104 105 function isequal_ (a,b) result(iseq) 106 type(distKind), intent(in) :: a, b 107 logical :: iseq 108 109 iseq = (a%i == b%i) 110 end function isequal_ 111! 112 subroutine newDistribution_(Blocksize,Group,this,dist_type,name) 113 !........................................ 114 ! Constructor 115 !........................................ 116 type (Dist), intent(inout) :: this 117 integer, intent(in) :: Blocksize 118 integer, intent(in) :: Group 119 type (distKind), intent(in) :: dist_type 120 character(len=*), intent(in), optional :: name 121 122 integer :: error 123 124 call init(this) 125 126 if (dist_type == TYPE_BLOCK_CYCLIC) then 127 call newDistribution(Blocksize,Group,this%data%bdist,name) 128 else if (dist_type == TYPE_PEXSI) then 129 call newDistribution(Blocksize,Group,this%data%pdist,name) 130 else if (dist_type == TYPE_NULL) then 131 ! do nothing 132 else 133 ! do nothing 134 end if 135 136 this%data%dist_type = dist_type 137 138 if (present(name)) then 139 this%data%name = trim(name) 140 else 141 this%data%name = "(Abstract Distribution from BS and Group)" 142 endif 143 call tag_new_object(this) 144 145 end subroutine newDistribution_ 146 147!----------------------------------------------------------- 148 function num_local_elements_(this,nels,Node) result(nl) 149 type(Dist), intent(in) :: this 150 integer, intent(in) :: nels 151 integer, intent(in) :: Node 152 integer :: nl 153 154 if (dist_type(this) == TYPE_BLOCK_CYCLIC) then 155 nl = num_local_elements(this%data%bdist,nels,Node) 156 else if (dist_type(this) == TYPE_PEXSI) then 157 nl = num_local_elements(this%data%pdist,nels,Node) 158 else if (dist_type(this) == TYPE_NULL) then 159 nl = 0 160 else 161 nl = 0 162 end if 163 164 end function num_local_elements_ 165 166 function index_local_to_global_(this,il,Node) result(ig) 167 type(Dist), intent(in) :: this 168 integer, intent(in) :: il 169 integer, intent(in) :: Node 170 integer :: ig 171 172 if (dist_type(this) == TYPE_BLOCK_CYCLIC) then 173 ig = index_local_to_global(this%data%bdist,il,Node) 174 else if (dist_type(this) == TYPE_PEXSI) then 175 ig = index_local_to_global(this%data%pdist,il,Node) 176 else if (dist_type(this) == TYPE_NULL) then 177 ig = 0 178 else 179 ig = 0 180 end if 181 182 end function index_local_to_global_ 183 184 function index_global_to_local_(this,ig,Node) result(il) 185 type(Dist), intent(in) :: this 186 integer, intent(in) :: ig 187 integer, intent(in), optional :: Node 188 integer :: il 189 190 if (dist_type(this) == TYPE_BLOCK_CYCLIC) then 191 il = index_global_to_local(this%data%bdist,ig,Node) 192 else if (dist_type(this) == TYPE_PEXSI) then 193 il = index_global_to_local(this%data%pdist,ig,Node) 194 else if (dist_type(this) == TYPE_NULL) then 195 il = 0 196 else 197 il = 0 198 end if 199 200 end function index_global_to_local_ 201 202 function node_handling_element_(this,ig) result(proc) 203 type(Dist), intent(in) :: this 204 integer, intent(in) :: ig 205 integer :: proc 206 207 if (dist_type(this) == TYPE_BLOCK_CYCLIC) then 208 proc = node_handling_element(this%data%bdist,ig) 209 else if (dist_type(this) == TYPE_PEXSI) then 210 proc = node_handling_element(this%data%pdist,ig) 211 else if (dist_type(this) == TYPE_NULL) then 212 proc = MPI_UNDEFINED 213 else 214 proc = MPI_UNDEFINED 215 end if 216 217 end function node_handling_element_ 218 219!----------------------------------------------------------- 220 function group_(this) result(g) 221 type(Dist), intent(in) :: this 222 integer :: g 223 224 if (dist_type(this) == TYPE_BLOCK_CYCLIC) then 225 g = this%data%bdist%data%group 226 else if (dist_type(this) == TYPE_PEXSI) then 227 g = this%data%pdist%data%group 228 else if (dist_type(this) == TYPE_NULL) then 229 g = MPI_GROUP_NULL 230 else 231 g = MPI_GROUP_NULL 232 end if 233 234 end function group_ 235 236 237end module class_Dist 238