1! --- 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt . 6! See Docs/Contributors.txt for a list of contributors. 7! --- 8 9! This should be used to specify the orbital distribution, in the spirit 10! of the Scalapack descriptors 11! Closer to Siesta, it should mimic the "order-N arrays" of J. Gale 12! 13! subroutine GetNodeOrbs( NOrb, Node, Nodes, NOrbNode) 14! calculates the number of orbitals stored on the (local)? Node 15! Could be in principle that of any node. 16! NOTE: NOrb is not always no_u: nocc in savepsi... (ksv) 17! NOrbNode = nOrbPerNode(Node+1) 18! subroutine GlobalToLocalOrb( GOrb, Node, Nodes, LOrb) 19! ditto, returning zero if Gorb is not handled in local node 20! LOrb = nG2L(GOrb) 21! subroutine LocalToGlobalOrb( LOrb, Node, Nodes, GOrb) 22! easy in principle... but Julian has the option to have the 23! info from any node... 24! GOrb = nL2G(LOrb,Node+1) 25! It is used in iodm_netcdf, but it could be emulated by having 26! each node send also the info about the global counterparts. 27! subroutine WhichNodeOrb( GOrb, Nodes, Node) 28! returns the Node number that handles GOrb 29! Node = nNode(GOrb) 30 31 module class_OrbitalDistribution 32 33 implicit none 34 35 character(len=*), parameter :: mod_name="class_OrbitalDistribution" 36 37 type OrbitalDistribution_ 38 integer :: refCount = 0 39 character(len=36) :: id 40 !------------------------ 41 character(len=256) :: name = "null OrbitalDistribution" 42 !------------------------ 43 integer :: comm = -1 ! MPI communicator 44 integer :: node = -1 ! MPI rank in comm (my_proc) 45 integer :: nodes = 0 ! MPI size of comm (nprocs) 46 integer :: node_io = -1 ! Node capable of IO 47 !------------------------ 48 integer :: blocksize = 0 ! Flag to determine whether the dist is block-cyclic 49 ! 50 integer :: isrcproc = 0 ! Processor holding the first element (unused) 51 !------------------------------------------ 52 ! 53 ! If the distribution is not block-cyclic, 54 ! we need explicit information 55 ! 56 ! 57 ! Just for sanity: number of global elements (see function num_local_elements) 58 ! This will have to be set at the time of filling in the arrays 59 ! for non-block-cyclic distributions 60 integer :: n = -1 61 ! Number of elements held by each proc 62 integer, pointer :: nroc_proc(:) => null() 63 ! Global index from Local index (in my_proc) 64 integer, pointer :: nl2g(:) => null() 65 ! Local index (in my_proc) from global index (zero if not held by my_proc) 66 integer, pointer :: ng2l(:) => null() 67 ! Processor number from global index (base 0) 68 integer, pointer :: ng2p(:) => null() 69 !------------------------------------------ 70 end type OrbitalDistribution_ 71 72 type OrbitalDistribution 73 type(OrbitalDistribution_), pointer :: data => null() 74 end type OrbitalDistribution 75 76 public :: newDistribution 77 public :: dist_comm, dist_node, dist_nodes 78 public :: num_local_elements, node_handling_element 79 public :: index_local_to_global, index_global_to_local 80 public :: print_type 81 82 interface newDistribution 83 module procedure newBlockCyclicDistribution 84 end interface 85 86 interface print_type 87 module procedure printOrbitalDistribution 88 end interface 89 90 interface dist_comm 91 module procedure comm_ 92 end interface 93 interface dist_node 94 module procedure dist_node_ 95 end interface 96 interface dist_nodes 97 module procedure dist_nodes_ 98 end interface 99 100!==================================== 101#define TYPE_NAME OrbitalDistribution 102#include "basic_type.inc" 103!==================================== 104 105 subroutine delete_Data(spdata) 106 type(OrbitalDistribution_) :: spdata 107 if ( associated(spdata%nroc_proc) ) & 108 deallocate( spdata%nroc_proc) 109 if ( associated(spdata%nl2g) ) & 110 deallocate( spdata%nl2g) 111 if ( associated(spdata%ng2l) ) & 112 deallocate( spdata%ng2l) 113 if ( associated(spdata%ng2p) ) & 114 deallocate( spdata%ng2p) 115 end subroutine delete_Data 116 117 subroutine newBlockCyclicDistribution(Blocksize,Comm,this,name) 118 !........................................ 119 ! Constructor 120 !........................................ 121 type (OrbitalDistribution), intent(inout) :: this 122 integer, intent(in) :: Blocksize 123 integer, intent(in) :: Comm 124 character(len=*), intent(in), optional :: name 125 126 integer :: error 127 128 call init(this) 129 130 this%data%blocksize = Blocksize 131 this%data%comm = Comm 132 133#ifdef MPI 134 call MPI_Comm_Rank( Comm, this%data%node, error ) 135 call MPI_Comm_Size( Comm, this%data%nodes, error ) 136#else 137 this%data%node = 0 138 this%data%nodes = 1 139#endif 140 this%data%node_io = 0 141 142 if (present(name)) then 143 this%data%name = trim(name) 144 else 145 this%data%name = "(Distribution from BlockSize and Comm)" 146 endif 147 call tag_new_object(this) 148 149 end subroutine newBlockCyclicDistribution 150 151!----------------------------------------------------------- 152 function num_local_elements(this,nels,Node) result(nl) 153#ifdef MPI 154 use mpi_siesta, only: MPI_COMM_Self 155#endif 156 type(OrbitalDistribution), intent(in) :: this 157 integer, intent(in) :: nels 158 integer, intent(in), optional :: Node 159 integer :: nl 160 161 integer :: lNode 162 integer :: Remainder, MinPerNode, RemainderBlocks 163 164 if ( present(Node) ) then 165 lNode = Node 166 else 167 lNode = this%data%Node 168 end if 169 170 if (this%data%blocksize == 0) then 171 ! Not a block-cyclic distribution 172 if (nels /= this%data%n) then 173 call die("Cannot figure out no_l if nels/=n") 174 nl = -1 175 endif 176 if (.not. associated(this%data%nroc_proc)) then 177 call die("Dist arrays not setup") 178 nl = -1 179 endif 180 nl = this%data%nroc_proc(lNode) 181#ifdef MPI 182 else if ( this%data%Comm == MPI_Comm_Self ) then ! globalized local distribution 183 184 ! No matter the node, it is the same orbital 185 ! The interface for creating such an orbital distribution 186 ! is 187 ! call newDistribution(nrows_g(sp),MPI_Comm_Self,dit) 188 ! TODO Talk with Alberto about adding: "lproc = this%data%nodes - MOD(proc,this%data%nodes)" 189 if ( nels /= this%data%blocksize ) & 190 call die('Contact Nick Papior Andersen, nickpapior@gmail.com nel') 191 nl = this%data%blocksize 192#endif 193 else ! block-cyclic distribution 194 195 ! Use Julian's code for now, although this is 196 ! really Scalapack's numroc (copied at the end) 197 ! nl = numroc(nels,this%data%blocksize,Node,this%data%isrcproc,this%data%nodes) 198 199 MinPerNode = nels / (this%data%nodes*this%data%blocksize) 200 Remainder = nels - MinPerNode * this%data%nodes*this%data%blocksize 201 RemainderBlocks = Remainder / this%data%blocksize 202 Remainder = Remainder - RemainderBlocks * this%data%blocksize 203 nl = MinPerNode*this%data%blocksize 204 if (lNode.lt.RemainderBlocks) nl = nl + this%data%blocksize 205 if (lNode.eq.RemainderBlocks) nl = nl + Remainder 206 207 endif 208 end function num_local_elements 209 210 function dist_node_(this) result(node) 211 type(OrbitalDistribution), intent(in) :: this 212 integer :: node 213 node = this%data%node 214 end function dist_node_ 215 function dist_nodes_(this) result(nodes) 216 type(OrbitalDistribution), intent(in) :: this 217 integer :: nodes 218 nodes = this%data%nodes 219 end function dist_nodes_ 220 221 function index_local_to_global(this,il,Node) result(ig) 222#ifdef MPI 223 use mpi_siesta, only: MPI_COMM_Self 224#endif 225 type(OrbitalDistribution), intent(in) :: this 226 integer, intent(in) :: il 227 integer, intent(in), optional :: Node 228 integer :: ig 229 230 integer :: lNode 231 integer :: LBlock, LEle 232 233 lNode = this%data%node 234 if ( present(Node) ) lNode = Node 235 236 if (this%data%blocksize == 0) then 237 ! Not a block-cyclic distribution 238 if (lNode /= this%data%node) then 239 call die("Cannot figure out ig if Node/=my_proc") 240 ig = -1 241 endif 242 if (.not. associated(this%data%nl2g)) then 243 call die("Dist arrays not setup") 244 ig = -1 245 endif 246 ig = this%data%nl2g(il) 247#ifdef MPI 248 else if ( this%data%Comm == MPI_Comm_Self ) then ! globalized local distribution 249 250 ! No matter the node, it is the same orbital 251 ! The interface for creating such an orbital distribution 252 ! is 253 ! call newDistribution(nrows_g(sp),MPI_Comm_Self,dit) 254 if ( il > this%data%blocksize ) & 255 call die('Contact Nick Papior Andersen, nickpapior@gmail.com l2g') 256 257 ig = il 258#endif 259 else ! block-cyclic distribution 260 261 ! Use Julian's code for now, although this is 262 ! really Scalapack's indxl2g (copied at the end) 263 ! ig = indxl2g(il,this%data%blocksize,Node,this%data%isrcproc,this%data%nodes) 264 265 ! Find local block number 266 LBlock = ((il -1)/this%data%blocksize) 267 ! Substract local base line to find element number within the block 268 LEle = il - LBlock*this%data%blocksize 269 ! Calculate global index 270 ig = (LBlock*this%data%nodes + lNode)*this%data%blocksize + LEle 271 272 endif 273 end function index_local_to_global 274 275 276 function index_global_to_local(this,ig,Node) result(il) 277#ifdef MPI 278 use mpi_siesta, only: MPI_COMM_Self 279#endif 280 type(OrbitalDistribution), intent(in) :: this 281 integer, intent(in) :: ig 282 ! In case Node is not supplied we expect it to request its 283 ! own index 284 integer, intent(in), optional :: Node 285 integer :: il 286 287 integer :: lNode 288 integer :: LBlock, LEle, GBlock, OrbCheck 289 290 lNode = this%data%node 291 if ( present(Node) ) lNode = Node 292 293 if (this%data%blocksize == 0) then 294 ! Not a block-cyclic distribution 295 if (lNode /= this%data%node) then 296 call die("Cannot figure out il if Node/=my_proc") 297 il = -1 298 endif 299 if (.not. associated(this%data%ng2l)) then 300 call die("Dist arrays not setup") 301 il = -1 302 endif 303 il = this%data%ng2l(ig) 304 ! Alternatively: Use an extended ng2p, in which "local" elements are negative. 305 ! if (ng2p(ig)) < 0 then il = -that, else 0 306 ! Example: 0 0 1 1 2 2 -1 -2 4 4 5 5 0 0 1 1 2 2 -3 -4 5 5 ..... (nb=2, np=6) 307#ifdef MPI 308 else if ( this%data%Comm == MPI_Comm_Self ) then ! globalized local distribution 309 310 ! No matter the node, it is the same orbital 311 ! The interface for creating such an orbital distribution 312 ! is 313 ! call newDistribution(nrows_g(sp),MPI_Comm_Self,dit) 314 if ( ig > this%data%blocksize ) & 315 call die('Contact Nick Papior Andersen, nickpapior@gmail.com g2l') 316 317 il = ig 318#endif 319 else ! block-cyclic distribution 320 321 ! Use Julian's code for now, although this is 322 ! really a combination of Scalapack's indxg2p and indxg2l (copied at the end) 323 ! owner = indxg2p(ig,this%data%blocksize,Node,this%data%isrcproc,this%data%nodes) 324 ! if (owner == Node) then 325 ! il = indxg2l(ig,this%data%blocksize,Node,this%data%isrcproc,this%data%nodes) 326 ! else 327 ! il = 0 328 ! endif 329 330 ! Find global block number 331 GBlock = ((ig -1)/this%data%blocksize) 332 ! Substract global base line to find element number within the block 333 LEle = ig - GBlock*this%data%blocksize 334 ! Find the block number on the local node 335 LBlock = ((GBlock - lNode)/this%data%nodes) 336 ! Generate the local orbital pointer based on the local block number 337 il = LEle + LBlock*this%data%blocksize 338 ! Check that this is consistent - if it is not then this 339 ! local orbital is not on this node and so we return 0 340 ! to indicate this. 341 OrbCheck = (LBlock*this%data%nodes + lNode)*this%data%blocksize + LEle 342 if (OrbCheck.ne.ig) il = 0 343 endif 344 345 end function index_global_to_local 346 347 function node_handling_element(this,ig) result(proc) 348#ifdef MPI 349 use mpi_siesta, only: MPI_COMM_Self 350#endif 351 type(OrbitalDistribution), intent(in) :: this 352 integer, intent(in) :: ig 353 integer :: proc 354 355 integer :: GBlock 356 357 if (this%data%blocksize == 0) then 358 ! Not a block-cyclic distribution 359 if (.not. associated(this%data%ng2p)) then 360 call die("Dist arrays not setup") 361 proc = -1 362 endif 363 proc = this%data%ng2p(ig) 364 ! Alternatively: Use an extended ng2p, in which "local" elements are negative. 365 ! if (ng2p(ig)) < 0 then il = -that, else 0 366 ! Example: 0 0 1 1 2 2 -1 -2 4 4 5 5 0 0 1 1 2 2 -3 -4 5 5 ..... (nb=2, np=6) 367#ifdef MPI 368 else if ( this%data%Comm == MPI_Comm_Self ) then ! globalized local distribution 369 370 ! No matter the node, it is the same orbital 371 ! The interface for creating such an orbital distribution 372 ! is 373 ! call newDistribution(nrows_g(sp),MPI_Comm_Self,dit) 374 if ( ig > this%data%blocksize ) & 375 call die('Contact Nick Papior Andersen, nickpapior@gmail.com nhe') 376 377 proc = this%data%node 378#endif 379 else ! block-cyclic distribution 380 381 ! Use Julian's code for now, although this is 382 ! really Scalapack's indxg2p (copied at the end) 383 ! proc = indxg2p(ig,this%data%blocksize,dummy_Node, & 384 ! this%data%isrcproc,this%data%nodes) 385 386 ! Find global block number 387 GBlock = ((ig -1)/this%data%blocksize) 388 ! Find the Node number that has this block 389 proc = mod(GBlock,this%data%nodes) 390 391 endif 392 393 end function node_handling_element 394 395 ! Returns the communicator assigned to this 396 ! distribution 397 function comm_(this) result(comm) 398 type(OrbitalDistribution), intent(in) :: this 399 integer :: comm 400 comm = this%data%comm 401 end function comm_ 402 403 subroutine printOrbitalDistribution(this) 404 type(OrbitalDistribution), intent(in) :: this 405 406 if (.not. initialized(this) ) then 407 print "(a)", "Orbital distribution Not Associated" 408 RETURN 409 endif 410 411 print "(a,/,t4,5(a,i0),a)", & 412 " <orb-dist:"//trim(this%data%name), & 413 " comm=", this%data%comm, & 414 " node/nodes=", this%data%node,' / ',this%data%nodes, & 415 " blocksize=",this%data%blocksize, & 416 ", refcount: ", refcount(this), ">" 417 end subroutine printOrbitalDistribution 418 419!======================================================================== 420 INTEGER FUNCTION NUMROC( N, NB, IPROC, ISRCPROC, NPROCS ) 421! 422! -- ScaLAPACK tools routine (version 1.7) -- 423! University of Tennessee, Knoxville, Oak Ridge National Laboratory, 424! and University of California, Berkeley. 425! May 1, 1997 426! 427! .. Scalar Arguments .. 428 INTEGER IPROC, ISRCPROC, N, NB, NPROCS 429! .. 430! 431! Purpose 432! ======= 433! 434! NUMROC computes the NUMber of Rows Or Columns of a distributed 435! matrix owned by the process indicated by IPROC. 436! 437! Arguments 438! ========= 439! 440! N (global input) INTEGER 441! The number of rows/columns in distributed matrix. 442! 443! NB (global input) INTEGER 444! Block size, size of the blocks the distributed matrix is 445! split into. 446! 447! IPROC (local input) INTEGER 448! The coordinate of the process whose local array row or 449! column is to be determined. 450! 451! ISRCPROC (global input) INTEGER 452! The coordinate of the process that possesses the first 453! row or column of the distributed matrix. 454! 455! NPROCS (global input) INTEGER 456! The total number processes over which the matrix is 457! distributed. 458! 459! ===================================================================== 460! 461! .. Local Scalars .. 462 INTEGER EXTRABLKS, MYDIST, NBLOCKS 463! .. 464! .. Intrinsic Functions .. 465 INTRINSIC MOD 466! .. 467! .. Executable Statements .. 468! 469! Figure PROC's distance from source process 470! 471 MYDIST = MOD( NPROCS+IPROC-ISRCPROC, NPROCS ) 472! 473! Figure the total number of whole NB blocks N is split up into 474! 475 NBLOCKS = N / NB 476! 477! Figure the minimum number of rows/cols a process can have 478! 479 NUMROC = (NBLOCKS/NPROCS) * NB 480! 481! See if there are any extra blocks 482! 483 EXTRABLKS = MOD( NBLOCKS, NPROCS ) 484! 485! If I have an extra block 486! 487 IF( MYDIST.LT.EXTRABLKS ) THEN 488 NUMROC = NUMROC + NB 489! 490! If I have last block, it may be a partial block 491! 492 ELSE IF( MYDIST.EQ.EXTRABLKS ) THEN 493 NUMROC = NUMROC + MOD( N, NB ) 494 END IF 495! 496 RETURN 497! 498! End of NUMROC 499! 500 END function numroc 501!========================================================================== 502 INTEGER FUNCTION INDXL2G( INDXLOC, NB, IPROC, ISRCPROC, NPROCS ) 503! 504! -- ScaLAPACK tools routine (version 1.7) -- 505! University of Tennessee, Knoxville, Oak Ridge National Laboratory, 506! and University of California, Berkeley. 507! May 1, 1997 508! 509! .. Scalar Arguments .. 510 INTEGER INDXLOC, IPROC, ISRCPROC, NB, NPROCS 511! .. 512! 513! Purpose 514! ======= 515! 516! INDXL2G computes the global index of a distributed matrix entry 517! pointed to by the local index INDXLOC of the process indicated by 518! IPROC. 519! 520! Arguments 521! ========= 522! 523! INDXLOC (global input) INTEGER 524! The local index of the distributed matrix entry. 525! 526! NB (global input) INTEGER 527! Block size, size of the blocks the distributed matrix is 528! split into. 529! 530! IPROC (local input) INTEGER 531! The coordinate of the process whose local array row or 532! column is to be determined. 533! 534! ISRCPROC (global input) INTEGER 535! The coordinate of the process that possesses the first 536! row/column of the distributed matrix. 537! 538! NPROCS (global input) INTEGER 539! The total number processes over which the distributed 540! matrix is distributed. 541! 542! ===================================================================== 543! 544! .. Intrinsic Functions .. 545 INTRINSIC MOD 546! .. 547! .. Executable Statements .. 548! 549 INDXL2G = NPROCS*NB*((INDXLOC-1)/NB) + MOD(INDXLOC-1,NB) + & 550 MOD(NPROCS+IPROC-ISRCPROC, NPROCS)*NB + 1 551! 552 RETURN 553! 554! End of INDXL2G 555! 556 END function indxl2g 557 558!=========================================================== 559 INTEGER FUNCTION INDXG2P( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS ) 560! 561! -- ScaLAPACK tools routine (version 1.7) -- 562! University of Tennessee, Knoxville, Oak Ridge National Laboratory, 563! and University of California, Berkeley. 564! May 1, 1997 565! 566! .. Scalar Arguments .. 567 INTEGER INDXGLOB, IPROC, ISRCPROC, NB, NPROCS 568! .. 569! 570! Purpose 571! ======= 572! 573! INDXG2P computes the process coordinate which posseses the entry of a 574! distributed matrix specified by a global index INDXGLOB. 575! 576! Arguments 577! ========= 578! 579! INDXGLOB (global input) INTEGER 580! The global index of the element. 581! 582! NB (global input) INTEGER 583! Block size, size of the blocks the distributed matrix is 584! split into. 585! 586! IPROC (local dummy) INTEGER 587! Dummy argument in this case in order to unify the calling 588! sequence of the tool-routines. 589! 590! ISRCPROC (global input) INTEGER 591! The coordinate of the process that possesses the first 592! row/column of the distributed matrix. 593! 594! NPROCS (global input) INTEGER 595! The total number processes over which the matrix is 596! distributed. 597! 598! ===================================================================== 599! 600! .. Intrinsic Functions .. 601 INTRINSIC MOD 602! .. 603! .. Executable Statements .. 604! 605 INDXG2P = MOD( ISRCPROC + (INDXGLOB - 1) / NB, NPROCS ) 606! 607 RETURN 608! 609! End of INDXG2P 610! 611 END function indxg2p 612 613!================================================================ 614 615 INTEGER FUNCTION INDXG2L( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS ) 616! 617! -- ScaLAPACK tools routine (version 1.7) -- 618! University of Tennessee, Knoxville, Oak Ridge National Laboratory, 619! and University of California, Berkeley. 620! May 1, 1997 621! 622! .. Scalar Arguments .. 623 INTEGER INDXGLOB, IPROC, ISRCPROC, NB, NPROCS 624! .. 625! 626! Purpose 627! ======= 628! 629! INDXG2L computes the local index of a distributed matrix entry 630! pointed to by the global index INDXGLOB. 631! 632! Arguments 633! ========= 634! 635! INDXGLOB (global input) INTEGER 636! The global index of the distributed matrix entry. 637! 638! NB (global input) INTEGER 639! Block size, size of the blocks the distributed matrix is 640! split into. 641! 642! IPROC (local dummy) INTEGER 643! Dummy argument in this case in order to unify the calling 644! sequence of the tool-routines. 645! 646! ISRCPROC (local dummy) INTEGER 647! Dummy argument in this case in order to unify the calling 648! sequence of the tool-routines. 649! 650! NPROCS (global input) INTEGER 651! The total number processes over which the distributed 652! matrix is distributed. 653! 654! ===================================================================== 655! 656! .. Intrinsic Functions .. 657 INTRINSIC MOD 658! .. 659! .. Executable Statements .. 660! 661 INDXG2L = NB*((INDXGLOB-1)/(NB*NPROCS))+MOD(INDXGLOB-1,NB)+1 662! 663 RETURN 664! 665! End of INDXG2L 666! 667 END function indxg2l 668 669end module class_OrbitalDistribution 670