communicator_base.F90 coverage: 0.00 %func 0.00 %block
1) module Communicator_Base_module
2)
3) use PFLOTRAN_Constants_module
4)
5) implicit none
6)
7) private
8)
9) #include "petsc/finclude/petscvec.h"
10) #include "petsc/finclude/petscvec.h90"
11)
12) type, abstract, public :: communicator_type
13) contains
14) procedure(SetDM), public, deferred :: SetDM
15) procedure(VecToVec), public, deferred :: GlobalToLocal
16) procedure(VecToVec), public, deferred :: LocalToGlobal
17) procedure(VecToVec), public, deferred :: LocalToLocal
18) procedure(VecToVec), public, deferred :: GlobalToNatural
19) procedure(VecToVec), public, deferred :: NaturalToGlobal
20) procedure(MapArray), public, deferred :: AONaturalToPetsc
21) procedure(BaseDestroy), public, deferred :: Destroy
22) end type communicator_type
23)
24) abstract interface
25)
26) #ifdef SIMPLIFY
27) subroutine SetDM(this)
28) import communicator_type
29) implicit none
30) class(communicator_type) :: this
31) #else
32) subroutine SetDM(this,dm_ptr)
33) use DM_Kludge_module
34) import communicator_type
35) implicit none
36) class(communicator_type) :: this
37) type(dm_ptr_type) :: dm_ptr
38) #endif
39) end subroutine
40)
41) subroutine VecToVec(this,source,destination)
42) import communicator_type
43) implicit none
44) class(communicator_type) :: this
45) Vec :: source
46) Vec :: destination
47) end subroutine VecToVec
48)
49) subroutine MapArray(this,array)
50) import communicator_type
51) implicit none
52) class(communicator_type) :: this
53) PetscInt :: array(:)
54) end subroutine MapArray
55)
56) subroutine BaseDestroy(this)
57) import communicator_type
58) implicit none
59) class(communicator_type) :: this
60) end subroutine BaseDestroy
61)
62) end interface
63)
64) public :: CommCreateProcessorGroups
65)
66) contains
67)
68) ! ************************************************************************** !
69)
70) subroutine CommCreateProcessorGroups(option,num_groups)
71) !
72) ! Splits MPI_COMM_WORLD into N separate
73) ! processor groups
74) !
75) ! Author: Glenn Hammond
76) ! Date: 08/11/09
77) !
78)
79) use Option_module
80)
81) implicit none
82)
83) #include "petsc/finclude/petscsys.h"
84)
85) type(option_type) :: option
86) PetscInt :: num_groups
87)
88) PetscInt :: local_commsize
89) PetscInt :: offset, delta, remainder
90) PetscInt :: igroup
91) PetscMPIInt :: mycolor_mpi, mykey_mpi
92) PetscErrorCode :: ierr
93)
94) local_commsize = option%global_commsize / num_groups
95) remainder = option%global_commsize - num_groups * local_commsize
96) offset = 0
97) do igroup = 1, num_groups
98) delta = local_commsize
99) if (igroup < remainder) delta = delta + 1
100) if (option%global_rank >= offset .and. &
101) option%global_rank < offset + delta) exit
102) offset = offset + delta
103) enddo
104) mycolor_mpi = igroup
105) option%mygroup_id = igroup
106) mykey_mpi = option%global_rank - offset
107) call MPI_Comm_split(MPI_COMM_WORLD,mycolor_mpi,mykey_mpi,option%mycomm,ierr)
108) call MPI_Comm_group(option%mycomm,option%mygroup,ierr)
109)
110) PETSC_COMM_WORLD = option%mycomm
111) call PetscInitialize(PETSC_NULL_CHARACTER, ierr);CHKERRQ(ierr)
112) call MPI_Comm_rank(option%mycomm,option%myrank, ierr)
113) call MPI_Comm_size(option%mycomm,option%mycommsize,ierr)
114)
115) end subroutine CommCreateProcessorGroups
116)
117) end module Communicator_Base_module