Line data Source code
1 : #define VERIFY_(A) if(.not.HCO_ESMF_VRFY(A,subname,__LINE__)) stop -1
2 : #define ASSERT_(A) if(.not.HCO_ESMF_ASRT(A,subname,__LINE__)) stop -1
3 : !------------------------------------------------------------------------------
4 : ! Harmonized Emissions Component (HEMCO) !
5 : !------------------------------------------------------------------------------
6 : !BOP
7 : !
8 : ! !MODULE: hco_esmf_grid
9 : !
10 : ! !DESCRIPTION: Module HCO\_ESMF\_GRID defines functions for the regridding and
11 : ! control of the intermediate HEMCO grid in the HEMCO to CAM interface.
12 : !\\
13 : !\\
14 : ! !INTERFACE:
15 : !
16 : module hco_esmf_grid
17 : !
18 : ! !USES:
19 : !
20 : ! ESMF function wrappers
21 : use hco_esmf_wrappers
22 :
23 : ! ESMF types
24 : use ESMF, only: ESMF_Mesh, ESMF_DistGrid, ESMF_Grid
25 : use ESMF, only: ESMF_Field
26 : use ESMF, only: ESMF_RouteHandle
27 : use ESMF, only: ESMF_SUCCESS
28 :
29 : use ESMF, only: ESMF_LogWrite, ESMF_LOGMSG_INFO
30 :
31 : ! MPI
32 : use mpi, only: MPI_PROC_NULL, MPI_SUCCESS, MPI_INTEGER
33 :
34 : ! Floating point type.
35 : ! FIXME: May change to HEMCO precision later down the line
36 : use shr_kind_mod, only: r8 => shr_kind_r8
37 : use ESMF, only: ESMF_KIND_I4, ESMF_KIND_R8
38 :
39 : implicit none
40 : private
41 : save
42 : !
43 : ! !PUBLIC MEMBER FUNCTIONS:
44 : !
45 : public :: HCO_Grid_Init
46 : public :: HCO_Grid_UpdateRegrid
47 :
48 : public :: HCO_Grid_HCO2CAM_2D ! Regrid HEMCO to CAM mesh on 2D field
49 : public :: HCO_Grid_HCO2CAM_3D ! ...on 3D field
50 : public :: HCO_Grid_CAM2HCO_2D ! Regrid CAM mesh to HEMCO on 2D field
51 : public :: HCO_Grid_CAM2HCO_3D ! ...on 3D field
52 :
53 : !
54 : ! !PRIVATE MEMBER FUNCTIONS:
55 : !
56 : private :: HCO_Grid_ESMF_CreateCAM
57 : private :: HCO_Grid_ESMF_CreateCAMField ! Create field on CAM physics mesh
58 : private :: HCO_Grid_ESMF_CreateHCO ! Create HEMCO lat-lon grid in ESMF
59 : private :: HCO_Grid_ESMF_CreateHCOField ! Create field on HEMCO ll grid
60 :
61 : private :: HCO_ESMF_Set2DHCO ! Set ESMF field with 2D HEMCO data
62 : private :: HCO_ESMF_Set3DHCO ! Set ESMF field with 3D HEMCO data
63 : private :: HCO_ESMF_Set2DCAM ! Set ESMF field with 2D CAM mesh data
64 : private :: HCO_ESMF_Set3DCAM ! Set ESMF field with 3D CAM mesh data
65 :
66 : private :: HCO_ESMF_Get1DField ! Retrieve 1D ESMF field pointer
67 : private :: HCO_ESMF_Get2DField ! Retrieve 2D ESMF field pointer
68 : private :: HCO_ESMF_Get3DField ! Retrieve 3D ESMF field pointer
69 :
70 : ! !REMARKS:
71 : !
72 : ! Horizontal grid specifications are based on GEOS-Chem definitions.
73 : ! Vertical grid specifications are copied from hyai, hybi from CAM.
74 : !
75 : ! Notes:
76 : ! (i) In GEOS-Chem, level 1 is bottom-of-atmos, so 'bottom2top' lev_sequence.
77 : ! (ii)
78 : !
79 : ! The CreateHCO routine replaces the edyn_esmf routines for create_geo and geo2phys
80 : ! grid as they only differ by a grid staggering configuration (center and corner).
81 : ! Thus they are merged into one routine and write to two grids, HCO_Grid and HCO2CAM_Grid.
82 : ! (hplin, 2/20/20)
83 : !
84 : ! The HEMCO grid is initialized in ESMF with 1 periodic dim on the longitude.
85 : ! It remains to be seen whether this is necessary and validation of the regridder's
86 : ! robustness is necessary (hplin, 2/23/20)
87 : !
88 : ! Only 3-D regridding implemented for now. (hplin, 2/24/20)
89 : !
90 : ! ---------------------- BELOW ARE SIDE REMARKS ---------------------------------
91 : ! This module was written by hplin on a gloomy day (as usual) in Cambridge/Boston
92 : ! where I patiently copied the following modules and translated the terminology
93 : ! between them:
94 : ! GEOS-Chem: state_grid_mod, hco_types_mod, state_met_mod, pressure_mod
95 : ! CAM Ionosphere Interface: edyn_mpi, edyn_geogrid, ionosphere_interface, hycoef, ref_pres
96 : ! The result is a mash-up of CAM and GEOS-Chem coding terminology, but in the end
97 : ! it is mostly GEOS-Chem. Conversions are done at initialization.
98 : !
99 : ! The vertical grid specification conversions were hand-written
100 : ! with napkin-based calculations, so buyer beware. They need to be carefully
101 : ! validated before final shipping.
102 : !
103 : ! Horizontal grid specifications were mostly borrowed from GEOS-Chem definitions.
104 : ! We use a global ??? grid for testing at this point and they need to be
105 : ! un-hardcoded in the future.
106 : !
107 : ! Thanks to Thibaud Fritz for supplying CESM-GC dev code which helped in the
108 : ! conversions from CAM to GC speak.
109 : !
110 : ! -- hplin, 2/11/20
111 : !
112 : ! !REVISION HISTORY:
113 : ! 11 Feb 2020 - H.P. Lin - First crack.
114 : ! 17 Feb 2020 - H.P. Lin - Start of work on regrid routines.
115 : ! 15 Apr 2020 - H.P. Lin - Regrid routines mostly finalized
116 : !EOP
117 : !------------------------------------------------------------------------------
118 : !BOC
119 : !
120 : ! !PUBLIC TYPES:
121 : !
122 : ! Global grid parameters.
123 : ! this may be refactored into some other structure that isn't global later.
124 : ! for now this will do (hplin, 2/11/20) -- and I am confident this will
125 : ! be the way for at least a few more years, because who touches working code? ;)
126 : integer, public, protected :: IM ! # of lons
127 : integer, public, protected :: JM ! # of lats
128 : integer, public, protected :: LM ! # of levs
129 :
130 : ! Computed parameters for compatibility with GEOS-Chem
131 : real(r8), public, protected:: DX ! Delta X [deg long]
132 : real(r8), public, protected:: DY ! Delta X [deg lat]
133 :
134 : ! Horizontal Coordinates
135 : real(r8), public, pointer :: &
136 : XMid (:,:), & ! Longitude centers [deg]
137 : XEdge(:,:), & ! Longitude edges [deg]
138 : YMid (:,:), & ! Latitude centers [deg]
139 : YEdge(:,:), & ! Latitude edges [deg]
140 : YEdge_R(:,:), & ! Latitude edges R [rad]
141 : YSin (:,:) ! SIN( lat edges ) [1]
142 :
143 : ! Shadow variables of geo-"meteorological fields" required by HEMCO
144 : !
145 : ! Hybrid Grid Coordinate Definition: (dsa, bmy, 8/27/02, 2/2/12)
146 : ! ============================================================================
147 : !
148 : ! The pressure at the bottom edge of grid box (I,J,L) is defined as follows:
149 : ! Pedge(I,J,L) = Ap(L) + [ Bp(L) * Psurface(I,J) ]
150 : ! where
151 : ! Psurface(I,J) is the "true" surface pressure at lon,lat (I,J)
152 : ! Ap(L) has the same units as surface pressure [hPa]
153 : ! Bp(L) is a unitless constant given at level edges
154 : !
155 : ! Note: I don't know what to do about PEDGE, surface pressure and the like.
156 : ! they likely need a regrid through ESMF from State%PSDry or something.
157 : ! This is only available in the GridComp Run, not worrying about it here.
158 : ! (hplin, 2/11/20)
159 : real(r8), public, pointer :: &
160 : AREA_M2(:,:), & ! Area of grid box [m^2]
161 : Ap (:), & ! "hyai" Hybrid-sigma Ap value [Pa]
162 : Bp (:) ! "hybi" Hybrid-sigma Bp value [Pa]
163 :
164 : ! MPI Descriptors.
165 : ! Ported mostly from edyn_geogrid and edyn_mpi
166 : ! -- What everyone knows --
167 : integer, public, protected :: HCO_mpicom
168 : integer, public, protected :: nPET ! Number of PETs
169 : integer, public, protected :: nPET_lon, nPET_lat ! # of PETs over lon, lat
170 :
171 : integer, public, allocatable :: HCO_petTable(:,:)! 2D table of tasks (dim'l nPET_lon+2, ..lat+2)
172 : ! extra left and right used for halos
173 : integer, public, allocatable :: HCO_petMap(:,:,:)! PETmap for ESMF
174 :
175 : ! -- Private to MPI process --
176 : ! Note L dimension (levs) not distributed
177 : integer, public, protected :: my_IM, my_JM ! # of lons, levs in this task
178 : integer, public, protected :: my_IS, my_IE ! First and last lons
179 : integer, public, protected :: my_JS, my_JE ! First and last lats
180 :
181 : integer, public, protected :: my_ID ! my task ID in HCO_Task
182 : integer, public, protected :: my_ID_I, my_ID_J ! mytidi, mytidj coord for current task
183 :
184 : integer, public, protected :: my_CE ! # of CAM ncols in this task
185 :
186 : type HCO_Task
187 : integer :: ID ! identifier
188 :
189 : integer :: ID_I ! task coord in longitude dim'l of task table
190 : integer :: ID_J ! task coord in latitude dim'l of task table
191 : integer :: IM ! # of lons on this task
192 : integer :: JM ! # of lats on this task
193 : integer :: IS, IE ! start and end longitude dim'l index
194 : integer :: JS, JE ! start and end latitude dim'l index
195 : end type HCO_Task
196 : type(HCO_Task), allocatable:: HCO_Tasks(:) ! HCO_Tasks(nPET) avail to all tasks
197 : !
198 : ! !PRIVATE TYPES:
199 : !
200 : ! ESMF grid and meshes for regridding
201 : type(ESMF_Grid) :: HCO_Grid
202 : type(ESMF_Mesh) :: CAM_PhysMesh ! Copy of CAM physics mesh decomposition
203 : type(ESMF_DistGrid) :: CAM_DistGrid ! DE-local allocation descriptor DistGrid (2D)
204 :
205 : ! ESMF fields for mapping between HEMCO to CAM fields
206 : type(ESMF_Field) :: CAM_2DFld, CAM_3DFld
207 : type(ESMF_Field) :: HCO_2DFld, HCO_3DFld
208 :
209 : ! Used to generate regridding weights
210 : integer :: cam_last_atm_id ! Last CAM atmospheric ID
211 :
212 : ! Regridding weight route handles
213 : type(ESMF_RouteHandle) :: HCO2CAM_RouteHandle_3D, &
214 : HCO2CAM_RouteHandle_2D, &
215 : CAM2HCO_RouteHandle_3D, &
216 : CAM2HCO_RouteHandle_2D
217 : contains
218 : !EOC
219 : !------------------------------------------------------------------------------
220 : ! Harmonized Emissions Component (HEMCO) !
221 : !------------------------------------------------------------------------------
222 : !BOP
223 : !
224 : ! !IROUTINE: HCO_Grid_Init
225 : !
226 : ! !DESCRIPTION: Subroutine HCO\_Grid\_Init initializes the HEMCO-CAM interface
227 : ! grid descriptions and MPI distribution.
228 : !\\
229 : !\\
230 : ! !INTERFACE:
231 : !
232 0 : subroutine HCO_Grid_Init( IM_in, JM_in, nPET_in, RC )
233 : !
234 : ! !USES:
235 : !
236 : ! MPI Properties from CAM
237 : ! Even though CAM's principle is that only spmd_utils uses MPI,
238 : ! ionos code uses MPI very liberally. Unfortunately we have to follow
239 : ! this example as spmd_utils does not provide many of the relevant
240 : ! functionality like the communicator split. But keep this in mind.
241 : use cam_logfile, only: iulog
242 : use spmd_utils, only: CAM_mpicom => mpicom, CAM_npes => npes
243 : use spmd_utils, only: MPI_SUCCESS
244 : use spmd_utils, only: iam
245 : use spmd_utils, only: masterproc
246 :
247 : ! Physical constants
248 : use shr_const_mod, only: pi => shr_const_pi
249 : use shr_const_mod, only: Re => shr_const_rearth
250 :
251 : ! Grid specifications and information from CAM
252 : use hycoef, only: ps0, hyai, hybi ! Vertical specs
253 : use ppgrid, only: pver ! # of levs
254 :
255 : !
256 : ! !INPUT PARAMETERS:
257 : !
258 : integer, intent(in) :: IM_in, JM_in ! # lon, lat, lev global
259 : integer, intent(in) :: nPET_in ! # of PETs to distribute to?
260 : integer, intent(inout) :: RC ! Return code
261 : !
262 : ! !REMARKS:
263 : !
264 : ! !REVISION HISTORY:
265 : ! 11 Feb 2020 - H.P. Lin - Initial version
266 : !EOP
267 : !------------------------------------------------------------------------------
268 : !BOC
269 : !
270 : ! !LOCAL VARIABLES:
271 : !
272 : character(len=*), parameter :: subname = 'HCO_Grid_Init'
273 : integer :: I, J, L, N
274 : real(r8) :: SIN_N, SIN_S, PI_180
275 :
276 : ! MPI stuff
277 : integer :: color
278 : integer :: lons_per_task, lons_overflow
279 : integer :: lats_per_task, lats_overflow
280 : integer :: lon_beg, lon_end, lat_beg, lat_end, task_cnt
281 :
282 : ! Send and receive buffers
283 0 : integer, allocatable :: itasks_send(:,:), itasks_recv(:,:)
284 :
285 : ! Reset CAM atmospheric ID because we know nothing about it (hplin, 2/20/20)
286 0 : cam_last_atm_id = -999
287 :
288 : ! Some physical constants...
289 0 : PI_180 = pi / 180.0_r8
290 :
291 : ! Accept external dimensions.
292 0 : IM = IM_in
293 0 : JM = JM_in
294 0 : nPET = nPET_in
295 :
296 : !-----------------------------------------------------------------------
297 : ! Compute vertical grid parameters
298 : !-----------------------------------------------------------------------
299 : ! Can be directly retrieved from hyai, hybi
300 : ! Although they need to be flipped to be passed from CAM (from tfritz)
301 : !
302 : ! Note: In GEOS-Chem, Ap, Bp are defined in hPa, 1
303 : ! but in HEMCO, they are defined as Pa, 1 (see hcoi_gc_main_mod.F90 :2813)
304 : ! So you have to be especially wary of the units.
305 :
306 : ! For now, use the CAM vertical grid verbatim
307 0 : LM = pver
308 :
309 : ! Ap, Bp has LM+1 edges for LM levels
310 0 : allocate(Ap(LM + 1), STAT=RC) ! LM levels, LM+1 edges
311 0 : allocate(Bp(LM + 1), STAT=RC)
312 0 : ASSERT_(RC==0)
313 :
314 : ! Allocate PEDGE information
315 :
316 : ! G-C def: Pedge(I,J,L) = Ap(L) + [ Bp(L) * Psurface(I,J) ]
317 : ! CAM def: Pifce( L) = hyai(k)*ps0 + [ hybi(k) * ps ]
318 : ! w.r.t. ps0 = base state srfc prs; ps = ref srfc prs.
319 : !
320 : ! Note that the vertical has to be flipped and this will need to be done
321 : ! everywhere else within HEMCO_CESM, too.
322 0 : do L = 1, (LM+1)
323 0 : Ap(L) = hyai(LM+2-L) * ps0
324 0 : Bp(L) = hybi(LM+2-L)
325 : enddo
326 :
327 : !-----------------------------------------------------------------------
328 : ! Compute horizontal grid parameters
329 : !-----------------------------------------------------------------------
330 : ! Notes: long range (i) goes from -180.0_r8 to +180.0_r8
331 : ! lat range (j) goes from - 90.0_r8 to + 90.0_r8
332 :
333 0 : allocate(XMid (IM, JM ), STAT=RC)
334 0 : allocate(XEdge(IM+1, JM ), STAT=RC)
335 0 : allocate(YMid (IM, JM ), STAT=RC)
336 0 : allocate(YEdge(IM, JM+1), STAT=RC)
337 0 : allocate(YEdge_R(IM, JM+1), STAT=RC)
338 0 : allocate(YSin (IM, JM+1), STAT=RC)
339 0 : allocate(AREA_M2(IM, JM ), STAT=RC)
340 0 : ASSERT_(RC==0)
341 :
342 : ! Compute DX, DY (lon, lat)
343 0 : DX = 360.0_r8 / real(IM, r8)
344 0 : DY = 180.0_r8 / real((JM - 1), r8)
345 :
346 : ! Loop over horizontal grid
347 : ! Note: Might require special handling at poles. FIXME. (hplin, 2/11/20)
348 0 : do J = 1, JM
349 0 : do I = 1, IM
350 : ! Longitude centers [deg]
351 0 : XMid(I, J) = (DX * (I-1)) - 180.0_r8
352 :
353 : ! Latitude centers [deg]
354 0 : YMid(I, J) = (DY * (J-1)) - 90.0_r8
355 :
356 : ! Note half-sized polar boxes for global grid, multiply DY by 1/4 at poles
357 0 : if(J == 1) then
358 0 : YMid(I, 1) = -90.0_r8 + (0.25_r8 * DY)
359 : endif
360 0 : if(J == JM) then
361 0 : YMid(I, JM) = 90.0_r8 - (0.25_r8 * DY)
362 : endif
363 :
364 : ! Edges [deg] (or called corners in CAM ionos speak)
365 0 : XEdge(I, J) = XMid(I, J) - DX * 0.5_r8
366 0 : YEdge(I, J) = YMid(I, J) - DY * 0.5_r8
367 0 : YEdge_R(I, J) = (PI_180 * YEdge(I, J))
368 0 : YSin (I, J) = SIN( YEdge_R(I, J) ) ! Needed for MAP_A2A regridding
369 :
370 : ! Compute the LAST edges
371 0 : if(I == IM) then
372 0 : XEdge(I+1,J) = XEdge(I, J) + DX
373 : endif
374 :
375 : ! Enforce half-sized polar boxes where northern edge of grid boxes
376 : ! along the SOUTH POLE to be -90 deg lat.
377 0 : if(J == 1) then
378 0 : YEdge(I, 1) = -90.0_r8
379 : endif
380 :
381 0 : if(J == JM) then
382 : ! Northern edge of grid boxes along the north pole to be +90 deg lat
383 0 : YEdge(I,J+1) = 90.0_r8
384 :
385 : ! Adjust for second-to-last lat edge
386 0 : YEdge(I,J ) = YEdge(I,J+1) - (DY * 0.5_r8)
387 0 : YEdge_R(I,J) = YEdge(I,J) * PI_180
388 0 : YSin(I, J) = SIN( YEdge_R(I, J) )
389 :
390 : ! Last latitude edge [radians]
391 0 : YEdge_R(I,J+1) = YEdge(I,J+1) * PI_180
392 0 : YSin(I,J+1) = SIN( YEdge_R(I,J+1) )
393 : endif
394 : enddo
395 : enddo
396 :
397 : ! Compute grid box areas after everything is populated...
398 0 : do J = 1, JM
399 0 : do I = 1, IM
400 : ! Sine of latitudes at N and S edges of grid box (I,J)
401 0 : SIN_N = SIN( YEdge_R(I,J+1) )
402 0 : SIN_S = SIN( YEdge_R(I,J ) )
403 :
404 : ! Grid box surface areas [m2]
405 0 : AREA_M2(I,J) = ( DX * PI_180 ) * ( Re**2 ) * (SIN_N - SIN_S)
406 : enddo
407 : enddo
408 :
409 : ! Output debug information on the global grid information
410 : ! Copied from gc_grid_mod.F90 and pressure_mod.F
411 0 : if(masterproc) then
412 0 : write( iulog, '(a)' )
413 0 : write( iulog, '(''%%%%%%%%%%%%%%% HEMCO GRID %%%%%%%%%%%%%%%'')' )
414 0 : write( iulog, '(a)' )
415 0 : write( iulog, *) 'DX', DX, 'DY', DY
416 0 : write( iulog, '(''Grid box longitude centers [degrees]: '')' )
417 0 : write( iulog, * ) size(XMid, 1), size(XMid, 2)
418 0 : write( iulog, '(8(f8.3,1x))' ) ( XMid(I,1), I=1,IM )
419 0 : write( iulog, '(a)' )
420 0 : write( iulog, '(''Grid box longitude edges [degrees]: '')' )
421 0 : write( iulog, * ) size(XEdge, 1), size(XEdge, 2)
422 0 : write( iulog, '(8(f8.3,1x))' ) ( XEdge(I,1), I=1,IM+1 )
423 0 : write( iulog, '(a)' )
424 0 : write( iulog, '(''Grid box latitude centers [degrees]: '')' )
425 0 : write( iulog, * ) size(YMid, 1), size(YMid, 2)
426 0 : write( iulog, '(8(f8.3,1x))' ) ( YMid(1,J), J=1,JM )
427 0 : write( iulog, '(a)' )
428 0 : write( iulog, '(''Grid box latitude edges [degrees]: '')' )
429 0 : write( iulog, * ) size(YEdge, 1), size(YEdge, 2)
430 0 : write( iulog, '(8(f8.3,1x))' ) ( YEdge(1,J), J=1,JM+1 )
431 0 : write( iulog, '(a)' )
432 0 : write( iulog, '(''SIN( grid box latitude edges )'')' )
433 0 : write( iulog, '(8(f8.3,1x))' ) ( YSin(1,J), J=1,JM+1 )
434 :
435 0 : write( iulog, '(a)' ) REPEAT( '=', 79 )
436 0 : write( iulog, '(a,/)' ) 'V E R T I C A L G R I D S E T U P'
437 0 : write( iulog, '( ''Ap '', /, 6(f11.6,1x) )' ) AP(1:LM+1)
438 0 : write( iulog, '(a)' )
439 0 : write( iulog, '( ''Bp '', /, 6(f11.6,1x) )' ) BP(1:LM+1)
440 0 : write( iulog, '(a)' ) REPEAT( '=', 79 )
441 : endif
442 :
443 : !-----------------------------------------------------------------------
444 : ! Distribute among parallelization in MPI 1: Compute distribution
445 : !-----------------------------------------------------------------------
446 :
447 : ! edyn_geogrid uses 1-D latitude decomposition, so nPET_lon = 1
448 : ! and nPET_lat = JM. From tfritz this may not work well with GEOS-Chem
449 : ! so we may need to attempt some other decomposition in the future.
450 : !
451 : ! The code below from edyn_geogrid may not be generic enough for that
452 : ! need, so we might do nPET_lon = 1 for now. (See code path below)
453 : ! (hplin, 2/13/20)
454 :
455 : ! It seems like ESMF conservative regridding will not work with DE < 2
456 : ! so the distribution must assign more than 1 lat and lon per PET.
457 : ! The code has been updated accordingly. (hplin, 2/22/20)
458 0 : do nPET_lon = 2, IM
459 0 : nPET_lat = nPET / nPET_lon
460 0 : if( nPET_lon * nPET_lat == nPET .and. nPET_lon .le. IM .and. nPET_lat .le. JM ) then
461 : ! Enforce more than 2-width lon and lat...
462 0 : if(int(IM / nPET_lon) > 1 .and. int(JM / nPET_lat) > 1) then
463 : exit
464 : endif
465 : endif
466 : enddo
467 : ! nPET_lon = 1
468 : ! nPET_lat = nPET
469 :
470 : ! Verify we have a correct decomposition
471 : ! Can't accept invalid decomp; also cannot accept IM, 1 (for sake of consistency)
472 0 : if( nPET_lon * nPET_lat /= nPET .or. nPET_lat == 1 ) then
473 : ! Fall back to same 1-D latitude decomposition like edyn_geogrid
474 0 : nPET_lon = 1
475 0 : nPET_lat = nPET
476 :
477 0 : if(masterproc) then
478 0 : write(iulog,*) "HEMCO: HCO_Grid_Init failed to find a secondary decomp."
479 0 : write(iulog,*) "Using 1-d latitude decomp. This may fail with ESMF regrid."
480 : endif
481 : endif
482 :
483 : ! Verify for the 1-D latitude decomposition case (edge case)
484 : ! if the number of CPUs is reasonably set. If nPET_lon = 1, then you cannot
485 : ! allow nPET_lat exceed JM or it will blow up.
486 0 : if(nPET_lon == 1 .and. nPET_lat == nPET) then
487 0 : if(nPET_lat > JM) then
488 0 : if(masterproc) then
489 0 : write(iulog,*) "HEMCO: Warning: Input nPET > JM", nPET, JM
490 0 : write(iulog,*) "HEMCO: I will use nPET = JM for now."
491 : endif
492 : endif
493 : endif
494 :
495 : ! Commit to the decomposition at this point
496 0 : if(masterproc) then
497 0 : write(iulog,*) "HEMCO: HCO_Grid_Init IM, JM, LM", IM, JM, LM
498 0 : write(iulog,*) "HEMCO: nPET_lon * nPET_lat = ", nPET_lon, nPET_lat, nPET
499 : endif
500 :
501 : ! Figure out beginning and ending coordinates for each task
502 : ! copied from edyn_geogrid
503 0 : lons_per_task = IM / nPET_lon
504 0 : lons_overflow = MOD(IM, nPET_lon)
505 0 : lats_per_task = JM / nPET_lat
506 0 : lats_overflow = MOD(JM, nPET_lat)
507 0 : lon_beg = 1
508 0 : lon_end = 0
509 0 : lat_beg = 1
510 0 : lat_end = 0
511 0 : task_cnt = 0
512 0 : jloop: do J = 0, nPET_lat - 1
513 0 : lat_beg = lat_end + 1
514 0 : lat_end = lat_beg + lats_per_task - 1
515 0 : if (J < lats_overflow) then
516 0 : lat_end = lat_end + 1
517 : endif
518 0 : lon_end = 0
519 0 : do I = 0, nPET_lon - 1
520 0 : lon_beg = lon_end + 1
521 0 : lon_end = lon_beg + lons_per_task - 1
522 0 : if (I < lons_overflow) then
523 0 : lon_end = lon_end + 1
524 : endif
525 0 : task_cnt = task_cnt+1
526 0 : if (task_cnt > iam) exit jloop ! This makes this loop CPU specific
527 : enddo
528 : enddo jloop
529 :
530 : !-----------------------------------------------------------------------
531 : ! Distribute among parallelization in MPI 2: Populate indices and task table
532 : !-----------------------------------------------------------------------
533 :
534 : ! Create communicator
535 : ! Color may be unnecessary if using all CAM processes for CAM_mpicom
536 : ! but we will retain this functionality for now incase needed (hplin, 2/12/20)
537 0 : color = iam / (nPET_lat * nPET_lon)
538 0 : call MPI_comm_split(CAM_mpicom, color, iam, HCO_mpicom, RC)
539 0 : ASSERT_(RC==MPI_SUCCESS)
540 :
541 : ! Distribute among MPI (mp_distribute_geo in edyn_geogrid)
542 : ! Merged all into this huge monolithic routine..
543 : ! (lon_beg, lon_end, lat_beg, lat_end, 1, LM, nPET_lon, nPET_lat)
544 : ! (lonndx0, lonndx1, latndx0, latndx1, levndx0, levndx1, ntaski_in, ntaskj_in)
545 :
546 : ! Get my indices!
547 0 : my_IS = lon_beg
548 0 : my_IE = lon_end
549 0 : my_JS = lat_beg
550 0 : my_JE = lat_end
551 :
552 : ! Allocate task info table
553 0 : allocate(HCO_Tasks(0:nPET-1), stat=RC)
554 :
555 : ! Allocate 2D table of TASKS (not i j coordinates)
556 0 : allocate(HCO_petTable(-1:nPET_lon, -1:nPET_lat), stat=RC)
557 0 : ASSERT_(RC==0)
558 :
559 : ! Allocate PET map for ESMF
560 0 : allocate(HCO_petMap(nPET_lon, nPET_lat, 1))
561 0 : ASSERT_(RC==0)
562 :
563 : ! 2D table of tasks communicates to MPI_PROC_NULL by default so talking
564 : ! to that PID has no effect in MPI comm
565 0 : HCO_petTable(:,:) = MPI_PROC_NULL
566 :
567 : ! Figure out ranks for the petTable, which is a table of I, J PETs
568 : ! with halo (hplin, 2/12/20)
569 0 : my_ID = iam
570 0 : N = 0
571 0 : do J = 0, nPET_lat-1
572 0 : do I = 0, nPET_lon-1
573 0 : HCO_petTable(I, J) = N
574 0 : HCO_petMap(I+1, J+1, 1) = N ! PETmap indices based off 1, so +1 here
575 0 : if(iam == N) then
576 0 : my_ID_I = I
577 0 : my_ID_J = J ! Found my place in the PET table
578 : endif
579 0 : N = N + 1 ! move on to the next rank
580 : enddo
581 :
582 : ! Tasks are periodic in longitude (from edyn_mpi) for haloing
583 : ! FIXME: Check if this is true in HCO distribution. Maybe not
584 0 : HCO_petTable(-1, J) = HCO_petTable(nPET_lon-1, J)
585 0 : HCO_petTable(nPET_lon, J) = HCO_petTable(0, J)
586 : enddo
587 :
588 : ! Print some debug information on the distribution
589 : if( .false. ) then
590 : ! FIXME Don't write to 6 once finished debugging, this literally
591 : ! writes to cesm.log
592 : write(6, "('HEMCO: MPIGrid mytid=',i4,' my_IM, my_JM=',2i4,' my_ID_I,J=',2i4, &
593 : ' lon0,1=',2i4,' lat0,1=',2i4)") &
594 : my_ID,my_IM,my_JM,my_ID_I,my_ID_J,my_IS,my_IE,my_JS,my_JE
595 :
596 : ! write(iulog,"(/,'nPET=',i3,' nPET_lon=',i2,' nPET_lat=',i2,' Task Table:')") &
597 : ! nPET,nPET_lon,nPET_lat
598 : ! do J=-1,nPET_lat
599 : ! write(iulog,"('J=',i3,' HCO_petTable(:,J)=',100i3)") J,HCO_petTable(:,J)
600 : ! enddo
601 : endif
602 :
603 : ! Each task should know its role now...
604 0 : my_IM = my_IE - my_IS + 1
605 0 : my_JM = my_JE - my_JS + 1
606 :
607 : ! Fill all PET info arrays with our information first
608 0 : do N = 0, nPET-1
609 0 : HCO_Tasks(N)%ID = iam ! identifier
610 0 : HCO_Tasks(N)%ID_I = my_ID_I ! task coord in longitude dim'l of task table
611 0 : HCO_Tasks(N)%ID_J = my_ID_J ! task coord in latitude dim'l of task table
612 0 : HCO_Tasks(N)%IM = my_IM ! # of lons on this task
613 0 : HCO_Tasks(N)%JM = my_JM ! # of lats on this task
614 0 : HCO_Tasks(N)%IS = my_IS
615 0 : HCO_Tasks(N)%IE = my_IE ! start and end longitude dim'l index
616 0 : HCO_Tasks(N)%JS = my_JS
617 0 : HCO_Tasks(N)%JE = my_JE ! start and end latitude dim'l index
618 : enddo
619 :
620 : ! For root task write out a debug output to make sure
621 0 : if(masterproc) then
622 0 : write(iulog,*) ">> HEMCO: Root task committing to sub-decomp"
623 0 : write(iulog,*) ">> my_IM,JM,IS,IE,JS,JE", my_IM, my_JM, my_IS, my_IE, my_JS, my_JE
624 : endif
625 :
626 : !-----------------------------------------------------------------------
627 : ! Distribute among parallelization in MPI 3: Distribute all-to-all task info
628 : !-----------------------------------------------------------------------
629 :
630 : ! After this, exchange task information between everyone so we are all
631 : ! on the same page. This was called from edynamo_init in the ionos code.
632 : ! We adapt the whole mp_exchange_tasks code here...
633 :
634 : ! Note: 9 here is the length of the HCO_Tasks(N) component.
635 :
636 : #define HCO_TASKS_ITEM_LENGTH 9
637 0 : allocate(itasks_send(HCO_TASKS_ITEM_LENGTH, 0:nPET-1), stat=RC)
638 0 : allocate(itasks_recv(HCO_TASKS_ITEM_LENGTH, 0:nPET-1), stat=RC)
639 0 : ASSERT_(RC==0)
640 :
641 : ! Fill my send PET info array
642 0 : do N = 0, nPET-1
643 0 : itasks_send(1,N) = iam ! %ID identifier
644 0 : itasks_send(2,N) = my_ID_I ! %ID_I task coord in longitude dim'l of task table
645 0 : itasks_send(3,N) = my_ID_J ! %ID_J task coord in latitude dim'l of task table
646 0 : itasks_send(4,N) = my_IM ! %IM # of lons on this task
647 0 : itasks_send(5,N) = my_JM ! %JM # of lats on this task
648 0 : itasks_send(6,N) = my_IS ! %IS
649 0 : itasks_send(7,N) = my_IE ! %IE start and end longitude dim'l index
650 0 : itasks_send(8,N) = my_JS ! %JS
651 0 : itasks_send(9,N) = my_JE ! %JE start and end latitude dim'l index
652 : enddo
653 :
654 : ! Send MPI all-to-all
655 : call mpi_alltoall(itasks_send, HCO_TASKS_ITEM_LENGTH, MPI_INTEGER, &
656 : itasks_recv, HCO_TASKS_ITEM_LENGTH, MPI_INTEGER, &
657 0 : CAM_mpicom, RC)
658 0 : ASSERT_(RC==MPI_SUCCESS)
659 :
660 : ! Unpack receiving data back
661 0 : do N = 0, nPET-1
662 0 : HCO_Tasks(N)%ID = itasks_recv(1,N)
663 0 : HCO_Tasks(N)%ID_I = itasks_recv(2,N)
664 0 : HCO_Tasks(N)%ID_J = itasks_recv(3,N)
665 0 : HCO_Tasks(N)%IM = itasks_recv(4,N)
666 0 : HCO_Tasks(N)%JM = itasks_recv(5,N)
667 0 : HCO_Tasks(N)%IS = itasks_recv(6,N)
668 0 : HCO_Tasks(N)%IE = itasks_recv(7,N)
669 0 : HCO_Tasks(N)%JS = itasks_recv(8,N)
670 0 : HCO_Tasks(N)%JE = itasks_recv(9,N)
671 :
672 : ! Debug output for masterproc
673 : ! if(masterproc) then
674 : ! write(iulog,*) "(mp) Task ", N
675 : ! write(iulog,*) "%ID ", HCO_Tasks(N)%ID
676 : ! write(iulog,*) "%ID_I", HCO_Tasks(N)%ID_I
677 : ! write(iulog,*) "%ID_J", HCO_Tasks(N)%ID_J
678 : ! write(iulog,*) "%IM ", HCO_Tasks(N)%IM
679 : ! write(iulog,*) "%JM ", HCO_Tasks(N)%JM
680 : ! write(iulog,*) "%IS ", HCO_Tasks(N)%IS
681 : ! write(iulog,*) "%IE ", HCO_Tasks(N)%IE
682 : ! write(iulog,*) "%JS ", HCO_Tasks(N)%JS
683 : ! write(iulog,*) "%JE ", HCO_Tasks(N)%JE
684 : ! endif
685 : enddo
686 :
687 : ! Reclaim space
688 0 : deallocate(itasks_send)
689 0 : deallocate(itasks_recv)
690 :
691 0 : if(masterproc) then
692 : ! Output information on the decomposition
693 0 : write(iulog,*) ">> HEMCO DEBUG - Committed HCO_Tasks decomposition"
694 0 : write(iulog,*) "nPET, nPET_lon, nPET_lat", nPET, nPET_lon, nPET_lat
695 0 : do N = 0, nPET-1
696 0 : write(iulog,*) "PET", N, " ID", HCO_Tasks(N)%ID, " ID_I", HCO_Tasks(N)%ID_I, " ID_J", HCO_Tasks(N)%ID_J, &
697 0 : "IM (IS,IE)", HCO_Tasks(N)%IM, HCO_Tasks(N)%IS, HCO_Tasks(N)%IE, " // JM (JS, JE)", HCO_Tasks(N)%JM, HCO_Tasks(N)%JS, HCO_Tasks(N)%JE
698 : enddo
699 : endif
700 :
701 : !
702 : ! Just remember that my_IM, my_JM ... are your keys to generating
703 : ! the relevant met fields and passing to HEMCO.
704 : !
705 : ! Only HCO_ESMF_Grid should be aware of the entire grid.
706 : ! Everyone else should be just doing work on its subset indices.
707 : !
708 0 : RC = ESMF_SUCCESS
709 :
710 0 : end subroutine HCO_Grid_Init
711 : !EOC
712 : !------------------------------------------------------------------------------
713 : ! Harmonized Emissions Component (HEMCO) !
714 : !------------------------------------------------------------------------------
715 : !BOP
716 : !
717 : ! !IROUTINE: HCO_Grid_UpdateRegrid
718 : !
719 : ! !DESCRIPTION: Subroutine HCO\_Grid\_UpdateRegrid initializes or updates the
720 : ! regridding information used in the HEMCO_CESM interface.
721 : !\\
722 : !\\
723 : ! !INTERFACE:
724 : !
725 0 : subroutine HCO_Grid_UpdateRegrid( RC )
726 : !
727 : ! !USES:
728 : !
729 : ! MPI Properties from CAM
730 : ! Even though CAM's principle is that only spmd_utils uses MPI,
731 : ! ionos code uses MPI very liberally. Unfortunately we have to follow
732 : ! this example as spmd_utils does not provide many of the relevant
733 : ! functionality like the communicator split. But keep this in mind.
734 : use cam_logfile, only: iulog
735 : use spmd_utils, only: CAM_mpicom => mpicom, CAM_npes => npes
736 : use spmd_utils, only: MPI_SUCCESS
737 : use spmd_utils, only: iam
738 : use spmd_utils, only: masterproc
739 :
740 : use phys_grid, only: get_ncols_p
741 : use ppgrid, only: begchunk, endchunk
742 :
743 : use cam_instance, only: atm_id
744 :
745 : ! ESMF
746 : use ESMF, only: ESMF_FieldRegridStore
747 : use ESMF, only: ESMF_REGRIDMETHOD_BILINEAR, ESMF_REGRIDMETHOD_CONSERVE
748 : use ESMF, only: ESMF_POLEMETHOD_ALLAVG, ESMF_POLEMETHOD_NONE
749 : use ESMF, only: ESMF_EXTRAPMETHOD_NEAREST_IDAVG
750 :
751 : use ESMF, only: ESMF_FieldGet
752 : !
753 : ! !OUTPUT PARAMETERS:
754 : !
755 : integer, intent(out) :: RC
756 : !
757 : ! !REMARKS:
758 : ! This field will ONLY update if it recognizes a change in the CAM instance
759 : ! information, as determined by cam_instance::atm_id which is saved in the
760 : ! module's cam_last_atm_id field.
761 : !
762 : ! This allows the function HCO_Grid_UpdateRegrid be called both in init and run
763 : ! without performance / memory repercussions (hopefully...)
764 : !
765 : ! Maybe we can use CONSERVE_2ND order for better accuracy in the future. To be tested.
766 : !
767 : ! !REVISION HISTORY:
768 : ! 20 Feb 2020 - H.P. Lin - Initial version
769 : ! 23 Feb 2020 - H.P. Lin - Finalized regridding handles, bilinear for CAM-HCO and
770 : ! conservative (1st order) for HCO-CAM.
771 : ! 29 Oct 2020 - H.P. Lin - Changed to 1st-order conservative based on feedback from T. Fritz
772 : ! about negative numerical artifacts coming out of ESMF.
773 : !EOP
774 : !------------------------------------------------------------------------------
775 : !BOC
776 : !
777 : ! !LOCAL VARIABLES:
778 : !
779 : character(len=*), parameter :: subname = 'HCO_Grid_UpdateRegrid'
780 : integer :: smm_srctermproc, smm_pipelinedep
781 :
782 : integer :: chnk
783 :
784 : ! Debug only
785 : integer :: lbnd_hco(3), ubnd_hco(3) ! 3-d bounds of HCO field
786 : integer :: lbnd_cam(2), ubnd_cam(2) ! 3-d bounds of CAM field
787 : real(r8), pointer :: fptr(:,:,:) ! debug
788 : real(r8), pointer :: fptr_cam(:,:) ! debug
789 :
790 : ! Assume success
791 0 : RC = ESMF_SUCCESS
792 :
793 : ! Parameters for ESMF RouteHandle (taken from ionos interface)
794 0 : smm_srctermproc = 0
795 : ! smm_pipelinedep = -1 ! Accept auto-tuning of the pipeline depth
796 0 : smm_pipelinedep = 16 ! Prescribe pipeline depth for BFB consistency (hplin, 6/7/23)
797 :
798 : ! Check if we need to update
799 0 : if(cam_last_atm_id == atm_id) then
800 0 : if(masterproc) then
801 0 : write(iulog,*) ">> UpdateRegrid received ", atm_id, " already set"
802 : endif
803 : return
804 : endif
805 :
806 0 : cam_last_atm_id = atm_id
807 :
808 : ! Create CAM physics mesh...
809 0 : call HCO_Grid_ESMF_CreateCAM(RC)
810 0 : ASSERT_(RC==ESMF_SUCCESS)
811 :
812 : ! How many columns in this task? my_CE
813 : ! Note that my_CE is LOCAL column index (1:my_CE) and is equivalent to "blksize"
814 : ! in the ionos coupling code. This was 144 in my testing.
815 : ! But this is NOT ppgrid::pcols, which was 16.
816 0 : my_CE = 0
817 0 : do chnk = begchunk, endchunk
818 0 : my_CE = my_CE + get_ncols_p(chnk)
819 : enddo
820 :
821 : ! Create HEMCO grid in ESMF format
822 0 : call HCO_Grid_ESMF_CreateHCO(RC)
823 0 : ASSERT_(RC==ESMF_SUCCESS)
824 :
825 : ! Create empty fields on the HEMCO grid and CAM physics mesh
826 : ! used for later regridding
827 :
828 : ! FIXME: Destroy fields before creating to prevent memory leak? ESMF_Destroy
829 : ! requires the field to be present (so no silent destruction) (hplin, 2/21/20)
830 :
831 0 : call HCO_Grid_ESMF_CreateCAMField(CAM_2DFld, CAM_PhysMesh, 'HCO_PHYS_2DFLD', 0, RC)
832 0 : call HCO_Grid_ESMF_CreateCAMField(CAM_3DFld, CAM_PhysMesh, 'HCO_PHYS_3DFLD', LM, RC)
833 :
834 0 : call HCO_Grid_ESMF_CreateHCOField(HCO_2DFld, HCO_Grid, 'HCO_2DFLD', 0, RC)
835 0 : call HCO_Grid_ESMF_CreateHCOField(HCO_3DFld, HCO_Grid, 'HCO_3DFLD', LM, RC)
836 0 : ASSERT_(RC==ESMF_SUCCESS)
837 :
838 : ! if(masterproc) then
839 : ! write(iulog,*) ">> HEMCO: HCO_PHYS and HCO_ fields initialized successfully"
840 : ! call ESMF_FieldGet(HCO_3DFld, localDE=0, farrayPtr=fptr, &
841 : ! computationalLBound=lbnd_hco, &
842 : ! computationalUBound=ubnd_hco, rc=RC)
843 : ! write(iulog,*) ">> HEMCO: Debug HCO Field: lbnd = (", lbnd_hco, "), ", my_IS, my_IE, "ubnd = (", ubnd_hco, ")", my_JS, my_JE
844 :
845 : ! call ESMF_FieldGet(CAM_3DFld, localDE=0, farrayPtr=fptr_cam, &
846 : ! computationalLBound=lbnd_cam, &
847 : ! computationalUBound=ubnd_cam, rc=RC)
848 : ! write(iulog,*) ">> HEMCO: Debug CAM Field: lbnd = (", lbnd_cam, "), ubnd = (", ubnd_cam, ")", my_CE
849 : ! endif
850 :
851 : ! CAM -> HCO 2-D
852 : call ESMF_FieldRegridStore( &
853 : srcField=CAM_2DFld, dstField=HCO_2DFld, &
854 : regridMethod=ESMF_REGRIDMETHOD_BILINEAR, &
855 : poleMethod=ESMF_POLEMETHOD_ALLAVG, &
856 : extrapMethod=ESMF_EXTRAPMETHOD_NEAREST_IDAVG, &
857 : routeHandle=CAM2HCO_RouteHandle_2D, &
858 : !srcTermProcessing=smm_srctermproc, &
859 0 : pipelineDepth=smm_pipelinedep, rc=RC)
860 0 : ASSERT_(RC==ESMF_SUCCESS)
861 :
862 0 : if(masterproc) then
863 0 : write(iulog,*) ">> After FieldRegridStore, CAM2D->HCO2D, pipeline", smm_pipelinedep
864 : endif
865 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
866 :
867 : ! Create and store ESMF route handles for regridding CAM <-> HCO
868 : ! CAM -> HCO 3-D
869 : call ESMF_FieldRegridStore( &
870 : srcField=CAM_3DFld, dstField=HCO_3DFld, &
871 : regridMethod=ESMF_REGRIDMETHOD_BILINEAR, &
872 : poleMethod=ESMF_POLEMETHOD_ALLAVG, &
873 : extrapMethod=ESMF_EXTRAPMETHOD_NEAREST_IDAVG, &
874 : routeHandle=CAM2HCO_RouteHandle_3D, &
875 : !srcTermProcessing=smm_srctermproc, &
876 0 : pipelineDepth=smm_pipelinedep, rc=RC)
877 0 : ASSERT_(RC==ESMF_SUCCESS)
878 :
879 0 : if(masterproc) then
880 0 : write(iulog,*) ">> After FieldRegridStore, CAM3D->HCO3D, pipeline", smm_pipelinedep
881 : endif
882 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
883 :
884 : ! HCO -> CAM 2-D
885 : call ESMF_FieldRegridStore( &
886 : srcField=HCO_2DFld, dstField=CAM_2DFld, &
887 : regridMethod=ESMF_REGRIDMETHOD_CONSERVE , &
888 : poleMethod=ESMF_POLEMETHOD_NONE, &
889 : routeHandle=HCO2CAM_RouteHandle_2D, &
890 : !srcTermProcessing=smm_srctermproc, &
891 0 : pipelineDepth=smm_pipelinedep, rc=RC)
892 0 : ASSERT_(RC==ESMF_SUCCESS)
893 :
894 0 : if(masterproc) then
895 0 : write(iulog,*) ">> After FieldRegridStore, HCO2D->CAM2D, pipeline", smm_pipelinedep
896 : endif
897 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
898 :
899 : ! HCO -> CAM 3-D
900 : ! 3-D conserv. regridding cannot be done on a stagger other than center
901 : ! (as of ESMF 8.0.0 in ESMF_FieldRegrid.F90::993)
902 : call ESMF_FieldRegridStore( &
903 : srcField=HCO_3DFld, dstField=CAM_3DFld, &
904 : regridMethod=ESMF_REGRIDMETHOD_CONSERVE , &
905 : poleMethod=ESMF_POLEMETHOD_NONE, &
906 : routeHandle=HCO2CAM_RouteHandle_3D, &
907 : !srcTermProcessing=smm_srctermproc, &
908 0 : pipelineDepth=smm_pipelinedep, rc=RC)
909 0 : ASSERT_(RC==ESMF_SUCCESS)
910 :
911 0 : if(masterproc) then
912 0 : write(iulog,*) ">> After FieldRegridStore, HCO3D->CAM3D, pipeline", smm_pipelinedep
913 : endif
914 :
915 0 : if(masterproc) then
916 0 : write(iulog,*) ">> HEMCO: FieldRegridStore four-way ready"
917 : endif
918 :
919 : ! Finished updating regrid routines
920 0 : end subroutine HCO_Grid_UpdateRegrid
921 : !EOC
922 : !------------------------------------------------------------------------------
923 : ! Harmonized Emissions Component (HEMCO) !
924 : !------------------------------------------------------------------------------
925 : !BOP
926 : !
927 : ! !IROUTINE: HCO_Grid_ESMF_CreateCAM
928 : !
929 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateCAM creates the physics mesh
930 : ! from CAM and stores in the ESMF state for regridding.
931 : !\\
932 : !\\
933 : ! !INTERFACE:
934 : !
935 0 : subroutine HCO_Grid_ESMF_CreateCAM( RC )
936 : !
937 : ! !USES:
938 : !
939 : ! MPI Properties from CAM
940 0 : use cam_logfile, only: iulog
941 : use spmd_utils, only: masterproc
942 :
943 : ! Phys constants
944 : use shr_const_mod, only: pi => shr_const_pi
945 :
946 : ! Grid properties in CAM
947 : use cam_instance, only: inst_name
948 : use phys_control, only: phys_getopts
949 : use phys_grid, only: get_ncols_p, get_gcol_p, get_rlon_all_p, get_rlat_all_p
950 : use ppgrid, only: begchunk, endchunk
951 : use ppgrid, only: pcols ! # of col max in chunks
952 :
953 : ! ESMF
954 : use ESMF, only: ESMF_DistGridCreate, ESMF_MeshCreate
955 : use ESMF, only: ESMF_MeshGet
956 : use ESMF, only: ESMF_FILEFORMAT_ESMFMESH
957 : use ESMF, only: ESMF_MeshIsCreated, ESMF_MeshDestroy
958 : !
959 : ! !OUTPUT PARAMETERS:
960 : !
961 : integer, intent(out) :: RC
962 : !
963 : ! !REMARKS:
964 : !
965 : ! !REVISION HISTORY:
966 : ! 11 Feb 2020 - H.P. Lin - Initial version
967 : !EOP
968 : !------------------------------------------------------------------------------
969 : !BOC
970 : !
971 : ! !LOCAL VARIABLES:
972 : !
973 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateCAM'
974 :
975 : ! For allocation of the distGrid and Mesh
976 : integer :: ncols
977 : integer :: chnk, col, dindex
978 0 : integer, allocatable :: decomp(:)
979 : integer :: col_total
980 : character(len=256) :: grid_file
981 :
982 : ! For verification of the mesh
983 : integer :: spatialDim
984 : integer :: numOwnedElements
985 0 : real(r8), pointer :: ownedElemCoords(:)
986 0 : real(r8), pointer :: latCAM(:), latMesh(:)
987 0 : real(r8), pointer :: lonCAM(:), lonMesh(:)
988 : real(r8) :: latCAM_R(pcols) ! array of chunk lat
989 : real(r8) :: lonCAM_R(pcols) ! array of chunk long
990 :
991 : integer :: i, c, n
992 : real(r8), parameter :: radtodeg = 180.0_r8/pi
993 :
994 : ! Assume success
995 0 : RC = ESMF_SUCCESS
996 :
997 : !-----------------------------------------------------------------------
998 : ! Compute the CAM_DistGrid and CAM_PhysMesh
999 : !-----------------------------------------------------------------------
1000 : ! Get the physics grid information
1001 0 : call phys_getopts(physics_grid_out=grid_file)
1002 :
1003 0 : if(masterproc) then
1004 0 : write(iulog,*) "physics_grid_out=", grid_file
1005 : endif
1006 :
1007 : ! Compute local decomposition (global variable in-module)
1008 0 : col_total = 0 ! Sum of columns on this PET in all chunks
1009 0 : do chnk = begchunk, endchunk
1010 0 : col_total = col_total + get_ncols_p(chnk)
1011 : enddo
1012 0 : allocate(decomp(col_total))
1013 0 : allocate(lonCAM(col_total))
1014 0 : allocate(latCAM(col_total)) ! the _R variants are already allocated to pcols
1015 :
1016 : ! ...and also attach to the loop computing lat and lons for CAM physics mesh
1017 : ! to be used later.
1018 0 : dindex = 0
1019 0 : do chnk = begchunk, endchunk
1020 0 : ncols = get_ncols_p(chnk)
1021 :
1022 : ! Get [rad] lat and lons
1023 0 : call get_rlon_all_p(chnk, ncols, lonCAM_R)
1024 0 : call get_rlat_all_p(chnk, ncols, latCAM_R)
1025 :
1026 0 : do col = 1, ncols
1027 0 : dindex = dindex + 1
1028 0 : decomp(dindex) = get_gcol_p(chnk, col)
1029 :
1030 0 : lonCAM(dindex) = lonCAM_R(col) * radtodeg
1031 0 : latCAM(dindex) = latCAM_R(col) * radtodeg
1032 : enddo
1033 : enddo
1034 :
1035 : ! Build the 2D field CAM DistGrid based on the physics decomposition
1036 0 : CAM_DistGrid = ESMF_DistGridCreate(arbSeqIndexList=decomp, rc=RC)
1037 0 : ASSERT_(RC==ESMF_SUCCESS)
1038 :
1039 : ! Release memory if any is being taken, to avoid leakage
1040 0 : if(ESMF_MeshIsCreated(CAM_PhysMesh)) then
1041 0 : call ESMF_MeshDestroy(CAM_PhysMesh)
1042 : endif
1043 :
1044 : ! Create the physics decomposition ESMF mesh
1045 : CAM_PhysMesh = ESMF_MeshCreate(trim(grid_file), ESMF_FILEFORMAT_ESMFMESH, &
1046 0 : elementDistGrid=CAM_DistGrid, rc=RC)
1047 0 : ASSERT_(RC==ESMF_SUCCESS)
1048 :
1049 : !-----------------------------------------------------------------------
1050 : ! Validate mesh coordinates against model physics column coords.
1051 : !-----------------------------------------------------------------------
1052 : ! (From edyn_esmf::edyn_create_physmesh)
1053 : call ESMF_MeshGet(CAM_PhysMesh, spatialDim=spatialDim, &
1054 0 : numOwnedElements=numOwnedElements, rc=RC)
1055 0 : ASSERT_(RC==ESMF_SUCCESS)
1056 :
1057 0 : if(numOwnedElements /= col_total) then
1058 0 : write(iulog,*) "HEMCO: ESMF_MeshGet numOwnedElements =", numOwnedElements, &
1059 0 : "col_total =", col_total, " MISMATCH! Aborting"
1060 0 : ASSERT_(.false.)
1061 : endif
1062 :
1063 : ! Coords for the CAM_PhysMesh
1064 0 : allocate(ownedElemCoords(spatialDim * numOwnedElements))
1065 0 : allocate(lonMesh(col_total), latMesh(col_total))
1066 :
1067 0 : call ESMF_MeshGet(CAM_PhysMesh, ownedElemCoords=ownedElemCoords, rc=RC)
1068 :
1069 0 : do n = 1, col_total
1070 0 : lonMesh(n) = ownedElemCoords(2*n-1)
1071 0 : latMesh(n) = ownedElemCoords(2*n)
1072 : enddo
1073 :
1074 : ! Error check coordinates
1075 0 : do n = 1, col_total
1076 0 : if(abs(lonMesh(n) - lonCAM(n)) > 0.000001_r8) then
1077 0 : if((abs(lonMesh(n) - lonCAM(n)) > 360.000001_r8) .or. &
1078 : (abs(lonMesh(n) - lonCAM(n)) < 359.99999_r8)) then
1079 0 : write(6,*) "HEMCO: ESMF_MeshGet VERIFY fail! n, lonMesh, lonCAM, delta"
1080 0 : write(6,*) n, lonMesh(n), lonCAM(n), abs(lonMesh(n)-lonCAM(n))
1081 0 : ASSERT_(.false.)
1082 : endif
1083 : endif
1084 :
1085 0 : if(abs(latMesh(n) - latCAM(n)) > 0.000001_r8) then
1086 0 : if(.not. ((abs(latCAM(n)) > 88.0_r8) .and. (abs(latMesh(n)) > 88.0_r8))) then
1087 0 : write(6,*) "HEMCO: ESMF_MeshGet VERIFY fail! n, latmesh, latCAM, delta"
1088 0 : write(6,*) n, latMesh(n), latCAM(n), abs(latMesh(n)-latCAM(n))
1089 0 : ASSERT_(.false.)
1090 : endif
1091 : endif
1092 : enddo
1093 :
1094 : ! Ready to go
1095 0 : if(masterproc) then
1096 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateCAM ok, dim'l = ", col_total
1097 : endif
1098 :
1099 : ! Free memory
1100 0 : deallocate(ownedElemCoords)
1101 0 : deallocate(lonCAM, lonMesh)
1102 0 : deallocate(latCAM, latMesh)
1103 0 : deallocate(decomp)
1104 :
1105 0 : end subroutine HCO_Grid_ESMF_CreateCAM
1106 : !EOC
1107 : !------------------------------------------------------------------------------
1108 : ! Harmonized Emissions Component (HEMCO) !
1109 : !------------------------------------------------------------------------------
1110 : !BOP
1111 : !
1112 : ! !IROUTINE: HCO_Grid_ESMF_CreateHCO
1113 : !
1114 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateHCO creates the HEMCO grid
1115 : ! in center and corner staggering modes in ESMF_Grid format,
1116 : ! and stores in the ESMF state for regridding.
1117 : !\\
1118 : !\\
1119 : ! !INTERFACE:
1120 : !
1121 0 : subroutine HCO_Grid_ESMF_CreateHCO( RC )
1122 : !
1123 : ! !USES:
1124 : !
1125 : ! MPI Properties from CAM
1126 0 : use cam_logfile, only: iulog
1127 : use spmd_utils, only: masterproc
1128 :
1129 : ! ESMF
1130 : use ESMF, only: ESMF_GridCreate1PeriDim, ESMF_INDEX_GLOBAL
1131 : use ESMF, only: ESMF_STAGGERLOC_CENTER, ESMF_STAGGERLOC_CORNER
1132 : use ESMF, only: ESMF_GridAddCoord, ESMF_GridGetCoord
1133 : !
1134 : ! !OUTPUT PARAMETERS:
1135 : !
1136 : integer, intent(out) :: RC
1137 : !
1138 : ! !REMARKS:
1139 : ! We initialize TWO coordinates here for the HEMCO grid. Both the center and corner
1140 : ! staggering information needs to be added to the grid, or ESMF_FieldRegridStore will
1141 : ! fail. It is a rather weird implementation, as ESMF_FieldRegridStore only accepts
1142 : ! the coordinate in CENTER staggering, but the conversion to mesh is in CORNER.
1143 : !
1144 : ! So both coordinates need to be specified even though (presumably?) the center field
1145 : ! is regridded once the route handle is generated. This remains to be seen but at this
1146 : ! point only the following "duplicated" code is the correct implementation for
1147 : ! conservative regridding handles to be generated properly.
1148 : !
1149 : ! Roughly 9 hours of debugging and reading the ESMF documentation and code
1150 : ! were wasted here. (hplin, 2/23/20)
1151 : !
1152 : ! !REVISION HISTORY:
1153 : ! 20 Feb 2020 - H.P. Lin - Initial version
1154 : ! 22 Feb 2020 - H.P. Lin - Support curvilinear in ESMF format (note Map_A2A)
1155 : !EOP
1156 : !------------------------------------------------------------------------------
1157 : !BOC
1158 : !
1159 : ! !LOCAL VARIABLES:
1160 : !
1161 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateHCO'
1162 :
1163 : integer :: i, j, n
1164 : integer :: lbnd(2), ubnd(2)
1165 0 : integer :: nlons_task(nPET_lon) ! # lons per task
1166 0 : integer :: nlats_task(nPET_lat) ! # lats per task
1167 0 : real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:)
1168 0 : real(ESMF_KIND_R8), pointer :: coordX_E(:,:), coordY_E(:,:)
1169 :
1170 : !-----------------------------------------------------------------------
1171 : ! Task distribution for ESMF grid
1172 : !-----------------------------------------------------------------------
1173 0 : do i = 1, nPET_lon
1174 0 : loop: do n = 0, nPET-1
1175 0 : if (HCO_Tasks(n)%ID_I == i-1) then
1176 0 : nlons_task(i) = HCO_Tasks(n)%IM
1177 0 : exit loop
1178 : endif
1179 : enddo loop
1180 : enddo
1181 :
1182 : ! Exclude periodic points for source grids
1183 : ! do n = 0, nPET-1
1184 : ! if (HCO_Tasks(n)%ID_I == nPET_lon-1) then ! Eastern edge
1185 : ! nlons_task(HCO_Tasks(n)%ID_I + 1) = HCO_Tasks(n)%IM - 1
1186 : ! ! overwrites %IM above...
1187 : ! endif
1188 : ! enddo
1189 :
1190 0 : do j = 1, nPET_lat
1191 0 : loop1: do n = 0, nPET-1
1192 0 : if (HCO_Tasks(n)%ID_J == j-1) then
1193 0 : nlats_task(j) = HCO_Tasks(n)%JM
1194 0 : exit loop1
1195 : endif
1196 : enddo loop1
1197 : enddo
1198 :
1199 : !-----------------------------------------------------------------------
1200 : ! Create source grids and allocate coordinates.
1201 : !-----------------------------------------------------------------------
1202 : ! Create pole-based 2D geographic source grid
1203 : HCO_Grid = ESMF_GridCreate1PeriDim( &
1204 : countsPerDEDim1=nlons_task, coordDep1=(/1,2/), &
1205 : countsPerDEDim2=nlats_task, coordDep2=(/1,2/), &
1206 : indexflag=ESMF_INDEX_GLOBAL, &
1207 : petmap=HCO_petMap, &
1208 : minIndex=(/1,1/), &
1209 0 : rc=RC)
1210 0 : ASSERT_(RC==ESMF_SUCCESS)
1211 :
1212 0 : call ESMF_GridAddCoord(HCO_Grid, staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1213 0 : ASSERT_(RC==ESMF_SUCCESS)
1214 :
1215 0 : call ESMF_GridAddCoord(HCO_Grid, staggerloc=ESMF_STAGGERLOC_CORNER, rc=RC)
1216 0 : ASSERT_(RC==ESMF_SUCCESS) ! why two grids? see remark above (hplin, 2/23/20)
1217 :
1218 : !-----------------------------------------------------------------------
1219 : ! Get pointer and set coordinates - CENTER HEMCO GRID
1220 : !-----------------------------------------------------------------------
1221 : call ESMF_GridGetCoord(HCO_Grid, coordDim=1, localDE=0, &
1222 : computationalLBound=lbnd, &
1223 : computationalUBound=ubnd, &
1224 : farrayPtr=coordX, &
1225 : staggerloc=ESMF_STAGGERLOC_CENTER, &
1226 0 : rc=RC)
1227 : ! Longitude range -180.0, 180.0 is XMid for center staggering
1228 :
1229 : ! Note: Compute bounds are not starting from 1 almost surely
1230 : ! and they should be on the same decomp as the global elems.
1231 : ! So they should be read through the global indices
1232 : ! and not offset ones (a la WRF) (hplin, 2/21/20)
1233 0 : do j = lbnd(2), ubnd(2)
1234 0 : do i = lbnd(1), ubnd(1)
1235 0 : coordX(i, j) = XMid(i, j)
1236 : enddo
1237 : enddo
1238 0 : ASSERT_(RC==ESMF_SUCCESS)
1239 :
1240 : ! Sanity checks here for coordinate staggering
1241 : ! write(iulog,*) ">> [X] IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE, &
1242 : ! "IM, JM, sizeX1,2 = ", &
1243 : ! my_IM, my_JM, size(coordX, 1), size(coordX, 2)
1244 : ! Off by one in periodic dim'n, don't assert I dim here
1245 : ! ASSERT_((ubnd(1) - lbnd(1))==(my_IE-my_IS))
1246 : ! ASSERT_((ubnd(2) - lbnd(2))==(my_JE-my_JS))
1247 :
1248 : call ESMF_GridGetCoord(HCO_Grid, coordDim=2, localDE=0, &
1249 : computationalLBound=lbnd, &
1250 : computationalUBound=ubnd, &
1251 : farrayPtr=coordY, &
1252 : staggerloc=ESMF_STAGGERLOC_CENTER, &
1253 0 : rc=RC)
1254 :
1255 : ! Sanity checks here for coordinate staggering
1256 : ! write(iulog,*) ">> [Y] IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE, &
1257 : ! "IM, JM, sizeY1,2 = ", &
1258 : ! my_IM, my_JM, size(coordY, 1), size(coordY, 2)
1259 : ! ASSERT_((ubnd(1) - lbnd(1))==(my_IE-my_IS))
1260 : ! ASSERT_((ubnd(2) - lbnd(2))==(my_JE-my_JS))
1261 :
1262 0 : do j = lbnd(2), ubnd(2)
1263 0 : do i = lbnd(1), ubnd(1)
1264 0 : coordY(i, j) = YMid(i, j)
1265 : enddo
1266 : enddo
1267 0 : ASSERT_(RC==ESMF_SUCCESS)
1268 :
1269 0 : if(masterproc) then
1270 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateHCO [Ctr] ok"
1271 0 : write(iulog,*) ">> lbnd,ubnd_lon = ", lbnd(1), ubnd(1)
1272 0 : write(iulog,*) ">> lbnd,ubnd_lat = ", lbnd(2), ubnd(2)
1273 0 : write(iulog,*) ">> IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE
1274 0 : write(iulog,*) ">> IM, JM, sizeX1,2 sizeY1,2 = ", &
1275 0 : my_IM, my_JM, size(coordX, 1), size(coordX, 2), size(coordY, 1), size(coordY, 2)
1276 : endif
1277 :
1278 : !-----------------------------------------------------------------------
1279 : ! Get pointer and set coordinates - CORNER HEMCO GRID
1280 : !-----------------------------------------------------------------------
1281 : call ESMF_GridGetCoord(HCO_Grid, coordDim=1, localDE=0, &
1282 : computationalLBound=lbnd, &
1283 : computationalUBound=ubnd, &
1284 : farrayPtr=coordX_E, &
1285 : staggerloc=ESMF_STAGGERLOC_CORNER, &
1286 0 : rc=RC)
1287 :
1288 : ! lbnd(1), ubnd(1) -> 1, 180; 181, 360 ... same as IS, IE ...
1289 :
1290 :
1291 : ! Note: Compute bounds are not starting from 1 almost surely
1292 : ! and they should be on the same decomp as the global elems.
1293 : ! So they should be read through the global indices
1294 : ! and not offset ones (a la WRF) (hplin, 2/21/20)
1295 :
1296 : ! It is a rectilinear grid - so it will not matter what j you pick
1297 : ! for the x-dimension edges, and vice-versa. This is a crude assumption
1298 : ! that is correct for rectilinear but should be revisited. The code
1299 : ! itself is capable of much more. (hplin, 5/29/20)
1300 0 : do j = lbnd(2), ubnd(2)
1301 0 : do i = lbnd(1), ubnd(1)
1302 0 : coordX_E(i, j) = XEdge(i, min(JM, j))
1303 : enddo
1304 : enddo
1305 0 : ASSERT_(RC==ESMF_SUCCESS)
1306 :
1307 : call ESMF_GridGetCoord(HCO_Grid, coordDim=2, localDE=0, &
1308 : computationalLBound=lbnd, &
1309 : computationalUBound=ubnd, &
1310 : farrayPtr=coordY_E, &
1311 : staggerloc=ESMF_STAGGERLOC_CORNER, &
1312 0 : rc=RC)
1313 :
1314 0 : do j = lbnd(2), ubnd(2)
1315 0 : do i = lbnd(1), ubnd(1)
1316 0 : coordY_E(i, j) = YEdge(min(IM, i), j)
1317 : enddo
1318 : enddo
1319 0 : ASSERT_(RC==ESMF_SUCCESS)
1320 :
1321 0 : if(masterproc) then
1322 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateHCO [Cnr] ok"
1323 0 : write(iulog,*) ">> lbnd,ubnd_lon = ", lbnd(1), ubnd(1)
1324 0 : write(iulog,*) ">> lbnd,ubnd_lat = ", lbnd(2), ubnd(2)
1325 0 : write(iulog,*) ">> IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE
1326 0 : write(iulog,*) ">> IM, JM, sizeX1,2 sizeY1,2 = ", &
1327 0 : my_IM, my_JM, size(coordX, 1), size(coordX, 2), size(coordY, 1), size(coordY, 2)
1328 : endif
1329 :
1330 0 : end subroutine HCO_Grid_ESMF_CreateHCO
1331 : !EOC
1332 : !------------------------------------------------------------------------------
1333 : ! Harmonized Emissions Component (HEMCO) !
1334 : !------------------------------------------------------------------------------
1335 : !BOP
1336 : !
1337 : ! !IROUTINE: HCO_Grid_ESMF_CreateCAMField
1338 : !
1339 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateCAMField creates an ESMF
1340 : ! 2D or 3D field on the mesh representation of CAM physics decomp.
1341 : !\\
1342 : !\\
1343 : ! !INTERFACE:
1344 : !
1345 0 : subroutine HCO_Grid_ESMF_CreateCAMField( field, mesh, name, nlev, RC )
1346 : !
1347 : ! !USES:
1348 : !
1349 : ! MPI Properties from CAM
1350 0 : use cam_logfile, only: iulog
1351 : use spmd_utils, only: masterproc
1352 :
1353 : use ESMF, only: ESMF_TYPEKIND_R8
1354 : use ESMF, only: ESMF_MESHLOC_ELEMENT
1355 : use ESMF, only: ESMF_ArraySpec, ESMF_ArraySpecSet
1356 : use ESMF, only: ESMF_FieldCreate
1357 : !
1358 : ! !INPUT PARAMETERS:
1359 : !
1360 : type(ESMF_Mesh), intent(in) :: mesh
1361 : character(len=*), intent(in) :: name
1362 : integer, intent(in) :: nlev ! nlev==0?2d:3d
1363 : !
1364 : ! !OUTPUT PARAMETERS:
1365 : !
1366 : type(ESMF_Field), intent(out) :: field
1367 : integer, intent(out) :: RC
1368 : !
1369 : ! !REMARKS:
1370 : ! If nlev == 0, field is 2D (i, j), otherwise 3D. 3rd dimension is ungridded.
1371 : !
1372 : ! !REVISION HISTORY:
1373 : ! 21 Feb 2020 - H.P. Lin - Initial version
1374 : !EOP
1375 : !------------------------------------------------------------------------------
1376 : !BOC
1377 : !
1378 : ! !LOCAL VARIABLES:
1379 : !
1380 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateCAMField'
1381 : type(ESMF_ArraySpec) :: arrayspec
1382 :
1383 0 : if(nlev > 0) then
1384 : ! 3D field (i,j,k) with nondistributed vertical
1385 0 : call ESMF_ArraySpecSet(arrayspec, 2, ESMF_TYPEKIND_R8, rc=RC)
1386 0 : ASSERT_(RC==ESMF_SUCCESS)
1387 :
1388 : field = ESMF_FieldCreate(mesh, arrayspec, &
1389 : name=name, &
1390 : ungriddedLBound=(/1/), &
1391 : ungriddedUBound=(/nlev/), &
1392 : gridToFieldMap=(/2/), & ! mapping between grid/field dims ??
1393 0 : meshloc=ESMF_MESHLOC_ELEMENT, rc=RC)
1394 0 : ASSERT_(RC==ESMF_SUCCESS)
1395 : else
1396 : ! 2D field (i,j)
1397 0 : call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=RC)
1398 0 : ASSERT_(RC==ESMF_SUCCESS)
1399 :
1400 : field = ESMF_FieldCreate(mesh, arrayspec, &
1401 : name=name, &
1402 0 : meshloc=ESMF_MESHLOC_ELEMENT, rc=RC)
1403 0 : ASSERT_(RC==ESMF_SUCCESS)
1404 : endif
1405 0 : end subroutine HCO_Grid_ESMF_CreateCAMField
1406 : !EOC
1407 : !------------------------------------------------------------------------------
1408 : ! Harmonized Emissions Component (HEMCO) !
1409 : !------------------------------------------------------------------------------
1410 : !BOP
1411 : !
1412 : ! !IROUTINE: HCO_Grid_ESMF_CreateHCOField
1413 : !
1414 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateHCO field creates an ESMF
1415 : ! 2D or 3D field on the HEMCO grid, excluding periodic points.
1416 : !\\
1417 : !\\
1418 : ! !INTERFACE:
1419 : !
1420 0 : subroutine HCO_Grid_ESMF_CreateHCOField( field, grid, name, nlev, RC )
1421 : !
1422 : ! !USES:
1423 : !
1424 : ! MPI Properties from CAM
1425 0 : use cam_logfile, only: iulog
1426 : use spmd_utils, only: masterproc
1427 :
1428 : use ESMF, only: ESMF_TYPEKIND_R8
1429 : use ESMF, only: ESMF_STAGGERLOC_CENTER
1430 : use ESMF, only: ESMF_ArraySpec, ESMF_ArraySpecSet
1431 : use ESMF, only: ESMF_FieldCreate
1432 :
1433 : use ESMF, only: ESMF_GridGet, ESMF_Array, ESMF_ArrayCreate, ESMF_INDEX_GLOBAL
1434 : !
1435 : ! !INPUT PARAMETERS:
1436 : !
1437 : type(ESMF_Grid), intent(in) :: grid
1438 : character(len=*), intent(in) :: name
1439 : integer, intent(in) :: nlev ! nlev==0?2d:3d
1440 : !
1441 : ! !OUTPUT PARAMETERS:
1442 : !
1443 : type(ESMF_Field), intent(out) :: field
1444 : integer, intent(out) :: RC
1445 : !
1446 : ! !REMARKS:
1447 : ! If nlev == 0, field is 2D (i, j), otherwise 3D. 3rd dimension is ungridded.
1448 : ! The grid input parameter accepts both HCO_Grid and HCO2CAM_Grid.
1449 : !
1450 : ! !REVISION HISTORY:
1451 : ! 21 Feb 2020 - H.P. Lin - Initial version
1452 : !EOP
1453 : !------------------------------------------------------------------------------
1454 : !BOC
1455 : !
1456 : ! !LOCAL VARIABLES:
1457 : !
1458 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateHCOField'
1459 : type(ESMF_ArraySpec) :: arrayspec
1460 :
1461 : type(ESMF_Array) :: array3D, array2D
1462 : type(ESMF_DistGrid) :: distgrid
1463 :
1464 : ! Get grid information from the HEMCO grid:
1465 : call ESMF_GridGet(grid, staggerloc=ESMF_STAGGERLOC_CENTER, &
1466 0 : distgrid=distgrid, rc=RC)
1467 0 : ASSERT_(RC==ESMF_SUCCESS)
1468 :
1469 0 : call ESMF_LogWrite("After ESMF_GridGet", ESMF_LOGMSG_INFO, rc=RC)
1470 :
1471 0 : if(nlev > 0) then
1472 : ! 3D field (i,j,k) with nondistributed vertical
1473 0 : call ESMF_ArraySpecSet(arrayspec, 3, ESMF_TYPEKIND_R8, rc=RC)
1474 0 : ASSERT_(RC==ESMF_SUCCESS)
1475 :
1476 0 : call ESMF_LogWrite("After array3D-ESMF_ArraySpecSet", ESMF_LOGMSG_INFO, rc=RC)
1477 :
1478 : array3D = ESMF_ArrayCreate(arrayspec=arrayspec, &
1479 : distgrid=distgrid, &
1480 : computationalEdgeUWidth=(/1,0/), &
1481 : undistLBound=(/1/), undistUBound=(/nlev/), &
1482 0 : indexflag=ESMF_INDEX_GLOBAL, rc=RC)
1483 0 : ASSERT_(RC==ESMF_SUCCESS)
1484 :
1485 0 : call ESMF_LogWrite("After array3D-ESMF_ArrayCreate", ESMF_LOGMSG_INFO, rc=RC)
1486 :
1487 : field = ESMF_FieldCreate(grid, array3D, & ! grid, arrayspec
1488 : name=name, &
1489 : ungriddedLBound=(/1/), &
1490 : ungriddedUBound=(/nlev/), &
1491 0 : staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1492 0 : ASSERT_(RC==ESMF_SUCCESS)
1493 :
1494 0 : call ESMF_LogWrite("After array3D-ESMF_FieldCreate", ESMF_LOGMSG_INFO, rc=RC)
1495 : else
1496 : ! 2D field (i,j)
1497 0 : call ESMF_ArraySpecSet(arrayspec, 2, ESMF_TYPEKIND_R8, rc=RC)
1498 0 : ASSERT_(RC==ESMF_SUCCESS)
1499 :
1500 : field = ESMF_FieldCreate(grid, arrayspec, &
1501 : name=name, &
1502 0 : staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1503 0 : ASSERT_(RC==ESMF_SUCCESS)
1504 : endif
1505 0 : end subroutine HCO_Grid_ESMF_CreateHCOField
1506 : !EOC
1507 : !------------------------------------------------------------------------------
1508 : ! Harmonized Emissions Component (HEMCO) !
1509 : !------------------------------------------------------------------------------
1510 : !BOP
1511 : !
1512 : ! !IROUTINE: HCO_Grid_HCO2CAM_3D
1513 : !
1514 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_3D regrids a HEMCO lat-lon grid
1515 : ! field (i,j,l) to CAM physics mesh (k,i).
1516 : !\\
1517 : !\\
1518 : ! !INTERFACE:
1519 : !
1520 0 : subroutine HCO_Grid_HCO2CAM_3D(hcoArray, camArray)
1521 : !
1522 : ! !USES:
1523 : !
1524 0 : use ESMF, only: ESMF_FieldRegrid
1525 : use ESMF, only: ESMF_TERMORDER_SRCSEQ, ESMF_TERMORDER_FREE
1526 :
1527 : ! MPI Properties from CAM
1528 : use cam_logfile, only: iulog
1529 : use spmd_utils, only: masterproc
1530 : !
1531 : ! !INPUT PARAMETERS:
1532 : !
1533 : real(r8), intent(in) :: hcoArray(my_IS:my_IE,my_JS:my_JE,1:LM)
1534 : !
1535 : ! !OUTPUT PARAMETERS:
1536 : !
1537 : real(r8), intent(inout) :: camArray(1:LM,1:my_CE) ! Col and lev start on 1
1538 : !
1539 : ! !REMARKS:
1540 : ! (1) There is no vertical regridding. Also, chunk and lev indices are assumed to
1541 : ! start at 1 always.
1542 : ! (2) The subroutine's inner workings abstract ESMF from the user, but this routine
1543 : ! needs to be called from within the gridded component (as suggested by Steve)
1544 : ! (3) Note that this automatically flips the vertical through HCO_ESMF_Get3DField.
1545 : ! In CAM, layer 1 is TOA. Layer 1 is ground in HEMCO and thus we flip the vertical
1546 : ! in regrid routines when retrieving the data.
1547 : !
1548 : ! !REVISION HISTORY:
1549 : ! 24 Feb 2020 - H.P. Lin - Initial version
1550 : ! 30 May 2020 - H.P. Lin - Remember to flip the vertical!
1551 : ! 06 Dec 2020 - H.P. Lin - Add a kludge to skip regridding (removed 12 Jan 2023)
1552 : !EOP
1553 : !------------------------------------------------------------------------------
1554 : !BOC
1555 : !
1556 : ! !LOCAL VARIABLES:
1557 : !
1558 : character(len=*), parameter :: subname = 'HCO_Grid_HCO2CAM_3D'
1559 : integer :: RC
1560 :
1561 0 : call HCO_ESMF_Set3DHCO(HCO_3DFld, hcoArray, my_IS, my_IE, my_JS, my_JE, 1, LM)
1562 :
1563 :
1564 : call ESMF_FieldRegrid(HCO_3DFld, CAM_3DFld, HCO2CAM_RouteHandle_3D, &
1565 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1566 0 : checkflag=.true., rc=RC)
1567 0 : ASSERT_(RC==ESMF_SUCCESS)
1568 :
1569 : ! (field_in, data_out, IS, IE, JS, JE)
1570 : ! For chunks, "I" is lev, "J" is chunk index, confusing, you are warned
1571 : ! (Physics "2D" fields on mesh are actually "3D" data)
1572 0 : call HCO_ESMF_Get2DField(CAM_3DFld, camArray, 1, LM, 1, my_CE, flip=.true.)
1573 :
1574 0 : end subroutine HCO_Grid_HCO2CAM_3D
1575 : !EOC
1576 : !------------------------------------------------------------------------------
1577 : ! Harmonized Emissions Component (HEMCO) !
1578 : !------------------------------------------------------------------------------
1579 : !BOP
1580 : !
1581 : ! !IROUTINE: HCO_Grid_CAM2HCO_3D
1582 : !
1583 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_3D regrids a HEMCO lat-lon grid
1584 : ! field (i,j,l) to CAM physics mesh (k,i).
1585 : !\\
1586 : !\\
1587 : ! !INTERFACE:
1588 : !
1589 0 : subroutine HCO_Grid_CAM2HCO_3D(camArray, hcoArray)
1590 : !
1591 : ! !USES:
1592 : !
1593 0 : use ESMF, only: ESMF_FieldRegrid, ESMF_TERMORDER_SRCSEQ
1594 :
1595 : ! MPI Properties from CAM
1596 : use cam_logfile, only: iulog
1597 : use spmd_utils, only: masterproc
1598 : !
1599 : ! !INPUT PARAMETERS:
1600 : !
1601 : real(r8), intent(in) :: camArray(1:LM,1:my_CE) ! Col and lev start on 1
1602 : !
1603 : ! !OUTPUT PARAMETERS:
1604 : !
1605 : real(r8), intent(inout) :: hcoArray(my_IS:my_IE,my_JS:my_JE,1:LM)
1606 : !
1607 : ! !REMARKS:
1608 : ! (1) There is no vertical regridding. Also, chunk and lev indices are assumed to
1609 : ! start at 1 always.
1610 : ! (2) The subroutine's inner workings abstract ESMF from the user, but this routine
1611 : ! needs to be called from within the gridded component (as suggested by Steve)
1612 : ! (3) Note that this automatically flips the vertical through HCO_ESMF_Get3DField.
1613 : ! In CAM, layer 1 is TOA. Layer 1 is ground in HEMCO and thus we flip the vertical
1614 : ! in regrid routines when retrieving the data.
1615 : !
1616 : ! !REVISION HISTORY:
1617 : ! 24 Feb 2020 - H.P. Lin - Initial version
1618 : ! 30 May 2020 - H.P. Lin - Remember to flip the vertical!
1619 : !EOP
1620 : !------------------------------------------------------------------------------
1621 : !BOC
1622 : !
1623 : ! !LOCAL VARIABLES:
1624 : !
1625 : character(len=*), parameter :: subname = 'HCO_Grid_CAM2HCO_3D'
1626 : integer :: RC
1627 :
1628 : integer :: J, K
1629 :
1630 : ! (field, data, KS, KE, CS, CE)
1631 0 : call HCO_ESMF_Set3DCAM(CAM_3DFld, camArray, 1, LM, 1, my_CE)
1632 :
1633 : call ESMF_FieldRegrid(CAM_3DFld, HCO_3DFld, CAM2HCO_RouteHandle_3D, &
1634 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1635 0 : checkflag=.true., rc=RC)
1636 0 : ASSERT_(RC==ESMF_SUCCESS)
1637 :
1638 0 : call HCO_ESMF_Get3DField(HCO_3DFld, hcoArray, my_IS, my_IE, my_JS, my_JE, 1, LM, flip=.true.)
1639 :
1640 0 : end subroutine HCO_Grid_CAM2HCO_3D
1641 : !EOC
1642 : !------------------------------------------------------------------------------
1643 : ! Harmonized Emissions Component (HEMCO) !
1644 : !------------------------------------------------------------------------------
1645 : !BOP
1646 : !
1647 : ! !IROUTINE: HCO_Grid_HCO2CAM_2D
1648 : !
1649 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_2D regrids a HEMCO lat-lon grid
1650 : ! field (i,j) to CAM physics mesh (i).
1651 : !\\
1652 : !\\
1653 : ! !INTERFACE:
1654 : !
1655 0 : subroutine HCO_Grid_HCO2CAM_2D(hcoArray, camArray)
1656 : !
1657 : ! !USES:
1658 : !
1659 0 : use ESMF, only: ESMF_FieldRegrid
1660 : use ESMF, only: ESMF_TERMORDER_SRCSEQ
1661 :
1662 : ! MPI Properties from CAM
1663 : use cam_logfile, only: iulog
1664 : use spmd_utils, only: masterproc
1665 : !
1666 : ! !INPUT PARAMETERS:
1667 : !
1668 : real(r8), intent(in) :: hcoArray(my_IS:my_IE,my_JS:my_JE)
1669 : !
1670 : ! !OUTPUT PARAMETERS:
1671 : !
1672 : real(r8), intent(inout) :: camArray(1:my_CE) ! Col and lev start on 1
1673 : !
1674 : ! !REMARKS:
1675 : ! (1) Used for 2-D fields, which are 1-D in CAM.
1676 : !
1677 : ! !REVISION HISTORY:
1678 : ! 31 Mar 2020 - H.P. Lin - Initial version
1679 : !EOP
1680 : !------------------------------------------------------------------------------
1681 : !BOC
1682 : !
1683 : ! !LOCAL VARIABLES:
1684 : !
1685 : character(len=*), parameter :: subname = 'HCO_Grid_HCO2CAM_2D'
1686 : integer :: RC
1687 :
1688 0 : call HCO_ESMF_Set2DHCO(HCO_2DFld, hcoArray, my_IS, my_IE, my_JS, my_JE)
1689 :
1690 : call ESMF_FieldRegrid(HCO_2DFld, CAM_2DFld, HCO2CAM_RouteHandle_2D, &
1691 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1692 0 : checkflag=.true., rc=RC)
1693 0 : ASSERT_(RC==ESMF_SUCCESS)
1694 :
1695 : ! HCO_ESMF_Get1DField(field_in, data_out, IS, IE)
1696 : ! (Physics "1D" fields on mesh are actually "2D" data)
1697 0 : call HCO_ESMF_Get1DField(CAM_2DFld, camArray, 1, my_CE)
1698 :
1699 0 : end subroutine HCO_Grid_HCO2CAM_2D
1700 : !EOC
1701 : !------------------------------------------------------------------------------
1702 : ! Harmonized Emissions Component (HEMCO) !
1703 : !------------------------------------------------------------------------------
1704 : !BOP
1705 : !
1706 : ! !IROUTINE: HCO_Grid_CAM2HCO_2D
1707 : !
1708 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_2D regrids a HEMCO lat-lon grid
1709 : ! field (i,j,l) to CAM physics mesh (k,i).
1710 : !\\
1711 : !\\
1712 : ! !INTERFACE:
1713 : !
1714 0 : subroutine HCO_Grid_CAM2HCO_2D(camArray, hcoArray)
1715 : !
1716 : ! !USES:
1717 : !
1718 0 : use ESMF, only: ESMF_FieldRegrid, ESMF_TERMORDER_SRCSEQ
1719 :
1720 : ! MPI Properties from CAM
1721 : use cam_logfile, only: iulog
1722 : use spmd_utils, only: masterproc
1723 : !
1724 : ! !INPUT PARAMETERS:
1725 : !
1726 : real(r8), intent(in) :: camArray(1:my_CE) ! Col and lev start on 1
1727 : !
1728 : ! !OUTPUT PARAMETERS:
1729 : !
1730 : real(r8), intent(inout) :: hcoArray(my_IS:my_IE,my_JS:my_JE)
1731 : !
1732 : ! !REVISION HISTORY:
1733 : ! 31 Mar 2020 - H.P. Lin - Initial version
1734 : !EOP
1735 : !------------------------------------------------------------------------------
1736 : !BOC
1737 : !
1738 : ! !LOCAL VARIABLES:
1739 : !
1740 : character(len=*), parameter :: subname = 'HCO_Grid_CAM2HCO_2D'
1741 : integer :: RC
1742 :
1743 : integer :: J
1744 :
1745 : ! HCO_ESMF_Set2DCAM(field, data, CS, CE)
1746 0 : call HCO_ESMF_Set2DCAM(CAM_2DFld, camArray, 1, my_CE)
1747 :
1748 : call ESMF_FieldRegrid(CAM_2DFld, HCO_2DFld, CAM2HCO_RouteHandle_2D, &
1749 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1750 0 : checkflag=.true., rc=RC)
1751 0 : ASSERT_(RC==ESMF_SUCCESS)
1752 :
1753 0 : call HCO_ESMF_Get2DField(HCO_2DFld, hcoArray, my_IS, my_IE, my_JS, my_JE)
1754 :
1755 : ! Kludge for periodic point
1756 : ! It seems like the last point for the PET in the x-edge direction is messed up,
1757 : ! because it is the "periodic" point wrapping around the globe. This part is
1758 : ! not smooth and needs to be manually extrapolated.
1759 : !
1760 : ! This is a kludge as we want to revisit the regrid mechanism later.
1761 : ! For now fix it by copying the edge, not IDAVG
1762 0 : if(my_IE .eq. IM) then
1763 0 : do J = my_JS, my_JE
1764 0 : if(hcoArray(my_IE, J) .le. 0.000001_r8) then
1765 0 : hcoArray(my_IE, J) = hcoArray(my_IE-1, J)
1766 : endif
1767 : enddo
1768 : endif
1769 :
1770 0 : end subroutine HCO_Grid_CAM2HCO_2D
1771 : !EOC
1772 : !------------------------------------------------------------------------------
1773 : ! Harmonized Emissions Component (HEMCO) !
1774 : !------------------------------------------------------------------------------
1775 : !BOP
1776 : !
1777 : ! !IROUTINE: HCO_ESMF_Set2DHCO
1778 : !
1779 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set2DHCO sets values of a ESMF field on
1780 : ! the HEMCO lat-lon grid. (Internal use)
1781 : !\\
1782 : !\\
1783 : ! !INTERFACE:
1784 : !
1785 0 : subroutine HCO_ESMF_Set2DHCO(field, data, IS, IE, JS, JE)
1786 : !
1787 : ! !USES:
1788 : !
1789 0 : use ESMF, only: ESMF_FieldGet
1790 : !
1791 : ! !INPUT PARAMETERS:
1792 : !
1793 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1794 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
1795 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
1796 : real(r8), intent(in) :: data(IS:IE, JS:JE) ! Data to write in
1797 : ! !REMARKS:
1798 : ! Note there might need to be handling for periodic points?
1799 : !
1800 : ! !REVISION HISTORY:
1801 : ! 24 Feb 2020 - H.P. Lin - Initial version
1802 : !EOP
1803 : !------------------------------------------------------------------------------
1804 : !BOC
1805 : !
1806 : ! !LOCAL VARIABLES:
1807 : !
1808 : character(len=*), parameter :: subname = 'HCO_ESMF_Set2DHCO'
1809 : integer :: I, J
1810 : integer :: RC
1811 : integer :: lbnd(2), ubnd(2)
1812 0 : real(r8), pointer :: fptr(:,:)
1813 :
1814 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1815 : computationalLBound=lbnd, &
1816 0 : computationalUBound=ubnd, rc=RC)
1817 0 : ASSERT_(RC==ESMF_SUCCESS)
1818 0 : fptr(:,:) = 0.0_r8 ! Arbitrary missval
1819 0 : do J = lbnd(2), ubnd(2)
1820 0 : do I = lbnd(1), ubnd(1)
1821 0 : fptr(I,J) = data(I,J)
1822 : enddo
1823 : enddo
1824 :
1825 0 : end subroutine HCO_ESMF_Set2DHCO
1826 : !EOC
1827 : !------------------------------------------------------------------------------
1828 : ! Harmonized Emissions Component (HEMCO) !
1829 : !------------------------------------------------------------------------------
1830 : !BOP
1831 : !
1832 : ! !IROUTINE: HCO_ESMF_Set3DHCO
1833 : !
1834 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set3DHCO sets values of a ESMF field on
1835 : ! the HEMCO lat-lon grid. (Internal use)
1836 : !\\
1837 : !\\
1838 : ! !INTERFACE:
1839 : !
1840 0 : subroutine HCO_ESMF_Set3DHCO(field, data, IS, IE, JS, JE, KS, KE)
1841 : !
1842 : ! !USES:
1843 : !
1844 0 : use ESMF, only: ESMF_FieldGet
1845 : !
1846 : ! !INPUT PARAMETERS:
1847 : !
1848 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1849 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
1850 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
1851 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
1852 : real(r8), intent(in) :: data(IS:IE, JS:JE, KS:KE) ! Data to write in
1853 : ! !REMARKS:
1854 : ! Note there might need to be handling for periodic points?
1855 : !
1856 : ! !REVISION HISTORY:
1857 : ! 24 Feb 2020 - H.P. Lin - Initial version
1858 : !EOP
1859 : !------------------------------------------------------------------------------
1860 : !BOC
1861 : !
1862 : ! !LOCAL VARIABLES:
1863 : !
1864 : character(len=*), parameter :: subname = 'HCO_ESMF_Set3DHCO'
1865 : integer :: I, J, K
1866 : integer :: RC
1867 : integer :: lbnd(3), ubnd(3)
1868 0 : real(r8), pointer :: fptr(:,:,:)
1869 :
1870 0 : call ESMF_LogWrite("In HCO_ESMF_Set3DHCO before ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1871 :
1872 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1873 : computationalLBound=lbnd, &
1874 0 : computationalUBound=ubnd, rc=RC)
1875 0 : ASSERT_(RC==ESMF_SUCCESS)
1876 :
1877 0 : call ESMF_LogWrite("In HCO_ESMF_Set3DHCO after ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1878 0 : fptr(:,:,:) = 0.0_r8 ! Arbitrary missval
1879 :
1880 : ! write(6,*) "In HCO_ESMF_Set3DHCO IS,IE,JS,JE,KS,KE ESMF bnds", lbnd(1),ubnd(1),lbnd(2),ubnd(2),lbnd(3),ubnd(3)
1881 0 : do K = lbnd(3), ubnd(3)
1882 0 : do J = lbnd(2), ubnd(2)
1883 0 : do I = lbnd(1), ubnd(1)
1884 0 : fptr(I,J,K) = data(I,J,K)
1885 : enddo
1886 : enddo
1887 : enddo
1888 :
1889 0 : end subroutine HCO_ESMF_Set3DHCO
1890 : !EOC
1891 : !------------------------------------------------------------------------------
1892 : ! Harmonized Emissions Component (HEMCO) !
1893 : !------------------------------------------------------------------------------
1894 : !BOP
1895 : !
1896 : ! !IROUTINE: HCO_ESMF_Set2DCAM
1897 : !
1898 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set2DCAM sets values of a ESMF field on
1899 : ! the physics mesh. (Internal use)
1900 : !\\
1901 : !\\
1902 : ! !INTERFACE:
1903 : !
1904 0 : subroutine HCO_ESMF_Set2DCAM(field, data, CS, CE)
1905 : !
1906 : ! !USES:
1907 : !
1908 0 : use ESMF, only: ESMF_FieldGet
1909 : !
1910 : ! !INPUT PARAMETERS:
1911 : !
1912 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1913 : integer, intent(in) :: CS, CE ! Start and end indices (chunk)
1914 : real(r8), intent(in) :: data(CS:CE) ! Data to write in
1915 : ! !REMARKS:
1916 : ! Remember that "2-D" fields in physics chunk is actually 3-D (k,i), and 2-D
1917 : ! is actually just column-level data. The target is a mesh so in ESMF it is
1918 : ! the same.
1919 : !
1920 : ! !REVISION HISTORY:
1921 : ! 23 Feb 2020 - H.P. Lin - Initial version
1922 : !EOP
1923 : !------------------------------------------------------------------------------
1924 : !BOC
1925 : !
1926 : ! !LOCAL VARIABLES:
1927 : !
1928 : character(len=*), parameter :: subname = 'HCO_ESMF_Set2DCAM'
1929 : integer :: I
1930 : integer :: RC
1931 : integer :: lbnd(1), ubnd(1)
1932 0 : real(r8), pointer :: fptr(:)
1933 :
1934 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1935 : computationalLBound=lbnd, &
1936 0 : computationalUBound=ubnd, rc=RC)
1937 0 : ASSERT_(RC==ESMF_SUCCESS)
1938 0 : fptr(:) = 0.0_r8 ! Arbitrary missval
1939 0 : do I = lbnd(1), ubnd(1)
1940 0 : fptr(i) = data(i)
1941 : enddo
1942 :
1943 0 : end subroutine HCO_ESMF_Set2DCAM
1944 : !EOC
1945 : !------------------------------------------------------------------------------
1946 : ! Harmonized Emissions Component (HEMCO) !
1947 : !------------------------------------------------------------------------------
1948 : !BOP
1949 : !
1950 : ! !IROUTINE: HCO_ESMF_Set3DCAM
1951 : !
1952 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set3DCAM sets values of a ESMF field on
1953 : ! the physics mesh. (Internal use)
1954 : !\\
1955 : !\\
1956 : ! !INTERFACE:
1957 : !
1958 0 : subroutine HCO_ESMF_Set3DCAM(field, data, KS, KE, CS, CE)
1959 : !
1960 : ! !USES:
1961 : !
1962 0 : use ESMF, only: ESMF_FieldGet
1963 : !
1964 : ! !INPUT PARAMETERS:
1965 : !
1966 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1967 : integer, intent(in) :: CS, CE ! Start and end indices (chunk)
1968 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
1969 : real(r8), intent(in) :: data(KS:KE, CS:CE) ! Data to write in
1970 : ! !REMARKS:
1971 : ! This is the actual 3-D data routine, pointers (i,k) on a mesh
1972 : !
1973 : ! !REVISION HISTORY:
1974 : ! 23 Feb 2020 - H.P. Lin - Initial version
1975 : !EOP
1976 : !------------------------------------------------------------------------------
1977 : !BOC
1978 : !
1979 : ! !LOCAL VARIABLES:
1980 : !
1981 : character(len=*), parameter :: subname = 'HCO_ESMF_Set3DCAM'
1982 : integer :: I, K
1983 : integer :: RC
1984 : integer :: lbnd(2), ubnd(2)
1985 0 : real(r8), pointer :: fptr(:,:)
1986 :
1987 : ! call ESMF_LogWrite("In HCO_ESMF_Set3DCAM before ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1988 :
1989 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1990 : computationalLBound=lbnd, &
1991 0 : computationalUBound=ubnd, rc=RC)
1992 0 : ASSERT_(RC==ESMF_SUCCESS)
1993 :
1994 :
1995 : ! call ESMF_LogWrite("In HCO_ESMF_Set3DCAM after ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1996 0 : fptr(:,:) = 0.0_r8 ! Arbitrary missval
1997 0 : do I = lbnd(2), ubnd(2)
1998 0 : do K = lbnd(1), ubnd(1)
1999 0 : fptr(k,i) = data(k,i)
2000 : enddo
2001 : enddo
2002 :
2003 0 : end subroutine HCO_ESMF_Set3DCAM
2004 : !EOC
2005 : !------------------------------------------------------------------------------
2006 : ! Harmonized Emissions Component (HEMCO) !
2007 : !------------------------------------------------------------------------------
2008 : !BOP
2009 : !
2010 : ! !IROUTINE: HCO_ESMF_Get1DField
2011 : !
2012 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get1DField gets a pointer to an 1-D ESMF
2013 : ! field.
2014 : !\\
2015 : !\\
2016 : ! !INTERFACE:
2017 : !
2018 0 : subroutine HCO_ESMF_Get1DField(field_in, data_out, IS, IE)
2019 : !
2020 : ! !USES:
2021 : !
2022 0 : use ESMF, only: ESMF_FieldGet
2023 : !
2024 : ! !INPUT PARAMETERS:
2025 : !
2026 : type(ESMF_Field), intent(in) :: field_in
2027 : integer, intent(in) :: IS, IE ! Start and end indices
2028 : !
2029 : ! !OUTPUT PARAMETERS:
2030 : !
2031 : real(r8), intent(out) :: data_out(IS:IE)
2032 : !
2033 : ! !REVISION HISTORY:
2034 : ! 23 Feb 2020 - H.P. Lin - Initial version
2035 : !EOP
2036 : !------------------------------------------------------------------------------
2037 : !BOC
2038 : !
2039 : ! !LOCAL VARIABLES:
2040 : !
2041 : character(len=*), parameter :: subname = 'HCO_ESMF_Get1DField'
2042 0 : real(r8), pointer :: fptr(:)
2043 : integer :: RC
2044 :
2045 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2046 0 : ASSERT_(RC==ESMF_SUCCESS)
2047 0 : data_out(:) = fptr(:)
2048 :
2049 0 : end subroutine HCO_ESMF_Get1DField
2050 : !EOC
2051 : !------------------------------------------------------------------------------
2052 : ! Harmonized Emissions Component (HEMCO) !
2053 : !------------------------------------------------------------------------------
2054 : !BOP
2055 : !
2056 : ! !IROUTINE: HCO_ESMF_Get2DField
2057 : !
2058 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get2DField gets a pointer to an 2-D ESMF
2059 : ! field.
2060 : !\\
2061 : !\\
2062 : ! !INTERFACE:
2063 : !
2064 0 : subroutine HCO_ESMF_Get2DField(field_in, data_out, IS, IE, JS, JE, flip)
2065 : !
2066 : ! !USES:
2067 : !
2068 0 : use ESMF, only: ESMF_FieldGet
2069 : !
2070 : ! !INPUT PARAMETERS:
2071 : !
2072 : type(ESMF_Field), intent(in) :: field_in
2073 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
2074 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
2075 : logical, optional, intent(in) :: flip ! Flip the FIRST dimension? Presence will flip
2076 : !
2077 : ! !OUTPUT PARAMETERS:
2078 : !
2079 : real(r8), intent(out) :: data_out(IS:IE, JS:JE)
2080 : !
2081 : ! !REMARKS:
2082 : ! If the flip argument is PROVIDED, then the FIRST dimension (IS:IE) will be flipped.
2083 : ! This is used when retrieving a CAM field regridded from HEMCO, when you have to flip
2084 : ! the vertical to fit correctly. The dimensions of a CAM 3d data is a 2d chunked array
2085 : ! with indices (k, i) so... the first dimension has to be flipped instead.
2086 : !
2087 : ! YES - THE FLIP ARGUMENT IS NOT READ. ITS MERE PRESENCE WILL TRIGGER .TRUE. FOR FLIPPING.
2088 : ! THIS IS INTENDED BEHAVIOR TO CIRCUMVENT AN INTEL COMPILER BUG.
2089 : !
2090 : ! !REVISION HISTORY:
2091 : ! 23 Feb 2020 - H.P. Lin - Initial version
2092 : ! 30 May 2020 - H.P. Lin - Add flip parameter
2093 : ! 03 Jan 2021 - H.P. Lin - I hate compiler bugs with PRESENT var
2094 : !EOP
2095 : !------------------------------------------------------------------------------
2096 : !BOC
2097 : !
2098 : ! !LOCAL VARIABLES:
2099 : !
2100 : character(len=*), parameter :: subname = 'HCO_ESMF_Get2DField'
2101 0 : real(r8), pointer :: fptr(:,:)
2102 : integer :: RC
2103 : logical :: flip1d = .false.
2104 :
2105 0 : if(present(flip)) then
2106 0 : flip1d = .true.
2107 : else
2108 0 : flip1d = .false.
2109 : endif
2110 :
2111 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2112 0 : ASSERT_(RC==ESMF_SUCCESS)
2113 :
2114 0 : if(flip1d) then
2115 0 : data_out(IS:IE:1, :) = fptr(IE:IS:-1, :)
2116 : else
2117 0 : data_out(:,:) = fptr(:,:)
2118 : endif
2119 :
2120 0 : end subroutine HCO_ESMF_Get2DField
2121 : !EOC
2122 : !------------------------------------------------------------------------------
2123 : ! Harmonized Emissions Component (HEMCO) !
2124 : !------------------------------------------------------------------------------
2125 : !BOP
2126 : !
2127 : ! !IROUTINE: HCO_ESMF_Get3DField
2128 : !
2129 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get3DField gets a pointer to an 3-D ESMF
2130 : ! field.
2131 : !\\
2132 : !\\
2133 : ! !INTERFACE:
2134 : !
2135 0 : subroutine HCO_ESMF_Get3DField(field_in, data_out, IS, IE, JS, JE, KS, KE, flip)
2136 : !
2137 : ! !USES:
2138 : !
2139 0 : use ESMF, only: ESMF_FieldGet
2140 : !
2141 : ! !INPUT PARAMETERS:
2142 : !
2143 : type(ESMF_Field), intent(in) :: field_in
2144 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
2145 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
2146 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
2147 : logical, optional, intent(in) :: flip ! Flip the 3rd dimension?
2148 : !
2149 : ! !OUTPUT PARAMETERS:
2150 : !
2151 : real(r8), intent(out) :: data_out(IS:IE, JS:JE, KS:KE)
2152 : !
2153 : ! !REMARKS:
2154 : ! If the flip argument is set to .true., then the third dimension (KS:KE) will be flipped.
2155 : ! This is used when retrieving a HEMCO field regridded from CAM, when you have to flip
2156 : ! the vertical to fit correctly.
2157 : !
2158 : ! YES - THE FLIP ARGUMENT IS NOT READ. ITS MERE PRESENCE WILL TRIGGER .TRUE. FOR FLIPPING.
2159 : ! THIS IS INTENDED BEHAVIOR TO CIRCUMVENT AN INTEL COMPILER BUG.
2160 : !
2161 : ! !REVISION HISTORY:
2162 : ! 23 Feb 2020 - H.P. Lin - Initial version
2163 : ! 30 May 2020 - H.P. Lin - Add flip parameter
2164 : !EOP
2165 : !------------------------------------------------------------------------------
2166 : !BOC
2167 : !
2168 : ! !LOCAL VARIABLES:
2169 : !
2170 : character(len=*), parameter :: subname = 'HCO_ESMF_Get3DField'
2171 0 : real(r8), pointer :: fptr(:,:,:)
2172 : integer :: RC
2173 : logical :: flip3d = .false.
2174 :
2175 0 : if(present(flip)) then
2176 0 : flip3d = .true.
2177 : else
2178 0 : flip3d = .false.
2179 : endif
2180 :
2181 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2182 0 : ASSERT_(RC==ESMF_SUCCESS)
2183 :
2184 0 : if(flip3d) then
2185 0 : data_out(:,:,KS:KE:1) = fptr(:,:,KE:KS:-1)
2186 : else
2187 0 : data_out(:,:,:) = fptr(:,:,:)
2188 : endif
2189 :
2190 0 : end subroutine HCO_ESMF_Get3DField
2191 : !EOC
2192 0 : end module hco_esmf_grid
|
