FV3DYCORE  Version 1.1.0
sorted_index.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the FV3 dynamical core.
5 !*
6 !* The FV3 dynamical core is free software: you can redistribute it
7 !* and/or modify it under the terms of the
8 !* GNU Lesser General Public License as published by the
9 !* Free Software Foundation, either version 3 of the License, or
10 !* (at your option) any later version.
11 !*
12 !* The FV3 dynamical core is distributed in the hope that it will be
13 !* useful, but WITHOUT ANYWARRANTY; without even the implied warranty
14 !* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
15 !* See the GNU General Public License for more details.
16 !*
17 !* You should have received a copy of the GNU Lesser General Public
18 !* License along with the FV3 dynamical core.
19 !* If not, see <http://www.gnu.org/licenses/>.
20 !***********************************************************************
29 
30  use fv_arrays_mod, only: r_grid
31 
32  implicit none
33  private
34  public :: sorted_inta, sorted_intb
35 
36 contains
37 
46  subroutine sorted_inta(isd, ied, jsd, jed, cubed_sphere, bgrid, iinta, jinta)
47  integer, intent(in) :: isd, ied, jsd, jed
48  real(kind=R_GRID), intent(in), dimension(isd:ied+1,jsd:jed+1,2) :: bgrid
49  logical, intent(in) :: cubed_sphere
50 
51  integer, intent(out), dimension(4,isd:ied,jsd:jed) :: iinta, jinta
52  !------------------------------------------------------------------!
53  ! local variables !
54  !------------------------------------------------------------------!
55  real, dimension(4) :: xsort, ysort
56  integer, dimension(4) :: isort, jsort
57  integer :: i, j
58  !------------------------------------------------------------------!
59  ! special treatment for cubed sphere !
60  !------------------------------------------------------------------!
61  if (cubed_sphere) then
62  !---------------------------------------------------------------!
63  ! get order of indices for line integral around a-grid cell !
64  !---------------------------------------------------------------!
65  do j=jsd,jed
66  do i=isd,ied
67  xsort(1)=bgrid(i ,j ,1); ysort(1)=bgrid(i ,j ,2); isort(1)=i ; jsort(1)=j
68  xsort(2)=bgrid(i ,j+1,1); ysort(2)=bgrid(i ,j+1,2); isort(2)=i ; jsort(2)=j+1
69  xsort(3)=bgrid(i+1,j+1,1); ysort(3)=bgrid(i+1,j+1,2); isort(3)=i+1; jsort(3)=j+1
70  xsort(4)=bgrid(i+1,j ,1); ysort(4)=bgrid(i+1,j ,2); isort(4)=i+1; jsort(4)=j
71  call sort_rectangle(iinta(1,i,j), jinta(1,i,j))
72  enddo
73  enddo
74  else
75  !---------------------------------------------------------------!
76  ! default behavior for other grids !
77  !---------------------------------------------------------------!
78  do j=jsd,jed
79  do i=isd,ied
80  iinta(i,j,1)=i ; jinta(i,j,1)=j
81  iinta(i,j,2)=i ; jinta(i,j,2)=j+1
82  iinta(i,j,3)=i+1; jinta(i,j,3)=j+1
83  iinta(i,j,4)=i+1; jinta(i,j,4)=j
84  enddo
85  enddo
86  endif
87 
88  contains
89  !------------------------------------------------------------------!
90  subroutine sort_rectangle(iind, jind)
91  integer, dimension(4), intent(inout) :: iind, jind
92  !----------------------------------------------------------------!
93  ! local variables !
94  !----------------------------------------------------------------!
95  real, dimension(4) :: xsorted, ysorted
96  integer, dimension(4) :: isorted, jsorted
97  integer :: l, ll, lll
98  !----------------------------------------------------------------!
99  ! sort in east west !
100  !----------------------------------------------------------------!
101  xsorted(:)=10.
102  ysorted(:)=10.
103  isorted(:)=0
104  jsorted(:)=0
105 
106  do l=1,4
107  do ll=1,4
108  if (xsort(l)<xsorted(ll)) then
109  do lll=3,ll,-1
110  xsorted(lll+1)=xsorted(lll)
111  ysorted(lll+1)=ysorted(lll)
112  isorted(lll+1)=isorted(lll)
113  jsorted(lll+1)=jsorted(lll)
114  enddo
115  xsorted(ll)=xsort(l)
116  ysorted(ll)=ysort(l)
117  isorted(ll)=isort(l)
118  jsorted(ll)=jsort(l)
119  exit
120  endif
121  enddo
122  enddo
123  !----------------------------------------------------------------!
124  ! sort in north south !
125  !----------------------------------------------------------------!
126  do l=1,4
127  xsort(l)=xsorted(l); ysort(l)=ysorted(l)
128  isort(l)=isorted(l); jsort(l)=jsorted(l)
129  enddo
130  xsorted(:)=10.
131  ysorted(:)=10.
132  isorted(:)=0
133  jsorted(:)=0
134 
135  do l=1,4
136  do ll=1,4
137  if (ysort(l)<ysorted(ll)) then
138  do lll=3,ll,-1
139  xsorted(lll+1)=xsorted(lll)
140  ysorted(lll+1)=ysorted(lll)
141  isorted(lll+1)=isorted(lll)
142  jsorted(lll+1)=jsorted(lll)
143  enddo
144  xsorted(ll)=xsort(l)
145  ysorted(ll)=ysort(l)
146  isorted(ll)=isort(l)
147  jsorted(ll)=jsort(l)
148  exit
149  endif
150  enddo
151  enddo
152  !----------------------------------------------------------------!
153  ! use first two grid point for start and orientation !
154  !----------------------------------------------------------------!
155  if ( isorted(1)==i .and. jsorted(1)==j ) then
156  if ( isorted(2)==i+1 .and. jsorted(2)==j+1 ) then
157  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
158  endif
159  if ( isorted(2)==i .and. jsorted(2)==j+1 ) then
160  iind(1)=i ; jind(1)=j
161  iind(2)=i ; jind(2)=j+1
162  iind(3)=i+1; jind(3)=j+1
163  iind(4)=i+1; jind(4)=j
164  elseif ( isorted(2)==i+1 .and. jsorted(2)==j ) then
165  iind(1)=i ; jind(1)=j
166  iind(2)=i+1; jind(2)=j
167  iind(3)=i+1; jind(3)=j+1
168  iind(4)=i ; jind(4)=j+1
169  endif
170 
171  elseif ( isorted(1)==i .and. jsorted(1)==j+1 ) then
172  if ( isorted(2)==i+1 .and. jsorted(2)==j ) then
173  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
174  endif
175  if ( isorted(2)==i+1 .and. jsorted(2)==j+1 ) then
176  iind(1)=i ; jind(1)=j+1
177  iind(2)=i+1; jind(2)=j+1
178  iind(3)=i+1; jind(3)=j
179  iind(4)=i ; jind(4)=j
180  elseif ( isorted(2)==i .and. jsorted(2)==j ) then
181  iind(1)=i ; jind(1)=j+1
182  iind(2)=i ; jind(2)=j
183  iind(3)=i+1; jind(3)=j
184  iind(4)=i+1; jind(4)=j+1
185  endif
186 
187  elseif ( isorted(1)==i+1 .and. jsorted(1)==j+1 ) then
188  if ( isorted(2)==i .and. jsorted(2)==j ) then
189  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
190  endif
191  if ( isorted(2)==i+1 .and. jsorted(2)==j ) then
192  iind(1)=i+1; jind(1)=j+1
193  iind(2)=i+1; jind(2)=j
194  iind(3)=i ; jind(3)=j
195  iind(4)=i ; jind(4)=j+1
196  elseif ( isorted(2)==i .and. jsorted(2)==j+1 ) then
197  iind(1)=i+1; jind(1)=j+1
198  iind(2)=i ; jind(2)=j+1
199  iind(3)=i ; jind(3)=j
200  iind(4)=i+1; jind(4)=j
201  endif
202 
203  elseif ( isorted(1)==i+1 .and. jsorted(1)==j ) then
204  if ( isorted(2)==i .and. jsorted(2)==j+1 ) then
205  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
206  endif
207  if ( isorted(2)==i .and. jsorted(2)==j ) then
208  iind(1)=i+1; jind(1)=j
209  iind(2)=i ; jind(2)=j
210  iind(3)=i ; jind(3)=j+1
211  iind(4)=i+1; jind(4)=j+1
212  elseif ( isorted(2)==i+1 .and. jsorted(2)==j+1 ) then
213  iind(1)=i+1; jind(1)=j
214  iind(2)=i+1; jind(2)=j+1
215  iind(3)=i ; jind(3)=j+1
216  iind(4)=i ; jind(4)=j
217  endif
218 
219  endif
220 
221  end subroutine sort_rectangle
222  !------------------------------------------------------------------!
223  end subroutine sorted_inta
224 
233  subroutine sorted_intb(isd, ied, jsd, jed, is, ie, js, je, npx, npy, &
234  cubed_sphere, agrid, iintb, jintb)
235  integer, intent(in) :: isd, ied, jsd, jed, is, ie, js, je, npx, npy
236  real(kind=R_GRID), intent(in), dimension(isd:ied,jsd:jed,2) :: agrid
237  logical, intent(in) :: cubed_sphere
238 
239  integer, dimension(4,is:ie+1,js:je+1), intent(out) :: iintb, jintb
240  !------------------------------------------------------------------!
241  ! local variables !
242  !------------------------------------------------------------------!
243  real, dimension(4) :: xsort, ysort, xsorted, ysorted
244  integer, dimension(4) :: isort, jsort, isorted, jsorted
245  integer :: i, j, l, ll, lll
246  !------------------------------------------------------------------!
247  ! special treatment for cubed sphere !
248  !------------------------------------------------------------------!
249  if (cubed_sphere) then
250  !---------------------------------------------------------------!
251  ! get order of indices for line integral around b-grid cell !
252  !---------------------------------------------------------------!
253  do j=js,je+1
254  do i=is,ie+1
255  xsort(1)=agrid(i ,j ,1); ysort(1)=agrid(i ,j ,2); isort(1)=i ; jsort(1)=j
256  xsort(2)=agrid(i ,j-1,1); ysort(2)=agrid(i ,j-1,2); isort(2)=i ; jsort(2)=j-1
257  xsort(3)=agrid(i-1,j-1,1); ysort(3)=agrid(i-1,j-1,2); isort(3)=i-1; jsort(3)=j-1
258  xsort(4)=agrid(i-1,j ,1); ysort(4)=agrid(i-1,j ,2); isort(4)=i-1; jsort(4)=j
259  call sort_rectangle(iintb(1,i,j), jintb(1,i,j))
260  enddo
261  enddo
262  !---------------------------------------------------------------!
263  ! take care of corner points !
264  !---------------------------------------------------------------!
265  if ( (is==1) .and. (js==1) ) then
266  i=1
267  j=1
268  xsort(1)=agrid(i ,j ,1); ysort(1)=agrid(i ,j ,2); isort(1)=i ; jsort(1)=j
269  xsort(2)=agrid(i ,j-1,1); ysort(2)=agrid(i ,j-1,2); isort(2)=i ; jsort(2)=j-1
270  xsort(3)=agrid(i-1,j ,1); ysort(3)=agrid(i-1,j ,2); isort(3)=i-1; jsort(3)=j
271  call sort_triangle()
272  iintb(4,i,j)=i-1; jintb(4,i,j)=j-1
273  endif
274 
275  if ( (ie+1==npx) .and. (js==1) ) then
276  i=npx
277  j=1
278  xsort(1)=agrid(i ,j ,1); ysort(1)=agrid(i ,j ,2); isort(1)=i ; jsort(1)=j
279  xsort(2)=agrid(i-1,j ,1); ysort(2)=agrid(i-1,j ,2); isort(2)=i-1; jsort(2)=j
280  xsort(3)=agrid(i-1,j-1,1); ysort(3)=agrid(i-1,j-1,2); isort(3)=i-1; jsort(3)=j-1
281  call sort_triangle()
282  iintb(4,i,j)=i; jintb(4,i,j)=j-1
283  endif
284 
285  if ( (ie+1==npx) .and. (je+1==npy) ) then
286  i=npx
287  j=npy
288  xsort(1)=agrid(i-1,j-1,1); ysort(1)=agrid(i-1,j-1,2); isort(1)=i-1; jsort(1)=j-1
289  xsort(2)=agrid(i ,j-1,1); ysort(2)=agrid(i ,j-1,2); isort(2)=i ; jsort(2)=j-1
290  xsort(3)=agrid(i-1,j ,1); ysort(3)=agrid(i-1,j ,2); isort(3)=i-1; jsort(3)=j
291  call sort_triangle()
292  iintb(4,i,j)=i; jintb(4,i,j)=j
293  endif
294 
295  if ( (is==1) .and. (je+1==npy) ) then
296  i=1
297  j=npy
298  xsort(1)=agrid(i ,j ,1); ysort(1)=agrid(i ,j ,2); isort(1)=i ; jsort(1)=j
299  xsort(2)=agrid(i-1,j-1,1); ysort(2)=agrid(i-1,j-1,2); isort(2)=i-1; jsort(2)=j-1
300  xsort(3)=agrid(i ,j-1,1); ysort(3)=agrid(i ,j-1,2); isort(3)=i ; jsort(3)=j-1
301  call sort_triangle()
302  iintb(4,i,j)=i-1; jintb(4,i,j)=j
303  endif
304  else
305  !---------------------------------------------------------------!
306  ! default behavior for other grids !
307  !---------------------------------------------------------------!
308  do j=js,je+1
309  do i=is,ie+1
310  iintb(1,i,j)=i ; jintb(1,i,j)=j
311  iintb(2,i,j)=i ; jintb(2,i,j)=j-1
312  iintb(3,i,j)=i-1; jintb(3,i,j)=j-1
313  iintb(4,i,j)=i-1; jintb(4,i,j)=j
314  enddo
315  enddo
316  endif
317 
318  contains
319  !------------------------------------------------------------------!
320  subroutine sort_rectangle(iind, jind)
321 
322  integer, dimension(4), intent(inout) :: iind, jind
323  !----------------------------------------------------------------!
324  ! local variables !
325  !----------------------------------------------------------------!
326  real, dimension(4) :: xsorted, ysorted
327  integer, dimension(4) :: isorted, jsorted
328  !----------------------------------------------------------------!
329  ! sort in east west !
330  !----------------------------------------------------------------!
331  xsorted(:)=10.
332  ysorted(:)=10.
333  isorted(:)=0
334  jsorted(:)=0
335 
336  do l=1,4
337  do ll=1,4
338  if (xsort(l)<xsorted(ll)) then
339  do lll=3,ll,-1
340  xsorted(lll+1)=xsorted(lll)
341  ysorted(lll+1)=ysorted(lll)
342  isorted(lll+1)=isorted(lll)
343  jsorted(lll+1)=jsorted(lll)
344  enddo
345  xsorted(ll)=xsort(l)
346  ysorted(ll)=ysort(l)
347  isorted(ll)=isort(l)
348  jsorted(ll)=jsort(l)
349  exit
350  endif
351  enddo
352  enddo
353  !----------------------------------------------------------------!
354  ! sort in north south !
355  !----------------------------------------------------------------!
356  do l=1,4
357  xsort(l)=xsorted(l); ysort(l)=ysorted(l)
358  isort(l)=isorted(l); jsort(l)=jsorted(l)
359  enddo
360  xsorted(:)=10.
361  ysorted(:)=10.
362  isorted(:)=0
363  jsorted(:)=0
364 
365  do l=1,4
366  do ll=1,4
367  if (ysort(l)<ysorted(ll)) then
368  do lll=3,ll,-1
369  xsorted(lll+1)=xsorted(lll)
370  ysorted(lll+1)=ysorted(lll)
371  isorted(lll+1)=isorted(lll)
372  jsorted(lll+1)=jsorted(lll)
373  enddo
374  xsorted(ll)=xsort(l)
375  ysorted(ll)=ysort(l)
376  isorted(ll)=isort(l)
377  jsorted(ll)=jsort(l)
378  exit
379  endif
380  enddo
381  enddo
382  !----------------------------------------------------------------!
383  ! use first two grid point for start and orientation !
384  !----------------------------------------------------------------!
385  if ( isorted(1)==i .and. jsorted(1)==j ) then
386  if ( isorted(2)==i-1 .and. jsorted(2)==j-1 ) then
387  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
388  endif
389  if ( isorted(2)==i .and. jsorted(2)==j-1 ) then
390  iind(1)=i ; jind(1)=j
391  iind(2)=i ; jind(2)=j-1
392  iind(3)=i-1; jind(3)=j-1
393  iind(4)=i-1; jind(4)=j
394  elseif ( isorted(2)==i-1 .and. jsorted(2)==j ) then
395  iind(1)=i ; jind(1)=j
396  iind(2)=i-1; jind(2)=j
397  iind(3)=i-1; jind(3)=j-1
398  iind(4)=i ; jind(4)=j-1
399  endif
400 
401  elseif ( isorted(1)==i .and. jsorted(1)==j-1 ) then
402  if ( isorted(2)==i-1 .and. jsorted(2)==j ) then
403  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
404  endif
405  if ( isorted(2)==i-1 .and. jsorted(2)==j-1 ) then
406  iind(1)=i ; jind(1)=j-1
407  iind(2)=i-1; jind(2)=j-1
408  iind(3)=i-1; jind(3)=j
409  iind(4)=i ; jind(4)=j
410  elseif ( isorted(2)==i .and. jsorted(2)==j ) then
411  iind(1)=i ; jind(1)=j-1
412  iind(2)=i ; jind(2)=j
413  iind(3)=i-1; jind(3)=j
414  iind(4)=i-1; jind(4)=j-1
415  endif
416 
417  elseif ( isorted(1)==i-1 .and. jsorted(1)==j-1 ) then
418  if ( isorted(2)==i .and. jsorted(2)==j ) then
419  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
420  endif
421  if ( isorted(2)==i-1 .and. jsorted(2)==j ) then
422  iind(1)=i-1; jind(1)=j-1
423  iind(2)=i-1; jind(2)=j
424  iind(3)=i ; jind(3)=j
425  iind(4)=i ; jind(4)=j-1
426  elseif ( isorted(2)==i .and. jsorted(2)==j-1 ) then
427  iind(1)=i-1; jind(1)=j-1
428  iind(2)=i ; jind(2)=j-1
429  iind(3)=i ; jind(3)=j
430  iind(4)=i-1; jind(4)=j
431  endif
432 
433  elseif ( isorted(1)==i-1 .and. jsorted(1)==j ) then
434  if ( isorted(2)==i .and. jsorted(2)==j-1 ) then
435  isorted(2)=isorted(3); jsorted(2)=jsorted(3)
436  endif
437  if ( isorted(2)==i .and. jsorted(2)==j ) then
438  iind(1)=i-1; jind(1)=j
439  iind(2)=i ; jind(2)=j
440  iind(3)=i ; jind(3)=j-1
441  iind(4)=i-1; jind(4)=j-1
442  elseif ( isorted(2)==i-1 .and. jsorted(2)==j-1 ) then
443  iind(1)=i-1; jind(1)=j
444  iind(2)=i-1; jind(2)=j-1
445  iind(3)=i ; jind(3)=j-1
446  iind(4)=i ; jind(4)=j
447  endif
448 
449  endif
450 
451  end subroutine sort_rectangle
452  !------------------------------------------------------------------!
453  subroutine sort_triangle()
455  xsorted(1:3)=10.
456  ysorted(1:3)=10.
457  isorted(1:3)=0
458  jsorted(1:3)=0
459  !----------------------------------------------------------------!
460  ! sort in east west !
461  !----------------------------------------------------------------!
462  do l=1,3
463  do ll=1,3
464  if (xsort(l)<xsorted(ll)) then
465  do lll=2,ll,-1
466  xsorted(lll+1)=xsorted(lll)
467  ysorted(lll+1)=ysorted(lll)
468  isorted(lll+1)=isorted(lll)
469  jsorted(lll+1)=jsorted(lll)
470  enddo
471  xsorted(ll)=xsort(l)
472  ysorted(ll)=ysort(l)
473  isorted(ll)=isort(l)
474  jsorted(ll)=jsort(l)
475  exit
476  endif
477  enddo
478  enddo
479  !----------------------------------------------------------------!
480  ! sort in north south !
481  !----------------------------------------------------------------!
482  do l=1,3
483  xsort(l)=xsorted(l); ysort(l)=ysorted(l)
484  isort(l)=isorted(l); jsort(l)=jsorted(l)
485  enddo
486  xsorted(1:3)=10.
487  ysorted(1:3)=10.
488  isorted(1:3)=0
489  jsorted(1:3)=0
490 
491  do l=1,3
492  do ll=1,3
493  if (ysort(l)<ysorted(ll)) then
494  do lll=2,ll,-1
495  xsorted(lll+1)=xsorted(lll)
496  ysorted(lll+1)=ysorted(lll)
497  isorted(lll+1)=isorted(lll)
498  jsorted(lll+1)=jsorted(lll)
499  enddo
500  xsorted(ll)=xsort(l)
501  ysorted(ll)=ysort(l)
502  isorted(ll)=isort(l)
503  jsorted(ll)=jsort(l)
504  exit
505  endif
506  enddo
507  enddo
508  !----------------------------------------------------------------!
509  ! use first two grid point for start and orientation !
510  !----------------------------------------------------------------!
511  iintb(1,i,j)=isorted(1) ; jintb(1,i,j)=jsorted(1)
512  iintb(2,i,j)=isorted(2) ; jintb(2,i,j)=jsorted(2)
513  iintb(3,i,j)=isorted(3) ; jintb(3,i,j)=jsorted(3)
514 
515  end subroutine sort_triangle
516  !------------------------------------------------------------------!
517  end subroutine sorted_intb
518 
519 end module sorted_index_mod
The module &#39;sorted_index&#39; sorts cell corner indices in lat-lon space to ensure the same order of oper...
subroutine, public sorted_intb(isd, ied, jsd, jed, is, ie, js, je, npx, npy, cubed_sphere, agrid, iintb, jintb)
The subroutine &#39;sorted_intb&#39; sorts cell corner indices in latlon space based on grid locations in ind...
integer, parameter, public r_grid
Definition: fv_arrays.F90:35
subroutine sort_triangle()
subroutine sort_rectangle(iind, jind)
The module &#39;fv_arrays&#39; contains the &#39;fv_atmos_type&#39; and associated datatypes.
Definition: fv_arrays.F90:24
subroutine, public sorted_inta(isd, ied, jsd, jed, cubed_sphere, bgrid, iinta, jinta)
The subroutine &#39;sorted_inta&#39; sorts cell corner indices in latlon space based on grid locations in ind...