1 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
---|
2 | ! Math and Computer Science Division, Argonne National Laboratory ! |
---|
3 | !----------------------------------------------------------------------- |
---|
4 | ! CVS mutils.F90,v 1.8 2005-11-18 23:15:38 rloy Exp |
---|
5 | ! CVS MCT_2_8_0 |
---|
6 | !BOP ------------------------------------------------------------------- |
---|
7 | ! |
---|
8 | ! !MODULE: mutils -- utilities for the sequential climate example |
---|
9 | ! |
---|
10 | ! !DESCRIPTION: |
---|
11 | ! |
---|
12 | ! !INTERFACE: |
---|
13 | ! |
---|
14 | module mutils |
---|
15 | |
---|
16 | ! module of utilties for the sequential climate example |
---|
17 | ! |
---|
18 | |
---|
19 | implicit none |
---|
20 | |
---|
21 | private |
---|
22 | ! except |
---|
23 | |
---|
24 | ! !PUBLIC TYPES: |
---|
25 | |
---|
26 | public get_index |
---|
27 | |
---|
28 | contains |
---|
29 | !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
---|
30 | ! Math and Computer Science Division, Argonne National Laboratory ! |
---|
31 | !BOP ------------------------------------------------------------------- |
---|
32 | ! |
---|
33 | ! !IROUTINE: get_index - get local index array and size |
---|
34 | ! for 3 standard decompositions of a grid. |
---|
35 | ! |
---|
36 | ! !DESCRIPTION: |
---|
37 | ! The routine get_index will return a local index array and size that can |
---|
38 | ! be passed to a GSMap_init routine for three possible decompositions: |
---|
39 | ! R - by row or latitude |
---|
40 | ! C - by column or longitude |
---|
41 | ! RC - row and column or checkerboard |
---|
42 | ! choice is determined by the value of ldecomp. |
---|
43 | ! |
---|
44 | ! !INTERFACE: |
---|
45 | |
---|
46 | subroutine get_index(ldecomp,nprocs,myproc,gnx,gny,gridbuf) |
---|
47 | ! !INPUT PARAMETERS: |
---|
48 | ! |
---|
49 | character(len=*),intent(inout) :: ldecomp ! decomp choice |
---|
50 | integer,intent(in) :: nprocs ! total number of MPI processes |
---|
51 | integer,intent(in) :: myproc ! my rank in local communicator |
---|
52 | integer,intent(in) :: gnx ! total points in X direction |
---|
53 | integer,intent(in) :: gny ! total points in Y direction |
---|
54 | |
---|
55 | ! !OUTPUT PARAMETERS: |
---|
56 | ! |
---|
57 | integer,dimension(:),pointer :: gridbuf ! local index array |
---|
58 | ! |
---|
59 | !EOP ___________________________________________________________________ |
---|
60 | |
---|
61 | integer :: npesx,npesy,ng,ny,n,i,j,nx,ig,jg,nseg,factor |
---|
62 | |
---|
63 | |
---|
64 | ! default decomp is R |
---|
65 | if((trim(ldecomp) .ne. 'R') .and. (ldecomp .ne. 'C') .and. (ldecomp .ne. 'RC')) then |
---|
66 | ldecomp = 'R' |
---|
67 | endif |
---|
68 | |
---|
69 | ! A 'by-row' or 'by-latitude' decomposition |
---|
70 | if(trim(ldecomp) .eq. 'R') then |
---|
71 | npesx=1 |
---|
72 | npesy=nprocs |
---|
73 | nx=gnx |
---|
74 | ny=gny/npesy |
---|
75 | allocate(gridbuf(nx*ny)) |
---|
76 | n=0 |
---|
77 | do j=1,ny |
---|
78 | do i=1,nx |
---|
79 | n=n+1 |
---|
80 | ig=i |
---|
81 | jg = j + myProc*ny |
---|
82 | ng =(jg-1)*gnx + ig |
---|
83 | gridbuf(n)=ng |
---|
84 | enddo |
---|
85 | enddo |
---|
86 | |
---|
87 | ! A 'by-column' or 'by-longitude' decomposition |
---|
88 | else if (ldecomp .eq. 'C') then |
---|
89 | npesx=nprocs |
---|
90 | npesy=1 |
---|
91 | nx=gnx/npesx |
---|
92 | ny=gny |
---|
93 | allocate(gridbuf(nx*ny)) |
---|
94 | n=0 |
---|
95 | do j=1,ny |
---|
96 | do i=1,nx |
---|
97 | n=n+1 |
---|
98 | ig=i + myProc*nx |
---|
99 | jg= j |
---|
100 | ng=(jg-1)*gnx + ig |
---|
101 | gridbuf(n)=ng |
---|
102 | enddo |
---|
103 | enddo |
---|
104 | |
---|
105 | ! A 'row-columen' or 'checkerboard' decomposition |
---|
106 | else if (ldecomp .eq. 'RC') then |
---|
107 | ! find the closest square |
---|
108 | factor=1 |
---|
109 | do i=2,INT(sqrt(FLOAT(nprocs))) |
---|
110 | if ( (nprocs/i) * i .eq. nprocs) then |
---|
111 | factor = i |
---|
112 | endif |
---|
113 | enddo |
---|
114 | npesx=factor |
---|
115 | npesy=nprocs/factor |
---|
116 | nx=gnx/npesx |
---|
117 | ny=gny/npesy |
---|
118 | ! write(6,*) 'RC',factor,npesy,nx,ny |
---|
119 | allocate(gridbuf(nx*ny)) |
---|
120 | n=0 |
---|
121 | do j=1,ny |
---|
122 | do i=1,nx |
---|
123 | n=n+1 |
---|
124 | ig=mod(myProc,npesx)*nx+i |
---|
125 | jg=(myProc/npesx)*ny+j |
---|
126 | ng=(jg-1)*gnx + ig |
---|
127 | gridbuf(n)=ng |
---|
128 | enddo |
---|
129 | enddo |
---|
130 | |
---|
131 | |
---|
132 | endif |
---|
133 | |
---|
134 | end subroutine get_index |
---|
135 | |
---|
136 | |
---|
137 | |
---|
138 | |
---|
139 | end module mutils |
---|