Oasis3-MCT
 All Classes Files Functions Variables Macros Pages
mod_oasis_grid.F90
Go to the documentation of this file.
1 
2 !> OASIS grid data and methods
3 
4 !> These interfaces support both grid data specified globally on the root task
5 !> as required in Oasis3 and grid data decomposed across tasks. If grid data
6 !> is decomposed across tasks, the optional partid argument must be specified
7 !> when it exists in the interface.
8 
10 !-----------------------------------------------------------------------
11 ! BOP
12 !
13 ! !MODULE: mod_prism_grid
14 ! !REMARKS:
15 !
16 ! **********************
17 ! THIS SHOULD BE CALLED BY A SINGLE PE ACCORDING TO THE OASIS3
18 ! STANDARD. THE DATA IS GLOBAL.
19 ! **********************
20 !
21 ! !REVISION HISTORY:
22 !
23 !
24 ! !PUBLIC MEMBER FUNCTIONS:
25 !
26 ! subroutine oasis_start_grids_writing(iwrite)
27 ! This subroutine initializes grid writing by receiving a
28 ! starting command from OASIS.
29 !
30 ! subroutine oasis_write_grid(cgrid, nx, ny, lon, lat, part_id)
31 ! This subroutine writes longitudes and latitudes for a model
32 ! grid. part_id is optional and indicates decomposed data is provided
33 !
34 ! subroutine oasis_write_corner(cgrid, nx, ny, nc, clon, clat, part_id)
35 ! This subroutine writes the longitudes and latitudes of the
36 ! grid cell corners. part_id is optional and indicates decomposed data
37 ! is provided
38 !
39 ! subroutine oasis_write_mask(cgrid, nx, ny, mask, part_id)
40 ! This subroutine writes the mask for a model grid. part_id is optional
41 ! and indicates decomposed data is provided
42 !
43 ! subroutine oasis_write_area(cgrid, nx, ny, area, part_id)
44 ! This subroutine writes the grid cell areas for a model grid. part_id
45 ! is optional and indicates decomposed data is provided
46 !
47 ! subroutine oasis_terminate_grids_writing()
48 ! This subroutine terminates grid writing by sending a flag
49 ! to OASIS, stating the all needed grid information was written.
50 !
51 
52 ! !USES:
53  USE mod_oasis_data
54  USE mod_oasis_io
55  USE mod_oasis_sys
56  USE mod_oasis_part
59  USE mod_oasis_timer
60  USE mct_mod
61 
62  implicit none
63 
64  private
65 
67  public oasis_write_grid
68  public oasis_write_angle
69  public oasis_write_corner
70  public oasis_write_mask
71  public oasis_write_area
73  public oasis_write2files
75 
76 #include "oasis_os.h"
77 
78  !> Generic interface to support writing 4 or 8 byte reals
79  interface oasis_write_grid
80 #ifndef __NO_4BYTE_REALS
81  module procedure oasis_write_grid_r4
82 #endif
83  module procedure oasis_write_grid_r8
84  end interface
85 
86  !> Generic interface to support writing 4 or 8 byte reals
88 #ifndef __NO_4BYTE_REALS
89  module procedure oasis_write_angle_r4
90 #endif
91  module procedure oasis_write_angle_r8
92  end interface
93 
94  !> Generic interface to support writing 4 or 8 byte reals
96 #ifndef __NO_4BYTE_REALS
97  module procedure oasis_write_corner_r4
98 #endif
99  module procedure oasis_write_corner_r8
100  end interface
101 
102  !> Generic interface to support writing 4 or 8 byte reals
104 #ifndef __NO_4BYTE_REALS
105  module procedure oasis_write_area_r4
106 #endif
107  module procedure oasis_write_area_r8
108  end interface
109 
110  !--- datatypes ---
111  public :: prism_grid_type !< Grid datatype
112 
113  integer(kind=ip_intwp_p),parameter :: mgrid = 100 !< maximum number of grids allowed
114  integer(kind=ip_intwp_p),save :: writing_grids_call=0
115 
116  !> Model grid data for creating mapping data and conserving fields
118  character(len=ic_med) :: gridname !< grid name
119  integer(kind=ip_i4_p) :: partid !< partition ID
120  integer(kind=ip_i4_p) :: nx !< global nx size
121  integer(kind=ip_i4_p) :: ny !< global ny size
122  integer(kind=ip_i4_p) :: nc !< number of corners per gridcell
123  logical :: grid_set !< flag to track user calls for grid
124  logical :: corner_set !< flag to track user calls for corner
125  logical :: angle_set !< flag to track user calls for angle
126  logical :: area_set !< flag to track user calls for area
127  logical :: mask_set !< flag to track user calls for mask
128  logical :: written !< flag to indicate grid has been written
129  logical :: terminated !< flag to indicate user grid calls complete
130  real(kind=ip_realwp_p),allocatable :: lon(:,:) !< user specified longitudes
131  real(kind=ip_realwp_p),allocatable :: lat(:,:) !< user specified latitudes
132  real(kind=ip_realwp_p),allocatable :: clon(:,:,:) !< user specified corner longitudes
133  real(kind=ip_realwp_p),allocatable :: clat(:,:,:) !< user specified corner latitudes
134  real(kind=ip_realwp_p),allocatable :: angle(:,:) !< user specified angle
135  real(kind=ip_realwp_p),allocatable :: area(:,:) !< user specified area
136  integer(kind=ip_i4_p) ,allocatable :: mask(:,:) !< user specified mask
137  end type prism_grid_type
138 
139  integer(kind=ip_intwp_p),public,save :: prism_ngrid = 0 !< counter for grids
140  type(prism_grid_type),public,save :: prism_grid(mgrid) !< array of grid datatypes
141 
142 
143 #ifdef use_netCDF
144 #include <netcdf.inc>
145 #endif
146 
147 !---------------------------------------------------------------------------
148 
149 CONTAINS
150 
151 !--------------------------------------------------------------------------
152 
153 !> Print grid information to log file.
154 
156 
157  implicit none
158 
159  !-------------------------------------------------
160  integer(kind=ip_intwp_p) :: n
161  character(len=*),parameter :: subname = '(oasis_print_grid_data)'
162  !-------------------------------------------------
163 
164  call oasis_debug_enter(subname)
165 
166  do n = 1,prism_ngrid
167  write(nulprt,*) ' '
168  write(nulprt,*) subname,trim(prism_grid(n)%gridname),' size',prism_grid(n)%nx,prism_grid(n)%ny
169  write(nulprt,*) subname,trim(prism_grid(n)%gridname),' set ',prism_grid(n)%grid_set, &
170  prism_grid(n)%corner_set, prism_grid(n)%angle_set, prism_grid(n)%area_set, prism_grid(n)%mask_set
171  if (prism_grid(n)%partid > 0 .and. prism_grid(n)%partid < prism_npart) then
172  write(nulprt,*) subname,'partid ',trim(prism_grid(n)%gridname),prism_grid(n)%partid, &
173  trim(prism_part(prism_grid(n)%partid)%partname)
174  else
175  write(nulprt,*) subname,'partid ',trim(prism_grid(n)%gridname),prism_grid(n)%partid
176  endif
177  if (prism_grid(n)%grid_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' lon ', &
178  size(prism_grid(n)%lon,dim=1),size(prism_grid(n)%lon,dim=2), &
179  minval(prism_grid(n)%lon),maxval(prism_grid(n)%lon)
180  if (prism_grid(n)%grid_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' lat ', &
181  size(prism_grid(n)%lat,dim=1),size(prism_grid(n)%lat,dim=2), &
182  minval(prism_grid(n)%lat),maxval(prism_grid(n)%lat)
183  if (prism_grid(n)%corner_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' clon ', &
184  size(prism_grid(n)%clon,dim=1),size(prism_grid(n)%clon,dim=2), &
185  minval(prism_grid(n)%clon),maxval(prism_grid(n)%clon)
186  if (prism_grid(n)%corner_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' clat ', &
187  size(prism_grid(n)%clat,dim=1),size(prism_grid(n)%clat,dim=2), &
188  minval(prism_grid(n)%clat),maxval(prism_grid(n)%clat)
189  if (prism_grid(n)%angle_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' angl ', &
190  size(prism_grid(n)%angle,dim=1),size(prism_grid(n)%angle,dim=2), &
191  minval(prism_grid(n)%angle),maxval(prism_grid(n)%angle)
192  if (prism_grid(n)%mask_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' mask ', &
193  size(prism_grid(n)%mask,dim=1),size(prism_grid(n)%mask,dim=2), &
194  minval(prism_grid(n)%mask),maxval(prism_grid(n)%mask)
195  if (prism_grid(n)%area_set) write(nulprt,*) subname,trim(prism_grid(n)%gridname),' area ', &
196  size(prism_grid(n)%area,dim=1),size(prism_grid(n)%area,dim=2), &
197  minval(prism_grid(n)%area),maxval(prism_grid(n)%area)
198  enddo
199 
200  call oasis_debug_exit(subname)
201 
202  END SUBROUTINE oasis_print_grid_data
203 
204 !--------------------------------------------------------------------------
205 
206 !> User interface to initialize grid writing
207 
208  SUBROUTINE oasis_start_grids_writing(iwrite)
209 
210  implicit none
211 
212  integer(kind=ip_intwp_p), intent (OUT) :: iwrite !< flag, obsolete
213 
214  !-------------------------------------------------
215  character(len=*),parameter :: subname = '(oasis_start_grids_writing)'
216  !-------------------------------------------------
217 
218  call oasis_debug_enter(subname)
219 
220  if (oasis_debug >= 15) then
221  write(nulprt,*) subname,' prism_ngrid = ',prism_ngrid
222  endif
223 
224  if (prism_ngrid == 0) then ! first call
225  prism_grid(:)%gridname = 'unSet'
226  prism_grid(:)%nx = -1
227  prism_grid(:)%ny = -1
228  prism_grid(:)%grid_set = .false.
229  prism_grid(:)%corner_set = .false.
230  prism_grid(:)%angle_set = .false.
231  prism_grid(:)%area_set = .false.
232  prism_grid(:)%mask_set = .false.
233  prism_grid(:)%written = .false.
234  prism_grid(:)%terminated = .false.
235  prism_grid(:)%partid = -1
236  endif
237  iwrite = 1 ! just set grids are needed always
238  writing_grids_call=1
239 
240  call oasis_debug_exit(subname)
241 
242  END SUBROUTINE oasis_start_grids_writing
243 
244 !--------------------------------------------------------------------------
245 
246 !> User interface to set latitudes and longitudes for 8 byte reals
247 
248  SUBROUTINE oasis_write_grid_r8(cgrid, nx, ny, lon, lat, partid)
249 
250  !-------------------------------------------------
251  ! Routine to create a new grids file or to add a grid description to an
252  ! existing grids file.
253  !-------------------------------------------------
254 
255  implicit none
256 
257  character(len=*), intent (in) :: cgrid !< grid name
258  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
259  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
260  real(kind=ip_double_p), intent (in) :: lon(:,:) !< longitudes
261  real(kind=ip_double_p), intent (in) :: lat(:,:) !< latitudes
262  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
263  !-------------------------------------------------
264  integer(kind=ip_intwp_p) :: GRIDID
265  integer(kind=ip_intwp_p) :: ierror
266  integer(kind=ip_intwp_p) :: lnx,lny
267  character(len=*),parameter :: subname = '(oasis_write_grid_r8)'
268  !-------------------------------------------------
269 
270  call oasis_debug_enter(subname)
271 
272  if (oasis_debug >= 15) then
273  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
274  endif
275 
276  call oasis_findgrid(cgrid,nx,ny,gridid)
277 
278  lnx = size(lon,dim=1)
279  lny = size(lon,dim=2)
280 
281  allocate(prism_grid(gridid)%lon(lnx,lny),stat=ierror)
282  IF (ierror /= 0) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
283  mpi_rank_local,' WARNING lon alloc'
284 
285  lnx = size(lat,dim=1)
286  lny = size(lat,dim=2)
287 
288  allocate(prism_grid(gridid)%lat(lnx,lny),stat=ierror)
289  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
290  mpi_rank_local,' WARNING lat alloc'
291 
292  prism_grid(gridid)%lon = lon
293  prism_grid(gridid)%lat = lat
294  prism_grid(gridid)%grid_set = .true.
295 
296  if (present(partid)) then
297  if (prism_grid(gridid)%partid > 0 .and. prism_grid(gridid)%partid /= partid) then
298  write(nulprt,*) subname,estr,'partid inconsistency',gridid,prism_grid(gridid)%partid,partid
299  call oasis_abort()
300  endif
301  prism_grid(gridid)%partid = partid
302  if (oasis_debug >= 15) then
303  write(nulprt,*) subname,' partid = ',trim(cgrid),partid
304  endif
305  endif
306 
307  call oasis_debug_exit(subname)
308 
309  END SUBROUTINE oasis_write_grid_r8
310 
311 !--------------------------------------------------------------------------
312 
313 !> User interface to set latitudes and longitudes for 4 byte reals
314 
315  SUBROUTINE oasis_write_grid_r4(cgrid, nx, ny, lon, lat, partid)
316 
317  !-------------------------------------------------
318  ! Routine to create a new grids file or to add a grid description to an
319  ! existing grids file.
320  !-------------------------------------------------
321 
322  implicit none
323 
324  character(len=*), intent (in) :: cgrid !< grid name
325  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
326  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
327  real(kind=ip_single_p), intent (in) :: lon(:,:) !< longitudes
328  real(kind=ip_single_p), intent (in) :: lat(:,:) !< latitudes
329  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
330  !-------------------------------------------------
331  real(kind=ip_double_p), allocatable :: lon8(:,:)
332  real(kind=ip_double_p), allocatable :: lat8(:,:)
333  integer(kind=ip_intwp_p) :: ierror
334  integer(kind=ip_intwp_p) :: lpartid
335  integer(kind=ip_intwp_p) :: lnx,lny
336  character(len=*),parameter :: subname = '(oasis_write_grid_r4)'
337  !-------------------------------------------------
338 
339  call oasis_debug_enter(subname)
340 
341  if (oasis_debug >= 15) then
342  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
343  endif
344 
345  lpartid = -1
346  if (present(partid)) then
347  lpartid = partid
348  endif
349  if (oasis_debug >= 15) then
350  write(nulprt,*) subname,' partid = ',trim(cgrid),lpartid
351  endif
352 
353  lnx = size(lon,dim=1)
354  lny = size(lon,dim=2)
355 
356  allocate(lon8(lnx,lny),stat=ierror)
357  IF (ierror /= 0) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
358  mpi_rank_local,' WARNING lon alloc'
359 
360  lnx = size(lat,dim=1)
361  lny = size(lat,dim=2)
362 
363  allocate(lat8(lnx,lny),stat=ierror)
364  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
365  mpi_rank_local,' WARNING lat alloc'
366 
367  lon8 = lon
368  lat8 = lat
369  call oasis_write_grid_r8(cgrid,nx,ny,lon8,lat8,partid=lpartid)
370  deallocate(lon8)
371  deallocate(lat8)
372 
373  call oasis_debug_exit(subname)
374 
375  END SUBROUTINE oasis_write_grid_r4
376 
377 !--------------------------------------------------------------------------
378 
379 !> User interface to set angle for 8 byte reals
380 
381  SUBROUTINE oasis_write_angle_r8(cgrid, nx, ny, angle, partid)
382 
383  !-------------------------------------------------
384  ! Routine to add angles to an existing grid file.
385  !-------------------------------------------------
386 
387  implicit none
388 
389  character(len=*), intent (in) :: cgrid !< grid name
390  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
391  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
392  real(kind=ip_double_p), intent (in) :: angle(:,:) !< angles
393  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
394  !-------------------------------------------------
395  integer(kind=ip_intwp_p) :: GRIDID
396  integer(kind=ip_intwp_p) :: ierror
397  integer(kind=ip_intwp_p) :: lnx,lny
398  character(len=*),parameter :: subname = '(oasis_write_angle_r8)'
399  !-------------------------------------------------
400 
401  call oasis_debug_enter(subname)
402 
403  if (oasis_debug >= 15) then
404  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
405  endif
406 
407  call oasis_findgrid(cgrid,nx,ny,gridid)
408 
409  lnx = size(angle,dim=1)
410  lny = size(angle,dim=2)
411 
412  allocate(prism_grid(gridid)%angle(lnx,lny),stat=ierror)
413  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
414  mpi_rank_local,' WARNING angle alloc'
415 
416  prism_grid(gridid)%angle = angle
417  prism_grid(gridid)%angle_set = .true.
418  if (present(partid)) then
419  if (prism_grid(gridid)%partid > 0 .and. prism_grid(gridid)%partid /= partid) then
420  write(nulprt,*) subname,estr,'partid inconsistency',gridid,prism_grid(gridid)%partid,partid
421  call oasis_abort()
422  endif
423  prism_grid(gridid)%partid = partid
424  if (oasis_debug >= 15) then
425  write(nulprt,*) subname,' partid = ',trim(cgrid),partid
426  endif
427  endif
428 
429  call oasis_debug_exit(subname)
430 
431  END SUBROUTINE oasis_write_angle_r8
432 
433 !--------------------------------------------------------------------------
434 
435 !> User interface to set angle for 4 byte reals
436 
437  SUBROUTINE oasis_write_angle_r4(cgrid, nx, ny, angle, partid)
438 
439  !-------------------------------------------------
440  ! Routine to add angles to an existing grid file.
441  !-------------------------------------------------
442 
443  implicit none
444 
445  character(len=*), intent (in) :: cgrid !< grid name
446  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
447  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
448  real(kind=ip_single_p), intent (in) :: angle(:,:) !< angles
449  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
450  !-------------------------------------------------
451  real(kind=ip_double_p),allocatable :: angle8(:,:)
452  integer(kind=ip_intwp_p) :: ierror
453  integer(kind=ip_intwp_p) :: lpartid
454  integer(kind=ip_intwp_p) :: lnx,lny
455  character(len=*),parameter :: subname = '(oasis_write_angle_r4)'
456  !-------------------------------------------------
457 
458  call oasis_debug_enter(subname)
459 
460  if (oasis_debug >= 15) then
461  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
462  endif
463 
464  lpartid = -1
465  if (present(partid)) then
466  lpartid = partid
467  endif
468  if (oasis_debug >= 15) then
469  write(nulprt,*) subname,' partid = ',trim(cgrid),lpartid
470  endif
471 
472  lnx = size(angle,dim=1)
473  lny = size(angle,dim=2)
474 
475  allocate(angle8(lnx,lny),stat=ierror)
476  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
477  mpi_rank_local,' WARNING angle8 alloc'
478 
479  angle8 = angle
480  call oasis_write_angle_r8(cgrid,nx,ny,angle8,partid=lpartid)
481 
482  deallocate(angle8)
483 
484  call oasis_debug_exit(subname)
485 
486  END SUBROUTINE oasis_write_angle_r4
487 
488 !--------------------------------------------------------------------------
489 
490 !> User interface to set corner latitudes and longitudes for 8 byte reals
491 
492  SUBROUTINE oasis_write_corner_r8(cgrid, nx, ny, nc, clon, clat, partid)
493 
494  !-------------------------------------------------
495  ! Routine to add longitudes and latitudes of grid cell corners to an
496  ! existing grids file.
497  !-------------------------------------------------
498 
499  implicit none
500 
501  character(len=*), intent (in) :: cgrid !< grid name
502  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
503  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
504  integer(kind=ip_intwp_p), intent (in) :: nc !< number of corners per cell
505  real(kind=ip_double_p), intent (in) :: clon(:,:,:) !< corner longitudes
506  real(kind=ip_double_p), intent (in) :: clat(:,:,:) !< corner latitudes
507  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
508  !-------------------------------------------------
509  integer(kind=ip_intwp_p) :: GRIDID
510  integer(kind=ip_intwp_p) :: ierror
511  integer(kind=ip_intwp_p) :: lnx,lny
512  character(len=*),parameter :: subname = '(oasis_write_corner_r8)'
513  !-------------------------------------------------
514 
515  call oasis_debug_enter(subname)
516 
517  if (oasis_debug >= 15) then
518  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
519  endif
520 
521  call oasis_findgrid(cgrid,nx,ny,gridid)
522 
523  lnx = size(clon,dim=1)
524  lny = size(clon,dim=2)
525 
526  allocate(prism_grid(gridid)%clon(lnx,lny,nc),stat=ierror)
527  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
528  mpi_rank_local,' WARNING clon alloc'
529 
530  lnx = size(clat,dim=1)
531  lny = size(clat,dim=2)
532 
533  allocate(prism_grid(gridid)%clat(lnx,lny,nc),stat=ierror)
534  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
535  mpi_rank_local,' WARNING clat alloc'
536 
537  prism_grid(gridid)%nc = nc
538  prism_grid(gridid)%clon = clon
539  prism_grid(gridid)%clat = clat
540  prism_grid(gridid)%corner_set = .true.
541  if (present(partid)) then
542  if (prism_grid(gridid)%partid > 0 .and. prism_grid(gridid)%partid /= partid) then
543  write(nulprt,*) subname,estr,'partid inconsistency',gridid,prism_grid(gridid)%partid,partid
544  call oasis_abort()
545  endif
546  prism_grid(gridid)%partid = partid
547  if (oasis_debug >= 15) then
548  write(nulprt,*) subname,' partid = ',trim(cgrid),partid
549  endif
550  endif
551 
552  call oasis_debug_exit(subname)
553 
554  END SUBROUTINE oasis_write_corner_r8
555 
556 !--------------------------------------------------------------------------
557 
558 !> User interface to set corner latitudes and longitudes for 4 byte reals
559 
560  SUBROUTINE oasis_write_corner_r4(cgrid, nx, ny, nc, clon, clat, partid)
561 
562  !-------------------------------------------------
563  ! Routine to add longitudes and latitudes of grid cell corners to an
564  ! existing grids file.
565  !-------------------------------------------------
566 
567  implicit none
568 
569  character(len=*), intent (in) :: cgrid !< grid name
570  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
571  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
572  integer(kind=ip_intwp_p), intent (in) :: nc !< number of corners per cell
573  real(kind=ip_single_p), intent (in) :: clon(:,:,:) !< corner longitudes
574  real(kind=ip_single_p), intent (in) :: clat(:,:,:) !< corner latitudes
575  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
576  !-------------------------------------------------
577  real(kind=ip_double_p), allocatable :: clon8(:,:,:),clat8(:,:,:)
578  integer(kind=ip_intwp_p) :: ierror
579  integer(kind=ip_intwp_p) :: lpartid
580  integer(kind=ip_intwp_p) :: lnx,lny
581  character(len=*),parameter :: subname = '(oasis_write_corner_r4)'
582  !-------------------------------------------------
583 
584  call oasis_debug_enter(subname)
585 
586  if (oasis_debug >= 15) then
587  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
588  endif
589 
590  lpartid = -1
591  if (present(partid)) then
592  lpartid = partid
593  endif
594  if (oasis_debug >= 15) then
595  write(nulprt,*) subname,' partid = ',trim(cgrid),lpartid
596  endif
597 
598  lnx = size(clon,dim=1)
599  lny = size(clon,dim=2)
600 
601  allocate(clon8(lnx,lny,nc),stat=ierror)
602  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
603  mpi_rank_local,' WARNING clon8 alloc'
604 
605  lnx = size(clat,dim=1)
606  lny = size(clat,dim=2)
607 
608  allocate(clat8(lnx,lny,nc),stat=ierror)
609  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
610  mpi_rank_local,' WARNING clat8 alloc'
611 
612  clon8 = clon
613  clat8 = clat
614  call oasis_write_corner_r8(cgrid,nx,ny,nc,clon8,clat8,partid=lpartid)
615 
616  deallocate(clon8)
617  deallocate(clat8)
618 
619  call oasis_debug_exit(subname)
620 
621  END SUBROUTINE oasis_write_corner_r4
622 
623 !--------------------------------------------------------------------------
624 
625 !> User interface to set integer mask values
626 
627  SUBROUTINE oasis_write_mask(cgrid, nx, ny, mask, partid)
628 
629  !-------------------------------------------------
630  ! Routine to create a new masks file or to add a land see mask to an
631  ! existing masks file.
632  !-------------------------------------------------
633 
634  implicit none
635 
636  character(len=*), intent (in) :: cgrid !< grid name
637  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
638  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
639  integer(kind=ip_intwp_p), intent (in) :: mask(:,:) !< mask
640  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
641  !-------------------------------------------------
642  integer(kind=ip_intwp_p) :: GRIDID
643  integer(kind=ip_intwp_p) :: ierror
644  integer(kind=ip_intwp_p) :: lnx,lny
645  character(len=*),parameter :: subname = '(oasis_write_mask)'
646  !-------------------------------------------------
647 
648  call oasis_debug_enter(subname)
649 
650  if (oasis_debug >= 15) then
651  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
652  endif
653 
654  call oasis_findgrid(cgrid,nx,ny,gridid)
655 
656  lnx = size(mask,dim=1)
657  lny = size(mask,dim=2)
658 
659  allocate(prism_grid(gridid)%mask(lnx,lny),stat=ierror)
660  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
661  mpi_rank_local,' WARNING mask alloc'
662 
663  prism_grid(gridid)%mask = mask
664  prism_grid(gridid)%mask_set = .true.
665  if (present(partid)) then
666  if (prism_grid(gridid)%partid > 0 .and. prism_grid(gridid)%partid /= partid) then
667  write(nulprt,*) subname,estr,'partid inconsistency',gridid,prism_grid(gridid)%partid,partid
668  call oasis_abort()
669  endif
670  prism_grid(gridid)%partid = partid
671  if (oasis_debug >= 15) then
672  write(nulprt,*) subname,' partid = ',trim(cgrid),partid
673  endif
674  endif
675 
676  call oasis_debug_exit(subname)
677 
678  END SUBROUTINE oasis_write_mask
679 
680 !--------------------------------------------------------------------------
681 
682 !> User interface to set area values for 8 byte reals
683 
684  SUBROUTINE oasis_write_area_r8(cgrid, nx, ny, area, partid)
685 
686  !-------------------------------------------------
687  ! Routine to create a new areas file or to add areas of a grid to an
688  ! existing areas file.
689  !-------------------------------------------------
690 
691  implicit none
692 
693  character(len=*), intent (in) :: cgrid !< grid name
694  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
695  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
696  real(kind=ip_double_p), intent (in) :: area(:,:) !< areas
697  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
698  !-------------------------------------------------
699  integer(kind=ip_intwp_p) :: GRIDID
700  integer(kind=ip_intwp_p) :: ierror
701  integer(kind=ip_intwp_p) :: lnx,lny
702  character(len=*),parameter :: subname = '(oasis_write_area_r8)'
703  !-------------------------------------------------
704 
705  call oasis_debug_enter(subname)
706 
707  if (oasis_debug >= 15) then
708  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
709  endif
710 
711  call oasis_findgrid(cgrid,nx,ny,gridid)
712 
713  lnx = size(area,dim=1)
714  lny = size(area,dim=2)
715 
716  allocate(prism_grid(gridid)%area(lnx,lny),stat=ierror)
717  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
718  mpi_rank_local,' WARNING area alloc'
719 
720  prism_grid(gridid)%area = area
721  prism_grid(gridid)%area_set = .true.
722  if (present(partid)) then
723  if (prism_grid(gridid)%partid > 0 .and. prism_grid(gridid)%partid /= partid) then
724  write(nulprt,*) subname,estr,'partid inconsistency',gridid,prism_grid(gridid)%partid,partid
725  call oasis_abort()
726  endif
727  prism_grid(gridid)%partid = partid
728  if (oasis_debug >= 15) then
729  write(nulprt,*) subname,' partid = ',trim(cgrid),partid
730  endif
731  endif
732 
733  call oasis_debug_exit(subname)
734 
735  END SUBROUTINE oasis_write_area_r8
736 
737 !--------------------------------------------------------------------------
738 
739 !> User interface to set area values for 4 byte reals
740 
741  SUBROUTINE oasis_write_area_r4(cgrid, nx, ny, area, partid)
742 
743  !-------------------------------------------------
744  ! Routine to create a new areas file or to add areas of a grid to an
745  ! existing areas file.
746  !-------------------------------------------------
747 
748  implicit none
749 
750  character(len=*), intent (in) :: cgrid !< grid name
751  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
752  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
753  real(kind=ip_single_p), intent (in) :: area(:,:) !< areas
754  integer(kind=ip_intwp_p), intent (in),optional :: partid !< partition id if nonglobal data
755  !-------------------------------------------------
756  real(kind=ip_double_p), allocatable :: area8(:,:)
757  integer(kind=ip_intwp_p) :: ierror
758  integer(kind=ip_intwp_p) :: lpartid
759  integer(kind=ip_intwp_p) :: lnx,lny
760  character(len=*),parameter :: subname = '(oasis_write_area_r4)'
761  !-------------------------------------------------
762 
763  call oasis_debug_enter(subname)
764 
765  if (oasis_debug >= 15) then
766  write(nulprt,*) subname,' size = ',trim(cgrid),nx,ny
767  endif
768 
769  lpartid = -1
770  if (present(partid)) then
771  lpartid = partid
772  endif
773  if (oasis_debug >= 15) then
774  write(nulprt,*) subname,' partid = ',trim(cgrid),lpartid
775  endif
776 
777  lnx = size(area,dim=1)
778  lny = size(area,dim=2)
779 
780  allocate(area8(lnx,lny),stat=ierror)
781  if (ierror /= 0) write(nulprt,*) subname,' model :',compid,' proc :',&
782  mpi_rank_local,' WARNING area8 alloc'
783 
784  area8 = area
785  call oasis_write_area_r8(cgrid,nx,ny,area8,partid=lpartid)
786 
787  deallocate(area8)
788 
789  call oasis_debug_exit(subname)
790 
791  END SUBROUTINE oasis_write_area_r4
792 
793 !--------------------------------------------------------------------------
794 
795 !> User interface to indicate user defined grids are done
796 
798  !-------------------------------------------------
799  ! Routine to terminate the grids writing.
800  !-------------------------------------------------
801 
802  implicit none
803  integer(kind=ip_i4_p) :: n
804  character(len=*),parameter :: subname = '(oasis_terminate_grids_writing)'
805 
806  call oasis_debug_enter(subname)
807 
808  if (oasis_debug >= 15) then
809  write(nulprt,*) subname,' prism_ngrid = ',prism_ngrid
810  endif
811 
812  do n = 1,prism_ngrid
813  prism_grid(n)%terminated = .true.
814  enddo
815 
816  call oasis_print_grid_data()
817 ! moved to prism_method_enddef for synchronization
818 ! call oasis_write2files()
819 
820  call oasis_debug_exit(subname)
821 
822  END SUBROUTINE oasis_terminate_grids_writing
823 
824 !--------------------------------------------------------------------------
825 
826 !> Interface that actually writes fields to grid files
827 
828  SUBROUTINE oasis_write2files()
829 
830  !-------------------------------------------------
831  !> Write fields to grid files.
832  !> Only write fields that have been buffered and
833  !> if prism_grid_terminate_grids_writing has been called.
834  !> This is called by all tasks from oasis_enddef.
835  !-------------------------------------------------
836 
837  implicit none
838 
839  !-------------------------------------------------
840  character(len=ic_med) :: filename ! grid filename
841  character(len=ic_med) :: fldname ! full field name
842  character(len=ic_med) :: cgrid ! grid name
843  logical :: exists ! check if file exists
844  integer(kind=ip_i4_p) :: m,n,n1,g,p ! counter
845  integer(kind=ip_i4_p) :: partid ! part id
846  integer(kind=ip_i4_p) :: taskid ! task id for writing
847  integer(kind=ip_i4_p) :: nx,ny,nc ! grid size
848  integer(kind=ip_i4_p) :: tnx,tny ! temporary grid size
849  logical :: partid_grid ! global or local grid
850  logical :: active_task ! determine which tasks are involved
851  logical :: write_task ! task for writing
852  real(kind=ip_realwp_p),allocatable :: rloc(:,:) ! local array
853  real(kind=ip_realwp_p),allocatable :: rglo(:,:) ! global array
854  real(kind=ip_realwp_p),allocatable :: r3glo(:,:,:) ! global array
855  integer(kind=ip_i4_p) ,allocatable :: iglo(:,:) ! global array
856  integer(kind=ip_intwp_p) :: gcnt
857  logical :: found
858  character(len=ic_med) ,pointer :: gname0(:),gname(:)
859  character(len=ic_lvar2) ,pointer :: pname0(:),pname(:)
860  logical, parameter :: local_timers_on = .false.
861  character(len=*),parameter :: undefined_partname = '(UnDeFiNeD_PArtnaME)'
862  character(len=*),parameter :: subname = '(oasis_write2files)'
863  !-------------------------------------------------
864 
865  call oasis_debug_enter(subname)
866  call oasis_timer_start('grid_write')
867 
868  call oasis_mpi_bcast(writing_grids_call,mpi_comm_local,subname//'writing_grids_call')
869  if (writing_grids_call .eq. 1) then
870 
871  call oasis_timer_start('grid_write_reducelists')
872  allocate(gname0(prism_ngrid))
873  allocate(pname0(prism_ngrid))
874  do n = 1,prism_ngrid
875  gname0(n) = prism_grid(n)%gridname
876  if (prism_grid(n)%partid > 0 .and. prism_grid(n)%partid <= prism_npart) then
877  pname0(n) = prism_part(prism_grid(n)%partid)%partname
878  elseif (prism_grid(n)%partid == -1) then
879  pname0(n) = undefined_partname
880  else
881  write(nulprt,*) subname,estr,'illegal partition id for grid ',trim(prism_grid(n)%gridname),prism_grid(n)%partid
882  endif
883  enddo
884 
885  call oasis_mpi_reducelists(gname0,mpi_comm_local,gcnt,gname,'write2files',fastcheck=.true., &
886  linp2=pname0,lout2=pname,spval2=undefined_partname)
887  deallocate(gname0)
888  deallocate(pname0)
889  call oasis_timer_stop('grid_write_reducelists')
890 
891  !-------------------------------------
892  !> * Check that a grid defined on a partitition is defined on all tasks on that partition.
893  !-------------------------------------
894 
895  do n = 1,gcnt
896  if (pname(n) /= undefined_partname) then
897  do p = 1,prism_npart
898  if (pname(n) == prism_part(p)%partname .and. prism_part(p)%mpicom /= mpi_comm_null) then
899  found = .false.
900  do g = 1,prism_ngrid
901  if (prism_grid(g)%gridname == gname(n)) found = .true.
902  enddo
903  if (.not. found) then
904  write(nulprt,*) subname,estr,'grid with partition not defined on all partition tasks: ',trim(gname(n))
905  call oasis_abort()
906  endif
907  endif
908  enddo
909  endif
910  enddo
911 
912  if (local_timers_on) call oasis_timer_start('grid_write_writefiles')
913 
914  !-------------------------------------
915  !> * Write grid information
916  !-------------------------------------
917 
918  do g = 1,gcnt
919  do n = 1,prism_ngrid
920  if (prism_grid(n)%terminated) then
921  if (prism_grid(n)%gridname == gname(g)) then
922 
923  cgrid = gname(g)
924  partid = prism_grid(n)%partid
925  prism_grid(n)%written = .true.
926  tnx = -1
927  tny = -1
928 
929  !-------------------------------------
930  !> * Determine which tasks are associated with the grid information
931  !-------------------------------------
932 
933  active_task = .false.
934  write_task = .false.
935  if (pname(g) == undefined_partname) then
936  partid_grid = .false.
937  if (mpi_rank_local == 0) active_task = .true.
938  if (mpi_rank_local == 0) write_task = .true.
939  else
940  partid_grid = .true.
941  if (partid > 0 .and. partid <= prism_npart) then
942  taskid = 0
943  if (prism_part(partid)%mpicom /= mpi_comm_null) active_task = .true.
944  if (prism_part(partid)%rank == taskid) write_task = .true.
945  elseif (partid == -1) then
946  active_task = .false.
947  write_task = .false.
948  else
949  write(nulprt,*) subname,estr,'illegal partid for grid:',trim(gname(g)),trim(pname(g)),partid
950  call oasis_abort()
951  endif
952  endif
953 
954  if (oasis_debug >= 15) then
955  write(nulprt,*) subname,' ',trim(gname(g)),':',trim(pname(g)),': partid_grid=', &
956  partid_grid,'active_task=',active_task,'write_task=',write_task
957  endif
958 
959  if (active_task) then
960 
961  nx = prism_grid(n)%nx
962  ny = prism_grid(n)%ny
963  nc = prism_grid(n)%nc
964 
965  allocate(rglo(nx,ny))
966 
967  !-------------------------------------
968  !> * Check that array sizes match for all fields
969  !-------------------------------------
970 
971  if (prism_grid(n)%grid_set) then
972  if (tnx <= 0 .or. tny <= 0) then
973  tnx = size(prism_grid(n)%lon,dim=1)
974  tny = size(prism_grid(n)%lon,dim=2)
975  endif
976  if (size(prism_grid(n)%lon,dim=1) /= tnx .or. &
977  size(prism_grid(n)%lon,dim=2) /= tny .or. &
978  size(prism_grid(n)%lat,dim=1) /= tnx .or. &
979  size(prism_grid(n)%lat,dim=2) /= tny ) then
980  write(nulprt,*) subname,estr,'inconsistent array size lon/lat ',tnx,tny, &
981  size(prism_grid(n)%lon,dim=1),size(prism_grid(n)%lon,dim=2), &
982  size(prism_grid(n)%lat,dim=1),size(prism_grid(n)%lat,dim=2)
983  call oasis_abort()
984  endif
985 
986  !-------------------------------------
987  !> * Gather longitudes if needed and write from root
988  !-------------------------------------
989 
990  filename = 'grids.nc'
991  fldname = trim(cgrid)//'.lon'
992  if (partid_grid) then
993  call oasis_grid_loc2glo(prism_grid(n)%lon,rglo,partid,0)
994  else
995  rglo = prism_grid(n)%lon
996  endif
997  if (write_task) call oasis_io_write_2dgridfld_fromroot(filename,fldname,rglo,nx,ny)
998 
999  !-------------------------------------
1000  !> * Gather latitudes if needed and write from root
1001  !-------------------------------------
1002 
1003  filename = 'grids.nc'
1004  fldname = trim(cgrid)//'.lat'
1005  if (partid_grid) then
1006  call oasis_grid_loc2glo(prism_grid(n)%lat,rglo,partid,0)
1007  else
1008  rglo = prism_grid(n)%lat
1009  endif
1010  if (write_task) call oasis_io_write_2dgridfld_fromroot(filename,fldname,rglo,nx,ny)
1011  endif ! grid_set
1012 
1013  if (prism_grid(n)%corner_set) then
1014  if (tnx <= 0 .or. tny <= 0) then
1015  tnx = size(prism_grid(n)%clon,dim=1)
1016  tny = size(prism_grid(n)%clon,dim=2)
1017  endif
1018  if (size(prism_grid(n)%clon,dim=1) /= tnx .or. &
1019  size(prism_grid(n)%clon,dim=2) /= tny .or. &
1020  size(prism_grid(n)%clat,dim=1) /= tnx .or. &
1021  size(prism_grid(n)%clat,dim=2) /= tny ) then
1022  write(nulprt,*) subname,estr,'inconsistent array size clon/clat ',tnx,tny, &
1023  size(prism_grid(n)%clon,dim=1),size(prism_grid(n)%clon,dim=2), &
1024  size(prism_grid(n)%clat,dim=1),size(prism_grid(n)%clat,dim=2)
1025  call oasis_abort()
1026  endif
1027 
1028  !-------------------------------------
1029  !> * Gather corner longitudes if needed and write from root
1030  !-------------------------------------
1031 
1032  allocate(r3glo(nx,ny,nc))
1033  filename = 'grids.nc'
1034  fldname = trim(cgrid)//'.clo'
1035  if (partid_grid) then
1036  allocate(rloc(tnx,tny))
1037  do n1 = 1,nc
1038  rloc(:,:) = prism_grid(n)%clon(:,:,n1)
1039  call oasis_grid_loc2glo(rloc,rglo,partid,0)
1040  r3glo(:,:,n1) = rglo(:,:)
1041  enddo
1042  deallocate(rloc)
1043  else
1044  r3glo = prism_grid(n)%clon
1045  endif
1046  if (write_task) call oasis_io_write_3dgridfld_fromroot(filename,fldname,r3glo,nx,ny,nc)
1047 
1048  !-------------------------------------
1049  !> * Gather corner latitudes if needed and write from root
1050  !-------------------------------------
1051 
1052  filename = 'grids.nc'
1053  fldname = trim(cgrid)//'.cla'
1054  if (partid_grid) then
1055  allocate(rloc(tnx,tny))
1056  do n1 = 1,nc
1057  rloc(:,:) = prism_grid(n)%clat(:,:,n1)
1058  call oasis_grid_loc2glo(rloc,rglo,partid,0)
1059  r3glo(:,:,n1) = rglo(:,:)
1060  enddo
1061  deallocate(rloc)
1062  else
1063  r3glo = prism_grid(n)%clat
1064  endif
1065  if (write_task) call oasis_io_write_3dgridfld_fromroot(filename,fldname,r3glo,nx,ny,nc)
1066  deallocate(r3glo)
1067  endif ! corner_set
1068 
1069  if (prism_grid(n)%area_set) then
1070  if (tnx <= 0 .or. tny <= 0) then
1071  tnx = size(prism_grid(n)%area,dim=1)
1072  tny = size(prism_grid(n)%area,dim=2)
1073  endif
1074  if (size(prism_grid(n)%area,dim=1) /= tnx .or. &
1075  size(prism_grid(n)%area,dim=2) /= tny ) then
1076  write(nulprt,*) subname,estr,'inconsistent array size area ',tnx,tny, &
1077  size(prism_grid(n)%area,dim=1),size(prism_grid(n)%area,dim=2)
1078  call oasis_abort()
1079  endif
1080 
1081  !-------------------------------------
1082  !> * Gather areas if needed and write from root
1083  !-------------------------------------
1084 
1085  filename = 'areas.nc'
1086  fldname = trim(cgrid)//'.srf'
1087  if (partid_grid) then
1088  call oasis_grid_loc2glo(prism_grid(n)%area,rglo,partid,0)
1089  else
1090  rglo = prism_grid(n)%area
1091  endif
1092  if (write_task) call oasis_io_write_2dgridfld_fromroot(filename,fldname,rglo,nx,ny)
1093  endif ! area_set
1094 
1095  if (prism_grid(n)%angle_set) then
1096  if (tnx <= 0 .or. tny <= 0) then
1097  tnx = size(prism_grid(n)%angle,dim=1)
1098  tny = size(prism_grid(n)%angle,dim=2)
1099  endif
1100  if (size(prism_grid(n)%angle,dim=1) /= tnx .or. &
1101  size(prism_grid(n)%angle,dim=2) /= tny ) then
1102  write(nulprt,*) subname,estr,'inconsistent array size angle ',tnx,tny, &
1103  size(prism_grid(n)%angle,dim=1),size(prism_grid(n)%angle,dim=2)
1104  call oasis_abort()
1105  endif
1106 
1107  !-------------------------------------
1108  !> * Gather angles if needed and write from root
1109  !-------------------------------------
1110 
1111  filename = 'grids.nc'
1112  fldname = trim(cgrid)//'.ang'
1113  if (partid_grid) then
1114  call oasis_grid_loc2glo(prism_grid(n)%lon,rglo,partid,0)
1115  else
1116  rglo = prism_grid(n)%angle
1117  endif
1118  if (write_task) call oasis_io_write_2dgridfld_fromroot(filename,fldname,rglo,nx,ny)
1119  endif ! angle_set
1120 
1121  if (prism_grid(n)%mask_set) then
1122  if (tnx <= 0 .or. tny <= 0) then
1123  tnx = size(prism_grid(n)%mask,dim=1)
1124  tny = size(prism_grid(n)%mask,dim=2)
1125  endif
1126  if (size(prism_grid(n)%mask,dim=1) /= tnx .or. &
1127  size(prism_grid(n)%mask,dim=2) /= tny ) then
1128  write(nulprt,*) subname,estr,'inconsistent array size mask ',tnx,tny, &
1129  size(prism_grid(n)%mask,dim=1),size(prism_grid(n)%mask,dim=2)
1130  call oasis_abort()
1131  endif
1132 
1133  !-------------------------------------
1134  !> * Gather masks if needed and write from root
1135  !-------------------------------------
1136 
1137  allocate(iglo(nx,ny))
1138  filename = 'masks.nc'
1139  fldname = trim(cgrid)//'.msk'
1140  if (partid_grid) then
1141  allocate(rloc(tnx,tny))
1142  rloc = prism_grid(n)%mask
1143  call oasis_grid_loc2glo(rloc,rglo,partid,0)
1144  iglo = nint(rglo)
1145  deallocate(rloc)
1146  else
1147  iglo = prism_grid(n)%mask
1148  endif
1149  if (write_task) call oasis_io_write_2dgridint_fromroot(filename,fldname,iglo,nx,ny)
1150  deallocate(iglo)
1151  endif ! mask_set
1152 
1153  deallocate(rglo)
1154 
1155  endif ! active_task
1156 
1157  endif ! grid_name = gname
1158  endif ! terminated
1159  enddo ! n = 1,prism_ngrid
1160  enddo ! g = 1,gcnt
1161 
1162  deallocate(gname,pname)
1163 
1164  if (local_timers_on) call oasis_timer_stop('grid_write_writefiles')
1165  endif ! writing_grids_call
1166  call oasis_timer_stop('grid_write')
1167  call oasis_debug_exit(subname)
1168 
1169  END SUBROUTINE oasis_write2files
1170 !--------------------------------------------------------------------------
1171 
1172 !> Local interface to find gridID for a specified grid name
1173 
1174  SUBROUTINE oasis_findgrid(cgrid,nx,ny,gridID)
1175  !-------------------------------------------------
1176  ! Routine that sets gridID, identifies existing
1177  ! grid with cgrid name or starts a new one
1178  !-------------------------------------------------
1179  implicit none
1180 
1181  character(len=*), intent (in) :: cgrid !< grid name
1182  integer(kind=ip_intwp_p), intent (in) :: nx !< global nx size
1183  integer(kind=ip_intwp_p), intent (in) :: ny !< global ny size
1184  integer(kind=ip_intwp_p), intent(out) :: gridID !< gridID matching cgrid
1185  !-------------------------------------------------
1186  integer(kind=ip_intwp_p) :: n
1187  character(len=*),parameter :: subname = '(oasis_findgrid)'
1188  !-------------------------------------------------
1189 
1190  call oasis_debug_enter(subname)
1191 
1192  gridid = -1
1193  do n = 1,prism_ngrid
1194  if (trim(cgrid) == trim(prism_grid(n)%gridname)) then
1195  gridid = n
1196  ! since grid is defined before, make sure nx,ny match
1197  if (nx /= prism_grid(gridid)%nx .or. ny /= prism_grid(gridid)%ny) then
1198  write(nulprt,*) subname,estr,'in predefined grid size = ',nx,ny, &
1199  prism_grid(gridid)%nx,prism_grid(gridid)%ny
1200  call oasis_abort()
1201  endif
1202  endif
1203  enddo
1204 
1205  if (gridid < 1) then
1206  prism_ngrid = prism_ngrid+1
1207  gridid = prism_ngrid
1208  endif
1209 
1210  prism_grid(gridid)%gridname = trim(cgrid)
1211  prism_grid(gridid)%nx = nx
1212  prism_grid(gridid)%ny = ny
1213 
1214  call oasis_debug_exit(subname)
1215 
1216  END SUBROUTINE oasis_findgrid
1217 
1218 !--------------------------------------------------------------------------
1219 
1220 !> Local routine that gathers the local array using partition information
1221 
1222  SUBROUTINE oasis_grid_loc2glo(aloc,aglo,partid,taskid)
1223  implicit none
1224 
1225  real(kind=ip_realwp_p),intent(in) :: aloc(:,:) !< local array
1226  real(kind=ip_realwp_p),intent(inout) :: aglo(:,:) !< global array
1227  integer(kind=ip_i4_p) ,intent(in) :: partid !< partition id for local data
1228  integer(kind=ip_i4_p) ,intent(in) :: taskid !< task id to gather data to
1229  !-------------------------------------------------
1230  type(mct_avect) :: avloc,avglo
1231  integer(kind=ip_i4_p) :: i,j,n
1232  integer(kind=ip_i4_p) :: lnx,lny,gnx,gny
1233  character(len=*),parameter :: subname = '(oasis_grid_loc2glo)'
1234  !-------------------------------------------------
1235 
1236  call oasis_debug_enter(subname)
1237 
1238  if (prism_part(partid)%mpicom /= mpi_comm_null) then
1239 
1240  lnx = size(aloc,dim=1)
1241  lny = size(aloc,dim=2)
1242  gnx = size(aglo,dim=1)
1243  gny = size(aglo,dim=2)
1244  call mct_avect_init(avloc,rlist='field',lsize=lnx*lny)
1245 
1246  n = 0
1247  do j = 1,lny
1248  do i = 1,lnx
1249  n = n + 1
1250  avloc%rattr(1,n) = aloc(i,j)
1251  enddo
1252  enddo
1253 
1254  call mct_avect_gather(avloc,avglo,prism_part(partid)%pgsmap,taskid,prism_part(partid)%mpicom)
1255 
1256  if (prism_part(partid)%rank == taskid) then
1257  n = 0
1258  do j = 1,gny
1259  do i = 1,gnx
1260  n = n + 1
1261  aglo(i,j) = avglo%rattr(1,n)
1262  enddo
1263  enddo
1264  call mct_avect_clean(avglo)
1265  endif
1266 
1267  call mct_avect_clean(avloc)
1268 
1269  endif
1270 
1271  call oasis_debug_exit(subname)
1272 
1273  END SUBROUTINE oasis_grid_loc2glo
1274 
1275 !--------------------------------------------------------------------------
1276 END MODULE mod_oasis_grid
1277 !--------------------------------------------------------------------------
1278 
1279 
1280 
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
OASIS partition data and methods.
System type methods.
Generic interface to support writing 4 or 8 byte reals.
subroutine oasis_findgrid(cgrid, nx, ny, gridID)
Local interface to find gridID for a specified grid name.
Generic overloaded interface into MPI max reduction.
Provides a generic and simpler interface into MPI calls for OASIS.
Generic overloaded interface into MPI broadcast.
Generic interface to support writing 4 or 8 byte reals.
subroutine, public oasis_timer_start(timer_label, barrier)
Start a timer.
subroutine, public oasis_timer_stop(timer_label)
Stop a timer.
Provides a common location for several OASIS variables.
Generic interface to support writing 4 or 8 byte reals.
subroutine, public oasis_mpi_reducelists(linp1, comm, cntout, lout1, callstr, fastcheck, fastcheckout, linp2, lout2, spval2, linp3, lout3, spval3)
Custom method for reducing MPI lists across pes for OASIS.
subroutine, public oasis_io_write_3dgridfld_fromroot(filename, fldname, fld, nx, ny, nc)
Write a 3d real array named field from the root task to a file.
subroutine, public oasis_debug_enter(string)
Used when a subroutine is entered, write info to log file at some debug level.
Generic interface to support writing 4 or 8 byte reals.
Generic overloaded interface into MPI min reduction.
subroutine oasis_grid_loc2glo(aloc, aglo, partid, taskid)
Local routine that gathers the local array using partition information.
subroutine, public oasis_terminate_grids_writing()
User interface to indicate user defined grids are done.
subroutine oasis_write_angle_r8(cgrid, nx, ny, angle, partid)
User interface to set angle for 8 byte reals.
subroutine, public oasis_mpi_chkerr(rcode, string)
Checks MPI error codes and aborts.
subroutine, public oasis_write2files()
Interface that actually writes fields to grid files.
Performance timer methods.
subroutine oasis_write_corner_r8(cgrid, nx, ny, nc, clon, clat, partid)
User interface to set corner latitudes and longitudes for 8 byte reals.
subroutine, public oasis_write_mask(cgrid, nx, ny, mask, partid)
User interface to set integer mask values.
subroutine, public oasis_abort(id_compid, cd_routine, cd_message)
OASIS abort method, publically available to users.
subroutine oasis_write_grid_r8(cgrid, nx, ny, lon, lat, partid)
User interface to set latitudes and longitudes for 8 byte reals.
subroutine, public oasis_print_grid_data()
Print grid information to log file.
subroutine, public oasis_io_write_2dgridfld_fromroot(filename, fldname, fld, nx, ny)
Write a real array named field from the root task to a file.
subroutine, public oasis_start_grids_writing(iwrite)
User interface to initialize grid writing.
Model grid data for creating mapping data and conserving fields.
subroutine, public oasis_mpi_barrier(comm, string)
Call MPI_BARRIER for a particular communicator.
OASIS grid data and methods.
subroutine oasis_write_area_r8(cgrid, nx, ny, area, partid)
User interface to set area values for 8 byte reals.
subroutine, public oasis_io_write_2dgridint_fromroot(filename, fldname, fld, nx, ny)
Write an integer array named field from the root task to a file.
subroutine oasis_write_area_r4(cgrid, nx, ny, area, partid)
User interface to set area values for 4 byte reals.
Provides reusable IO routines for OASIS.
Definition: mod_oasis_io.F90:4
subroutine oasis_write_grid_r4(cgrid, nx, ny, lon, lat, partid)
User interface to set latitudes and longitudes for 4 byte reals.
subroutine oasis_write_corner_r4(cgrid, nx, ny, nc, clon, clat, partid)
User interface to set corner latitudes and longitudes for 4 byte reals.
subroutine oasis_write_angle_r4(cgrid, nx, ny, angle, partid)
User interface to set angle for 4 byte reals.