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(f14.6,1x) )' ) AP(1:LM+1)
438 0 : write( iulog, '(a)' )
439 0 : write( iulog, '( ''Bp '', /, 6(f14.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 0 : 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,*) "Committed HCO_Tasks decomposition:"
694 0 : write(iulog,*) "nPET = ", nPET, ", nPET_lon = ", nPET_lon, ", nPET_lat = ", nPET_lat
695 0 : do N = 0, nPET-1
696 : ! more info to print if needed
697 : !write(iulog,*) "PET", N, " ID", HCO_Tasks(N)%ID, " ID_I", HCO_Tasks(N)%ID_I, " ID_J", HCO_Tasks(N)%ID_J, &
698 : ! "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
699 : enddo
700 : endif
701 :
702 : !
703 : ! Just remember that my_IM, my_JM ... are your keys to generating
704 : ! the relevant met fields and passing to HEMCO.
705 : !
706 : ! Only HCO_ESMF_Grid should be aware of the entire grid.
707 : ! Everyone else should be just doing work on its subset indices.
708 : !
709 0 : RC = ESMF_SUCCESS
710 :
711 0 : end subroutine HCO_Grid_Init
712 : !EOC
713 : !------------------------------------------------------------------------------
714 : ! Harmonized Emissions Component (HEMCO) !
715 : !------------------------------------------------------------------------------
716 : !BOP
717 : !
718 : ! !IROUTINE: HCO_Grid_UpdateRegrid
719 : !
720 : ! !DESCRIPTION: Subroutine HCO\_Grid\_UpdateRegrid initializes or updates the
721 : ! regridding information used in the HEMCO_CESM interface.
722 : !\\
723 : !\\
724 : ! !INTERFACE:
725 : !
726 0 : subroutine HCO_Grid_UpdateRegrid( RC )
727 : !
728 : ! !USES:
729 : !
730 : ! MPI Properties from CAM
731 : ! Even though CAM's principle is that only spmd_utils uses MPI,
732 : ! ionos code uses MPI very liberally. Unfortunately we have to follow
733 : ! this example as spmd_utils does not provide many of the relevant
734 : ! functionality like the communicator split. But keep this in mind.
735 : use cam_logfile, only: iulog
736 : use spmd_utils, only: CAM_mpicom => mpicom, CAM_npes => npes
737 : use spmd_utils, only: MPI_SUCCESS
738 : use spmd_utils, only: iam
739 : use spmd_utils, only: masterproc
740 :
741 : use phys_grid, only: get_ncols_p
742 : use ppgrid, only: begchunk, endchunk
743 :
744 : use cam_instance, only: atm_id
745 :
746 : ! ESMF
747 : use ESMF, only: ESMF_FieldRegridStore
748 : use ESMF, only: ESMF_REGRIDMETHOD_BILINEAR, ESMF_REGRIDMETHOD_CONSERVE
749 : use ESMF, only: ESMF_POLEMETHOD_ALLAVG, ESMF_POLEMETHOD_NONE
750 : use ESMF, only: ESMF_EXTRAPMETHOD_NEAREST_IDAVG
751 :
752 : use ESMF, only: ESMF_FieldGet
753 : !
754 : ! !OUTPUT PARAMETERS:
755 : !
756 : integer, intent(out) :: RC
757 : !
758 : ! !REMARKS:
759 : ! This field will ONLY update if it recognizes a change in the CAM instance
760 : ! information, as determined by cam_instance::atm_id which is saved in the
761 : ! module's cam_last_atm_id field.
762 : !
763 : ! This allows the function HCO_Grid_UpdateRegrid be called both in init and run
764 : ! without performance / memory repercussions (hopefully...)
765 : !
766 : ! Maybe we can use CONSERVE_2ND order for better accuracy in the future. To be tested.
767 : !
768 : ! !REVISION HISTORY:
769 : ! 20 Feb 2020 - H.P. Lin - Initial version
770 : ! 23 Feb 2020 - H.P. Lin - Finalized regridding handles, bilinear for CAM-HCO and
771 : ! conservative (1st order) for HCO-CAM.
772 : ! 29 Oct 2020 - H.P. Lin - Changed to 1st-order conservative based on feedback from T. Fritz
773 : ! about negative numerical artifacts coming out of ESMF.
774 : !EOP
775 : !------------------------------------------------------------------------------
776 : !BOC
777 : !
778 : ! !LOCAL VARIABLES:
779 : !
780 : character(len=*), parameter :: subname = 'HCO_Grid_UpdateRegrid'
781 : integer :: smm_srctermproc, smm_pipelinedep
782 :
783 : integer :: chnk
784 :
785 : ! Debug only
786 : integer :: lbnd_hco(3), ubnd_hco(3) ! 3-d bounds of HCO field
787 : integer :: lbnd_cam(2), ubnd_cam(2) ! 3-d bounds of CAM field
788 : real(r8), pointer :: fptr(:,:,:) ! debug
789 : real(r8), pointer :: fptr_cam(:,:) ! debug
790 :
791 : ! Assume success
792 0 : RC = ESMF_SUCCESS
793 :
794 : ! Parameters for ESMF RouteHandle (taken from ionos interface)
795 0 : smm_srctermproc = 0
796 : ! smm_pipelinedep = -1 ! Accept auto-tuning of the pipeline depth
797 0 : smm_pipelinedep = 16 ! Prescribe pipeline depth for BFB consistency (hplin, 6/7/23)
798 :
799 : ! Check if we need to update
800 0 : if(cam_last_atm_id == atm_id) then
801 0 : if(masterproc) then
802 0 : write(iulog,'(A,I2)') "HEMCO_CESM: UpdateRegrid is already in atmospheric grid #", atm_id
803 : endif
804 0 : return
805 : else
806 0 : if(masterproc) then
807 0 : write(iulog,'(A,I2)') "HEMCO_CESM: UpdateRegrid now updating for atmospheric grid #", atm_id
808 : endif
809 : endif
810 :
811 0 : cam_last_atm_id = atm_id
812 :
813 : ! Create CAM physics mesh...
814 0 : call HCO_Grid_ESMF_CreateCAM(RC)
815 0 : ASSERT_(RC==ESMF_SUCCESS)
816 :
817 : ! How many columns in this task? my_CE
818 : ! Note that my_CE is LOCAL column index (1:my_CE) and is equivalent to "blksize"
819 : ! in the ionos coupling code. This was 144 in my testing.
820 : ! But this is NOT ppgrid::pcols, which was 16.
821 0 : my_CE = 0
822 0 : do chnk = begchunk, endchunk
823 0 : my_CE = my_CE + get_ncols_p(chnk)
824 : enddo
825 :
826 : ! Create HEMCO grid in ESMF format
827 0 : call HCO_Grid_ESMF_CreateHCO(RC)
828 0 : ASSERT_(RC==ESMF_SUCCESS)
829 :
830 : ! Create empty fields on the HEMCO grid and CAM physics mesh
831 : ! used for later regridding
832 :
833 : ! FIXME: Destroy fields before creating to prevent memory leak? ESMF_Destroy
834 : ! requires the field to be present (so no silent destruction) (hplin, 2/21/20)
835 :
836 0 : call HCO_Grid_ESMF_CreateCAMField(CAM_2DFld, CAM_PhysMesh, 'HCO_PHYS_2DFLD', 0, RC)
837 0 : call HCO_Grid_ESMF_CreateCAMField(CAM_3DFld, CAM_PhysMesh, 'HCO_PHYS_3DFLD', LM, RC)
838 :
839 0 : call HCO_Grid_ESMF_CreateHCOField(HCO_2DFld, HCO_Grid, 'HCO_2DFLD', 0, RC)
840 0 : call HCO_Grid_ESMF_CreateHCOField(HCO_3DFld, HCO_Grid, 'HCO_3DFLD', LM, RC)
841 0 : ASSERT_(RC==ESMF_SUCCESS)
842 :
843 : ! if(masterproc) then
844 : ! write(iulog,*) ">> HEMCO: HCO_PHYS and HCO_ fields initialized successfully"
845 : ! call ESMF_FieldGet(HCO_3DFld, localDE=0, farrayPtr=fptr, &
846 : ! computationalLBound=lbnd_hco, &
847 : ! computationalUBound=ubnd_hco, rc=RC)
848 : ! write(iulog,*) ">> HEMCO: Debug HCO Field: lbnd = (", lbnd_hco, "), ", my_IS, my_IE, "ubnd = (", ubnd_hco, ")", my_JS, my_JE
849 :
850 : ! call ESMF_FieldGet(CAM_3DFld, localDE=0, farrayPtr=fptr_cam, &
851 : ! computationalLBound=lbnd_cam, &
852 : ! computationalUBound=ubnd_cam, rc=RC)
853 : ! write(iulog,*) ">> HEMCO: Debug CAM Field: lbnd = (", lbnd_cam, "), ubnd = (", ubnd_cam, ")", my_CE
854 : ! endif
855 :
856 : ! CAM -> HCO 2-D
857 : call ESMF_FieldRegridStore( &
858 : srcField=CAM_2DFld, dstField=HCO_2DFld, &
859 : regridMethod=ESMF_REGRIDMETHOD_BILINEAR, &
860 : poleMethod=ESMF_POLEMETHOD_ALLAVG, &
861 : extrapMethod=ESMF_EXTRAPMETHOD_NEAREST_IDAVG, &
862 : routeHandle=CAM2HCO_RouteHandle_2D, &
863 : !srcTermProcessing=smm_srctermproc, &
864 0 : pipelineDepth=smm_pipelinedep, rc=RC)
865 0 : ASSERT_(RC==ESMF_SUCCESS)
866 :
867 0 : if(masterproc) then
868 0 : write(iulog,*) ">> After FieldRegridStore, CAM2D->HCO2D, pipeline", smm_pipelinedep
869 : endif
870 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
871 :
872 : ! Create and store ESMF route handles for regridding CAM <-> HCO
873 : ! CAM -> HCO 3-D
874 : call ESMF_FieldRegridStore( &
875 : srcField=CAM_3DFld, dstField=HCO_3DFld, &
876 : regridMethod=ESMF_REGRIDMETHOD_BILINEAR, &
877 : poleMethod=ESMF_POLEMETHOD_ALLAVG, &
878 : extrapMethod=ESMF_EXTRAPMETHOD_NEAREST_IDAVG, &
879 : routeHandle=CAM2HCO_RouteHandle_3D, &
880 : !srcTermProcessing=smm_srctermproc, &
881 0 : pipelineDepth=smm_pipelinedep, rc=RC)
882 0 : ASSERT_(RC==ESMF_SUCCESS)
883 :
884 0 : if(masterproc) then
885 0 : write(iulog,*) ">> After FieldRegridStore, CAM3D->HCO3D, pipeline", smm_pipelinedep
886 : endif
887 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
888 :
889 : ! HCO -> CAM 2-D
890 : call ESMF_FieldRegridStore( &
891 : srcField=HCO_2DFld, dstField=CAM_2DFld, &
892 : regridMethod=ESMF_REGRIDMETHOD_CONSERVE , &
893 : poleMethod=ESMF_POLEMETHOD_NONE, &
894 : routeHandle=HCO2CAM_RouteHandle_2D, &
895 : !srcTermProcessing=smm_srctermproc, &
896 0 : pipelineDepth=smm_pipelinedep, rc=RC)
897 0 : ASSERT_(RC==ESMF_SUCCESS)
898 :
899 0 : if(masterproc) then
900 0 : write(iulog,*) ">> After FieldRegridStore, HCO2D->CAM2D, pipeline", smm_pipelinedep
901 : endif
902 : ! smm_pipelinedep = -1 ! Only replace with -1 to accept auto-tuning of smm pipeline depth.
903 :
904 : ! HCO -> CAM 3-D
905 : ! 3-D conserv. regridding cannot be done on a stagger other than center
906 : ! (as of ESMF 8.0.0 in ESMF_FieldRegrid.F90::993)
907 : call ESMF_FieldRegridStore( &
908 : srcField=HCO_3DFld, dstField=CAM_3DFld, &
909 : regridMethod=ESMF_REGRIDMETHOD_CONSERVE , &
910 : poleMethod=ESMF_POLEMETHOD_NONE, &
911 : routeHandle=HCO2CAM_RouteHandle_3D, &
912 : !srcTermProcessing=smm_srctermproc, &
913 0 : pipelineDepth=smm_pipelinedep, rc=RC)
914 0 : ASSERT_(RC==ESMF_SUCCESS)
915 :
916 0 : if(masterproc) then
917 0 : write(iulog,*) ">> After FieldRegridStore, HCO3D->CAM3D, pipeline", smm_pipelinedep
918 : endif
919 :
920 0 : if(masterproc) then
921 0 : write(iulog,*) "HEMCO_CESM: FieldRegridStore four-way complete", atm_id
922 : endif
923 :
924 : ! Finished updating regrid routines
925 0 : end subroutine HCO_Grid_UpdateRegrid
926 : !EOC
927 : !------------------------------------------------------------------------------
928 : ! Harmonized Emissions Component (HEMCO) !
929 : !------------------------------------------------------------------------------
930 : !BOP
931 : !
932 : ! !IROUTINE: HCO_Grid_ESMF_CreateCAM
933 : !
934 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateCAM creates the physics mesh
935 : ! from CAM and stores in the ESMF state for regridding.
936 : !\\
937 : !\\
938 : ! !INTERFACE:
939 : !
940 0 : subroutine HCO_Grid_ESMF_CreateCAM( RC )
941 : !
942 : ! !USES:
943 : !
944 : ! MPI Properties from CAM
945 0 : use cam_logfile, only: iulog
946 : use spmd_utils, only: masterproc
947 :
948 : ! Phys constants
949 : use shr_const_mod, only: pi => shr_const_pi
950 :
951 : ! Grid properties in CAM
952 : use cam_instance, only: inst_name
953 : use phys_control, only: phys_getopts
954 : use phys_grid, only: get_ncols_p, get_gcol_p, get_rlon_all_p, get_rlat_all_p
955 : use ppgrid, only: begchunk, endchunk
956 : use ppgrid, only: pcols ! # of col max in chunks
957 :
958 : ! ESMF
959 : use ESMF, only: ESMF_DistGridCreate, ESMF_MeshCreate
960 : use ESMF, only: ESMF_MeshGet
961 : use ESMF, only: ESMF_FILEFORMAT_ESMFMESH
962 : use ESMF, only: ESMF_MeshIsCreated, ESMF_MeshDestroy
963 : !
964 : ! !OUTPUT PARAMETERS:
965 : !
966 : integer, intent(out) :: RC
967 : !
968 : ! !REMARKS:
969 : !
970 : ! !REVISION HISTORY:
971 : ! 11 Feb 2020 - H.P. Lin - Initial version
972 : !EOP
973 : !------------------------------------------------------------------------------
974 : !BOC
975 : !
976 : ! !LOCAL VARIABLES:
977 : !
978 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateCAM'
979 :
980 : ! For allocation of the distGrid and Mesh
981 : integer :: ncols
982 : integer :: chnk, col, dindex
983 0 : integer, allocatable :: decomp(:)
984 : integer :: col_total
985 : character(len=256) :: grid_file
986 :
987 : ! For verification of the mesh
988 : integer :: spatialDim
989 : integer :: numOwnedElements
990 0 : real(r8), pointer :: ownedElemCoords(:)
991 0 : real(r8), pointer :: latCAM(:), latMesh(:)
992 0 : real(r8), pointer :: lonCAM(:), lonMesh(:)
993 : real(r8) :: latCAM_R(pcols) ! array of chunk lat
994 : real(r8) :: lonCAM_R(pcols) ! array of chunk long
995 :
996 : integer :: i, c, n
997 : real(r8), parameter :: radtodeg = 180.0_r8/pi
998 :
999 : ! Assume success
1000 0 : RC = ESMF_SUCCESS
1001 :
1002 : !-----------------------------------------------------------------------
1003 : ! Compute the CAM_DistGrid and CAM_PhysMesh
1004 : !-----------------------------------------------------------------------
1005 : ! Get the physics grid information
1006 0 : call phys_getopts(physics_grid_out=grid_file)
1007 :
1008 0 : if(masterproc) then
1009 0 : write(iulog,*) "physics_grid_out=", grid_file
1010 : endif
1011 :
1012 : ! Compute local decomposition (global variable in-module)
1013 0 : col_total = 0 ! Sum of columns on this PET in all chunks
1014 0 : do chnk = begchunk, endchunk
1015 0 : col_total = col_total + get_ncols_p(chnk)
1016 : enddo
1017 0 : allocate(decomp(col_total))
1018 0 : allocate(lonCAM(col_total))
1019 0 : allocate(latCAM(col_total)) ! the _R variants are already allocated to pcols
1020 :
1021 : ! ...and also attach to the loop computing lat and lons for CAM physics mesh
1022 : ! to be used later.
1023 0 : dindex = 0
1024 0 : do chnk = begchunk, endchunk
1025 0 : ncols = get_ncols_p(chnk)
1026 :
1027 : ! Get [rad] lat and lons
1028 0 : call get_rlon_all_p(chnk, ncols, lonCAM_R)
1029 0 : call get_rlat_all_p(chnk, ncols, latCAM_R)
1030 :
1031 0 : do col = 1, ncols
1032 0 : dindex = dindex + 1
1033 0 : decomp(dindex) = get_gcol_p(chnk, col)
1034 :
1035 0 : lonCAM(dindex) = lonCAM_R(col) * radtodeg
1036 0 : latCAM(dindex) = latCAM_R(col) * radtodeg
1037 : enddo
1038 : enddo
1039 :
1040 : ! Build the 2D field CAM DistGrid based on the physics decomposition
1041 0 : CAM_DistGrid = ESMF_DistGridCreate(arbSeqIndexList=decomp, rc=RC)
1042 0 : ASSERT_(RC==ESMF_SUCCESS)
1043 :
1044 : ! Release memory if any is being taken, to avoid leakage
1045 0 : if(ESMF_MeshIsCreated(CAM_PhysMesh)) then
1046 0 : call ESMF_MeshDestroy(CAM_PhysMesh)
1047 : endif
1048 :
1049 : ! Create the physics decomposition ESMF mesh
1050 : CAM_PhysMesh = ESMF_MeshCreate(trim(grid_file), ESMF_FILEFORMAT_ESMFMESH, &
1051 0 : elementDistGrid=CAM_DistGrid, rc=RC)
1052 0 : ASSERT_(RC==ESMF_SUCCESS)
1053 :
1054 : !-----------------------------------------------------------------------
1055 : ! Validate mesh coordinates against model physics column coords.
1056 : !-----------------------------------------------------------------------
1057 : ! (From edyn_esmf::edyn_create_physmesh)
1058 : call ESMF_MeshGet(CAM_PhysMesh, spatialDim=spatialDim, &
1059 0 : numOwnedElements=numOwnedElements, rc=RC)
1060 0 : ASSERT_(RC==ESMF_SUCCESS)
1061 :
1062 0 : if(numOwnedElements /= col_total) then
1063 0 : write(iulog,*) "HEMCO: ESMF_MeshGet numOwnedElements =", numOwnedElements, &
1064 0 : "col_total =", col_total, " MISMATCH! Aborting"
1065 0 : ASSERT_(.false.)
1066 : endif
1067 :
1068 : ! Coords for the CAM_PhysMesh
1069 0 : allocate(ownedElemCoords(spatialDim * numOwnedElements))
1070 0 : allocate(lonMesh(col_total), latMesh(col_total))
1071 :
1072 0 : call ESMF_MeshGet(CAM_PhysMesh, ownedElemCoords=ownedElemCoords, rc=RC)
1073 :
1074 0 : do n = 1, col_total
1075 0 : lonMesh(n) = ownedElemCoords(2*n-1)
1076 0 : latMesh(n) = ownedElemCoords(2*n)
1077 : enddo
1078 :
1079 : ! Error check coordinates
1080 0 : do n = 1, col_total
1081 0 : if(abs(lonMesh(n) - lonCAM(n)) > 0.000001_r8) then
1082 0 : if((abs(lonMesh(n) - lonCAM(n)) > 360.000001_r8) .or. &
1083 0 : (abs(lonMesh(n) - lonCAM(n)) < 359.99999_r8)) then
1084 0 : write(6,*) "HEMCO: ESMF_MeshGet VERIFY fail! n, lonMesh, lonCAM, delta"
1085 0 : write(6,*) n, lonMesh(n), lonCAM(n), abs(lonMesh(n)-lonCAM(n))
1086 0 : ASSERT_(.false.)
1087 : endif
1088 : endif
1089 :
1090 0 : if(abs(latMesh(n) - latCAM(n)) > 0.000001_r8) then
1091 0 : if(.not. ((abs(latCAM(n)) > 88.0_r8) .and. (abs(latMesh(n)) > 88.0_r8))) then
1092 0 : write(6,*) "HEMCO: ESMF_MeshGet VERIFY fail! n, latmesh, latCAM, delta"
1093 0 : write(6,*) n, latMesh(n), latCAM(n), abs(latMesh(n)-latCAM(n))
1094 0 : ASSERT_(.false.)
1095 : endif
1096 : endif
1097 : enddo
1098 :
1099 : ! Ready to go
1100 0 : if(masterproc) then
1101 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateCAM ok, dim'l = ", col_total
1102 : endif
1103 :
1104 : ! Free memory
1105 0 : deallocate(ownedElemCoords)
1106 0 : deallocate(lonCAM, lonMesh)
1107 0 : deallocate(latCAM, latMesh)
1108 0 : deallocate(decomp)
1109 :
1110 0 : end subroutine HCO_Grid_ESMF_CreateCAM
1111 : !EOC
1112 : !------------------------------------------------------------------------------
1113 : ! Harmonized Emissions Component (HEMCO) !
1114 : !------------------------------------------------------------------------------
1115 : !BOP
1116 : !
1117 : ! !IROUTINE: HCO_Grid_ESMF_CreateHCO
1118 : !
1119 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateHCO creates the HEMCO grid
1120 : ! in center and corner staggering modes in ESMF_Grid format,
1121 : ! and stores in the ESMF state for regridding.
1122 : !\\
1123 : !\\
1124 : ! !INTERFACE:
1125 : !
1126 0 : subroutine HCO_Grid_ESMF_CreateHCO( RC )
1127 : !
1128 : ! !USES:
1129 : !
1130 : ! MPI Properties from CAM
1131 0 : use cam_logfile, only: iulog
1132 : use spmd_utils, only: masterproc
1133 :
1134 : ! ESMF
1135 : use ESMF, only: ESMF_GridCreate1PeriDim, ESMF_INDEX_GLOBAL
1136 : use ESMF, only: ESMF_STAGGERLOC_CENTER, ESMF_STAGGERLOC_CORNER
1137 : use ESMF, only: ESMF_GridAddCoord, ESMF_GridGetCoord
1138 : !
1139 : ! !OUTPUT PARAMETERS:
1140 : !
1141 : integer, intent(out) :: RC
1142 : !
1143 : ! !REMARKS:
1144 : ! We initialize TWO coordinates here for the HEMCO grid. Both the center and corner
1145 : ! staggering information needs to be added to the grid, or ESMF_FieldRegridStore will
1146 : ! fail. It is a rather weird implementation, as ESMF_FieldRegridStore only accepts
1147 : ! the coordinate in CENTER staggering, but the conversion to mesh is in CORNER.
1148 : !
1149 : ! So both coordinates need to be specified even though (presumably?) the center field
1150 : ! is regridded once the route handle is generated. This remains to be seen but at this
1151 : ! point only the following "duplicated" code is the correct implementation for
1152 : ! conservative regridding handles to be generated properly.
1153 : !
1154 : ! Roughly 9 hours of debugging and reading the ESMF documentation and code
1155 : ! were wasted here. (hplin, 2/23/20)
1156 : !
1157 : ! !REVISION HISTORY:
1158 : ! 20 Feb 2020 - H.P. Lin - Initial version
1159 : ! 22 Feb 2020 - H.P. Lin - Support curvilinear in ESMF format (note Map_A2A)
1160 : !EOP
1161 : !------------------------------------------------------------------------------
1162 : !BOC
1163 : !
1164 : ! !LOCAL VARIABLES:
1165 : !
1166 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateHCO'
1167 :
1168 : integer :: i, j, n
1169 : integer :: lbnd(2), ubnd(2)
1170 0 : integer :: nlons_task(nPET_lon) ! # lons per task
1171 0 : integer :: nlats_task(nPET_lat) ! # lats per task
1172 0 : real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:)
1173 0 : real(ESMF_KIND_R8), pointer :: coordX_E(:,:), coordY_E(:,:)
1174 :
1175 : !-----------------------------------------------------------------------
1176 : ! Task distribution for ESMF grid
1177 : !-----------------------------------------------------------------------
1178 0 : do i = 1, nPET_lon
1179 0 : loop: do n = 0, nPET-1
1180 0 : if (HCO_Tasks(n)%ID_I == i-1) then
1181 0 : nlons_task(i) = HCO_Tasks(n)%IM
1182 0 : exit loop
1183 : endif
1184 : enddo loop
1185 : enddo
1186 :
1187 : ! Exclude periodic points for source grids
1188 : ! do n = 0, nPET-1
1189 : ! if (HCO_Tasks(n)%ID_I == nPET_lon-1) then ! Eastern edge
1190 : ! nlons_task(HCO_Tasks(n)%ID_I + 1) = HCO_Tasks(n)%IM - 1
1191 : ! ! overwrites %IM above...
1192 : ! endif
1193 : ! enddo
1194 :
1195 0 : do j = 1, nPET_lat
1196 0 : loop1: do n = 0, nPET-1
1197 0 : if (HCO_Tasks(n)%ID_J == j-1) then
1198 0 : nlats_task(j) = HCO_Tasks(n)%JM
1199 0 : exit loop1
1200 : endif
1201 : enddo loop1
1202 : enddo
1203 :
1204 : !-----------------------------------------------------------------------
1205 : ! Create source grids and allocate coordinates.
1206 : !-----------------------------------------------------------------------
1207 : ! Create pole-based 2D geographic source grid
1208 : HCO_Grid = ESMF_GridCreate1PeriDim( &
1209 0 : countsPerDEDim1=nlons_task, coordDep1=(/1,2/), &
1210 0 : countsPerDEDim2=nlats_task, coordDep2=(/1,2/), &
1211 : indexflag=ESMF_INDEX_GLOBAL, &
1212 : petmap=HCO_petMap, &
1213 : minIndex=(/1,1/), &
1214 0 : rc=RC)
1215 0 : ASSERT_(RC==ESMF_SUCCESS)
1216 :
1217 0 : call ESMF_GridAddCoord(HCO_Grid, staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1218 0 : ASSERT_(RC==ESMF_SUCCESS)
1219 :
1220 0 : call ESMF_GridAddCoord(HCO_Grid, staggerloc=ESMF_STAGGERLOC_CORNER, rc=RC)
1221 0 : ASSERT_(RC==ESMF_SUCCESS) ! why two grids? see remark above (hplin, 2/23/20)
1222 :
1223 : !-----------------------------------------------------------------------
1224 : ! Get pointer and set coordinates - CENTER HEMCO GRID
1225 : !-----------------------------------------------------------------------
1226 : call ESMF_GridGetCoord(HCO_Grid, coordDim=1, localDE=0, &
1227 : computationalLBound=lbnd, &
1228 : computationalUBound=ubnd, &
1229 : farrayPtr=coordX, &
1230 : staggerloc=ESMF_STAGGERLOC_CENTER, &
1231 0 : rc=RC)
1232 : ! Longitude range -180.0, 180.0 is XMid for center staggering
1233 :
1234 : ! Note: Compute bounds are not starting from 1 almost surely
1235 : ! and they should be on the same decomp as the global elems.
1236 : ! So they should be read through the global indices
1237 : ! and not offset ones (a la WRF) (hplin, 2/21/20)
1238 0 : do j = lbnd(2), ubnd(2)
1239 0 : do i = lbnd(1), ubnd(1)
1240 0 : coordX(i, j) = XMid(i, j)
1241 : enddo
1242 : enddo
1243 0 : ASSERT_(RC==ESMF_SUCCESS)
1244 :
1245 : ! Sanity checks here for coordinate staggering
1246 : ! write(iulog,*) ">> [X] IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE, &
1247 : ! "IM, JM, sizeX1,2 = ", &
1248 : ! my_IM, my_JM, size(coordX, 1), size(coordX, 2)
1249 : ! Off by one in periodic dim'n, don't assert I dim here
1250 : ! ASSERT_((ubnd(1) - lbnd(1))==(my_IE-my_IS))
1251 : ! ASSERT_((ubnd(2) - lbnd(2))==(my_JE-my_JS))
1252 :
1253 : call ESMF_GridGetCoord(HCO_Grid, coordDim=2, localDE=0, &
1254 : computationalLBound=lbnd, &
1255 : computationalUBound=ubnd, &
1256 : farrayPtr=coordY, &
1257 : staggerloc=ESMF_STAGGERLOC_CENTER, &
1258 0 : rc=RC)
1259 :
1260 : ! Sanity checks here for coordinate staggering
1261 : ! write(iulog,*) ">> [Y] IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE, &
1262 : ! "IM, JM, sizeY1,2 = ", &
1263 : ! my_IM, my_JM, size(coordY, 1), size(coordY, 2)
1264 : ! ASSERT_((ubnd(1) - lbnd(1))==(my_IE-my_IS))
1265 : ! ASSERT_((ubnd(2) - lbnd(2))==(my_JE-my_JS))
1266 :
1267 0 : do j = lbnd(2), ubnd(2)
1268 0 : do i = lbnd(1), ubnd(1)
1269 0 : coordY(i, j) = YMid(i, j)
1270 : enddo
1271 : enddo
1272 0 : ASSERT_(RC==ESMF_SUCCESS)
1273 :
1274 0 : if(masterproc) then
1275 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateHCO [Ctr] ok"
1276 0 : write(iulog,*) ">> lbnd,ubnd_lon = ", lbnd(1), ubnd(1)
1277 0 : write(iulog,*) ">> lbnd,ubnd_lat = ", lbnd(2), ubnd(2)
1278 0 : write(iulog,*) ">> IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE
1279 0 : write(iulog,*) ">> IM, JM, sizeX1,2 sizeY1,2 = ", &
1280 0 : my_IM, my_JM, size(coordX, 1), size(coordX, 2), size(coordY, 1), size(coordY, 2)
1281 : endif
1282 :
1283 : !-----------------------------------------------------------------------
1284 : ! Get pointer and set coordinates - CORNER HEMCO GRID
1285 : !-----------------------------------------------------------------------
1286 : call ESMF_GridGetCoord(HCO_Grid, coordDim=1, localDE=0, &
1287 : computationalLBound=lbnd, &
1288 : computationalUBound=ubnd, &
1289 : farrayPtr=coordX_E, &
1290 : staggerloc=ESMF_STAGGERLOC_CORNER, &
1291 0 : rc=RC)
1292 :
1293 : ! lbnd(1), ubnd(1) -> 1, 180; 181, 360 ... same as IS, IE ...
1294 :
1295 :
1296 : ! Note: Compute bounds are not starting from 1 almost surely
1297 : ! and they should be on the same decomp as the global elems.
1298 : ! So they should be read through the global indices
1299 : ! and not offset ones (a la WRF) (hplin, 2/21/20)
1300 :
1301 : ! It is a rectilinear grid - so it will not matter what j you pick
1302 : ! for the x-dimension edges, and vice-versa. This is a crude assumption
1303 : ! that is correct for rectilinear but should be revisited. The code
1304 : ! itself is capable of much more. (hplin, 5/29/20)
1305 0 : do j = lbnd(2), ubnd(2)
1306 0 : do i = lbnd(1), ubnd(1)
1307 0 : coordX_E(i, j) = XEdge(i, min(JM, j))
1308 : enddo
1309 : enddo
1310 0 : ASSERT_(RC==ESMF_SUCCESS)
1311 :
1312 : call ESMF_GridGetCoord(HCO_Grid, coordDim=2, localDE=0, &
1313 : computationalLBound=lbnd, &
1314 : computationalUBound=ubnd, &
1315 : farrayPtr=coordY_E, &
1316 : staggerloc=ESMF_STAGGERLOC_CORNER, &
1317 0 : rc=RC)
1318 :
1319 0 : do j = lbnd(2), ubnd(2)
1320 0 : do i = lbnd(1), ubnd(1)
1321 0 : coordY_E(i, j) = YEdge(min(IM, i), j)
1322 : enddo
1323 : enddo
1324 0 : ASSERT_(RC==ESMF_SUCCESS)
1325 :
1326 0 : if(masterproc) then
1327 0 : write(iulog,*) ">> HCO_Grid_ESMF_CreateHCO [Cnr] ok"
1328 0 : write(iulog,*) ">> lbnd,ubnd_lon = ", lbnd(1), ubnd(1)
1329 0 : write(iulog,*) ">> lbnd,ubnd_lat = ", lbnd(2), ubnd(2)
1330 0 : write(iulog,*) ">> IS, IE, JS, JE = ", my_IS, my_IE, my_JS, my_JE
1331 0 : write(iulog,*) ">> IM, JM, sizeX1,2 sizeY1,2 = ", &
1332 0 : my_IM, my_JM, size(coordX, 1), size(coordX, 2), size(coordY, 1), size(coordY, 2)
1333 : endif
1334 :
1335 0 : end subroutine HCO_Grid_ESMF_CreateHCO
1336 : !EOC
1337 : !------------------------------------------------------------------------------
1338 : ! Harmonized Emissions Component (HEMCO) !
1339 : !------------------------------------------------------------------------------
1340 : !BOP
1341 : !
1342 : ! !IROUTINE: HCO_Grid_ESMF_CreateCAMField
1343 : !
1344 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateCAMField creates an ESMF
1345 : ! 2D or 3D field on the mesh representation of CAM physics decomp.
1346 : !\\
1347 : !\\
1348 : ! !INTERFACE:
1349 : !
1350 0 : subroutine HCO_Grid_ESMF_CreateCAMField( field, mesh, name, nlev, RC )
1351 : !
1352 : ! !USES:
1353 : !
1354 : ! MPI Properties from CAM
1355 0 : use cam_logfile, only: iulog
1356 : use spmd_utils, only: masterproc
1357 :
1358 : use ESMF, only: ESMF_TYPEKIND_R8
1359 : use ESMF, only: ESMF_MESHLOC_ELEMENT
1360 : use ESMF, only: ESMF_ArraySpec, ESMF_ArraySpecSet
1361 : use ESMF, only: ESMF_FieldCreate
1362 : !
1363 : ! !INPUT PARAMETERS:
1364 : !
1365 : type(ESMF_Mesh), intent(in) :: mesh
1366 : character(len=*), intent(in) :: name
1367 : integer, intent(in) :: nlev ! nlev==0?2d:3d
1368 : !
1369 : ! !OUTPUT PARAMETERS:
1370 : !
1371 : type(ESMF_Field), intent(out) :: field
1372 : integer, intent(out) :: RC
1373 : !
1374 : ! !REMARKS:
1375 : ! If nlev == 0, field is 2D (i, j), otherwise 3D. 3rd dimension is ungridded.
1376 : !
1377 : ! !REVISION HISTORY:
1378 : ! 21 Feb 2020 - H.P. Lin - Initial version
1379 : !EOP
1380 : !------------------------------------------------------------------------------
1381 : !BOC
1382 : !
1383 : ! !LOCAL VARIABLES:
1384 : !
1385 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateCAMField'
1386 : type(ESMF_ArraySpec) :: arrayspec
1387 :
1388 0 : if(nlev > 0) then
1389 : ! 3D field (i,j,k) with nondistributed vertical
1390 0 : call ESMF_ArraySpecSet(arrayspec, 2, ESMF_TYPEKIND_R8, rc=RC)
1391 0 : ASSERT_(RC==ESMF_SUCCESS)
1392 :
1393 : field = ESMF_FieldCreate(mesh, arrayspec, &
1394 : name=name, &
1395 : ungriddedLBound=(/1/), &
1396 : ungriddedUBound=(/nlev/), &
1397 : gridToFieldMap=(/2/), & ! mapping between grid/field dims ??
1398 0 : meshloc=ESMF_MESHLOC_ELEMENT, rc=RC)
1399 0 : ASSERT_(RC==ESMF_SUCCESS)
1400 : else
1401 : ! 2D field (i,j)
1402 0 : call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=RC)
1403 0 : ASSERT_(RC==ESMF_SUCCESS)
1404 :
1405 : field = ESMF_FieldCreate(mesh, arrayspec, &
1406 : name=name, &
1407 0 : meshloc=ESMF_MESHLOC_ELEMENT, rc=RC)
1408 0 : ASSERT_(RC==ESMF_SUCCESS)
1409 : endif
1410 0 : end subroutine HCO_Grid_ESMF_CreateCAMField
1411 : !EOC
1412 : !------------------------------------------------------------------------------
1413 : ! Harmonized Emissions Component (HEMCO) !
1414 : !------------------------------------------------------------------------------
1415 : !BOP
1416 : !
1417 : ! !IROUTINE: HCO_Grid_ESMF_CreateHCOField
1418 : !
1419 : ! !DESCRIPTION: Subroutine HCO\_Grid\_ESMF\_CreateHCO field creates an ESMF
1420 : ! 2D or 3D field on the HEMCO grid, excluding periodic points.
1421 : !\\
1422 : !\\
1423 : ! !INTERFACE:
1424 : !
1425 0 : subroutine HCO_Grid_ESMF_CreateHCOField( field, grid, name, nlev, RC )
1426 : !
1427 : ! !USES:
1428 : !
1429 : ! MPI Properties from CAM
1430 0 : use cam_logfile, only: iulog
1431 : use spmd_utils, only: masterproc
1432 :
1433 : use ESMF, only: ESMF_TYPEKIND_R8
1434 : use ESMF, only: ESMF_STAGGERLOC_CENTER
1435 : use ESMF, only: ESMF_ArraySpec, ESMF_ArraySpecSet
1436 : use ESMF, only: ESMF_FieldCreate
1437 :
1438 : use ESMF, only: ESMF_GridGet, ESMF_Array, ESMF_ArrayCreate, ESMF_INDEX_GLOBAL
1439 : !
1440 : ! !INPUT PARAMETERS:
1441 : !
1442 : type(ESMF_Grid), intent(in) :: grid
1443 : character(len=*), intent(in) :: name
1444 : integer, intent(in) :: nlev ! nlev==0?2d:3d
1445 : !
1446 : ! !OUTPUT PARAMETERS:
1447 : !
1448 : type(ESMF_Field), intent(out) :: field
1449 : integer, intent(out) :: RC
1450 : !
1451 : ! !REMARKS:
1452 : ! If nlev == 0, field is 2D (i, j), otherwise 3D. 3rd dimension is ungridded.
1453 : ! The grid input parameter accepts both HCO_Grid and HCO2CAM_Grid.
1454 : !
1455 : ! !REVISION HISTORY:
1456 : ! 21 Feb 2020 - H.P. Lin - Initial version
1457 : !EOP
1458 : !------------------------------------------------------------------------------
1459 : !BOC
1460 : !
1461 : ! !LOCAL VARIABLES:
1462 : !
1463 : character(len=*), parameter :: subname = 'HCO_Grid_ESMF_CreateHCOField'
1464 : type(ESMF_ArraySpec) :: arrayspec
1465 :
1466 : type(ESMF_Array) :: array3D, array2D
1467 : type(ESMF_DistGrid) :: distgrid
1468 :
1469 : ! Get grid information from the HEMCO grid:
1470 : call ESMF_GridGet(grid, staggerloc=ESMF_STAGGERLOC_CENTER, &
1471 0 : distgrid=distgrid, rc=RC)
1472 0 : ASSERT_(RC==ESMF_SUCCESS)
1473 :
1474 0 : call ESMF_LogWrite("After ESMF_GridGet", ESMF_LOGMSG_INFO, rc=RC)
1475 :
1476 0 : if(nlev > 0) then
1477 : ! 3D field (i,j,k) with nondistributed vertical
1478 0 : call ESMF_ArraySpecSet(arrayspec, 3, ESMF_TYPEKIND_R8, rc=RC)
1479 0 : ASSERT_(RC==ESMF_SUCCESS)
1480 :
1481 0 : call ESMF_LogWrite("After array3D-ESMF_ArraySpecSet", ESMF_LOGMSG_INFO, rc=RC)
1482 :
1483 : array3D = ESMF_ArrayCreate(arrayspec=arrayspec, &
1484 : distgrid=distgrid, &
1485 : computationalEdgeUWidth=(/1,0/), &
1486 : undistLBound=(/1/), undistUBound=(/nlev/), &
1487 0 : indexflag=ESMF_INDEX_GLOBAL, rc=RC)
1488 0 : ASSERT_(RC==ESMF_SUCCESS)
1489 :
1490 0 : call ESMF_LogWrite("After array3D-ESMF_ArrayCreate", ESMF_LOGMSG_INFO, rc=RC)
1491 :
1492 : field = ESMF_FieldCreate(grid, array3D, & ! grid, arrayspec
1493 : name=name, &
1494 : ungriddedLBound=(/1/), &
1495 : ungriddedUBound=(/nlev/), &
1496 0 : staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1497 0 : ASSERT_(RC==ESMF_SUCCESS)
1498 :
1499 0 : call ESMF_LogWrite("After array3D-ESMF_FieldCreate", ESMF_LOGMSG_INFO, rc=RC)
1500 : else
1501 : ! 2D field (i,j)
1502 0 : call ESMF_ArraySpecSet(arrayspec, 2, ESMF_TYPEKIND_R8, rc=RC)
1503 0 : ASSERT_(RC==ESMF_SUCCESS)
1504 :
1505 : field = ESMF_FieldCreate(grid, arrayspec, &
1506 : name=name, &
1507 0 : staggerloc=ESMF_STAGGERLOC_CENTER, rc=RC)
1508 0 : ASSERT_(RC==ESMF_SUCCESS)
1509 : endif
1510 0 : end subroutine HCO_Grid_ESMF_CreateHCOField
1511 : !EOC
1512 : !------------------------------------------------------------------------------
1513 : ! Harmonized Emissions Component (HEMCO) !
1514 : !------------------------------------------------------------------------------
1515 : !BOP
1516 : !
1517 : ! !IROUTINE: HCO_Grid_HCO2CAM_3D
1518 : !
1519 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_3D regrids a HEMCO lat-lon grid
1520 : ! field (i,j,l) to CAM physics mesh (k,i).
1521 : !\\
1522 : !\\
1523 : ! !INTERFACE:
1524 : !
1525 0 : subroutine HCO_Grid_HCO2CAM_3D(hcoArray, camArray)
1526 : !
1527 : ! !USES:
1528 : !
1529 0 : use ESMF, only: ESMF_FieldRegrid
1530 : use ESMF, only: ESMF_TERMORDER_SRCSEQ, ESMF_TERMORDER_FREE
1531 :
1532 : ! MPI Properties from CAM
1533 : use cam_logfile, only: iulog
1534 : use spmd_utils, only: masterproc
1535 : !
1536 : ! !INPUT PARAMETERS:
1537 : !
1538 : real(r8), intent(in) :: hcoArray(my_IS:my_IE,my_JS:my_JE,1:LM)
1539 : !
1540 : ! !OUTPUT PARAMETERS:
1541 : !
1542 : real(r8), intent(inout) :: camArray(1:LM,1:my_CE) ! Col and lev start on 1
1543 : !
1544 : ! !REMARKS:
1545 : ! (1) There is no vertical regridding. Also, chunk and lev indices are assumed to
1546 : ! start at 1 always.
1547 : ! (2) The subroutine's inner workings abstract ESMF from the user, but this routine
1548 : ! needs to be called from within the gridded component (as suggested by Steve)
1549 : ! (3) Note that this automatically flips the vertical through HCO_ESMF_Get3DField.
1550 : ! In CAM, layer 1 is TOA. Layer 1 is ground in HEMCO and thus we flip the vertical
1551 : ! in regrid routines when retrieving the data.
1552 : !
1553 : ! !REVISION HISTORY:
1554 : ! 24 Feb 2020 - H.P. Lin - Initial version
1555 : ! 30 May 2020 - H.P. Lin - Remember to flip the vertical!
1556 : ! 06 Dec 2020 - H.P. Lin - Add a kludge to skip regridding (removed 12 Jan 2023)
1557 : !EOP
1558 : !------------------------------------------------------------------------------
1559 : !BOC
1560 : !
1561 : ! !LOCAL VARIABLES:
1562 : !
1563 : character(len=*), parameter :: subname = 'HCO_Grid_HCO2CAM_3D'
1564 : integer :: RC
1565 :
1566 0 : call HCO_ESMF_Set3DHCO(HCO_3DFld, hcoArray, my_IS, my_IE, my_JS, my_JE, 1, LM)
1567 :
1568 :
1569 : call ESMF_FieldRegrid(HCO_3DFld, CAM_3DFld, HCO2CAM_RouteHandle_3D, &
1570 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1571 0 : checkflag=.true., rc=RC)
1572 0 : ASSERT_(RC==ESMF_SUCCESS)
1573 :
1574 : ! (field_in, data_out, IS, IE, JS, JE)
1575 : ! For chunks, "I" is lev, "J" is chunk index, confusing, you are warned
1576 : ! (Physics "2D" fields on mesh are actually "3D" data)
1577 0 : call HCO_ESMF_Get2DField(CAM_3DFld, camArray, 1, LM, 1, my_CE, flip=.true.)
1578 :
1579 0 : end subroutine HCO_Grid_HCO2CAM_3D
1580 : !EOC
1581 : !------------------------------------------------------------------------------
1582 : ! Harmonized Emissions Component (HEMCO) !
1583 : !------------------------------------------------------------------------------
1584 : !BOP
1585 : !
1586 : ! !IROUTINE: HCO_Grid_CAM2HCO_3D
1587 : !
1588 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_3D regrids a HEMCO lat-lon grid
1589 : ! field (i,j,l) to CAM physics mesh (k,i).
1590 : !\\
1591 : !\\
1592 : ! !INTERFACE:
1593 : !
1594 0 : subroutine HCO_Grid_CAM2HCO_3D(camArray, hcoArray)
1595 : !
1596 : ! !USES:
1597 : !
1598 0 : use ESMF, only: ESMF_FieldRegrid, ESMF_TERMORDER_SRCSEQ
1599 :
1600 : ! MPI Properties from CAM
1601 : use cam_logfile, only: iulog
1602 : use spmd_utils, only: masterproc
1603 : !
1604 : ! !INPUT PARAMETERS:
1605 : !
1606 : real(r8), intent(in) :: camArray(1:LM,1:my_CE) ! Col and lev start on 1
1607 : !
1608 : ! !OUTPUT PARAMETERS:
1609 : !
1610 : real(r8), intent(inout) :: hcoArray(my_IS:my_IE,my_JS:my_JE,1:LM)
1611 : !
1612 : ! !REMARKS:
1613 : ! (1) There is no vertical regridding. Also, chunk and lev indices are assumed to
1614 : ! start at 1 always.
1615 : ! (2) The subroutine's inner workings abstract ESMF from the user, but this routine
1616 : ! needs to be called from within the gridded component (as suggested by Steve)
1617 : ! (3) Note that this automatically flips the vertical through HCO_ESMF_Get3DField.
1618 : ! In CAM, layer 1 is TOA. Layer 1 is ground in HEMCO and thus we flip the vertical
1619 : ! in regrid routines when retrieving the data.
1620 : !
1621 : ! !REVISION HISTORY:
1622 : ! 24 Feb 2020 - H.P. Lin - Initial version
1623 : ! 30 May 2020 - H.P. Lin - Remember to flip the vertical!
1624 : !EOP
1625 : !------------------------------------------------------------------------------
1626 : !BOC
1627 : !
1628 : ! !LOCAL VARIABLES:
1629 : !
1630 : character(len=*), parameter :: subname = 'HCO_Grid_CAM2HCO_3D'
1631 : integer :: RC
1632 :
1633 : integer :: J, K
1634 :
1635 : ! (field, data, KS, KE, CS, CE)
1636 0 : call HCO_ESMF_Set3DCAM(CAM_3DFld, camArray, 1, LM, 1, my_CE)
1637 :
1638 : call ESMF_FieldRegrid(CAM_3DFld, HCO_3DFld, CAM2HCO_RouteHandle_3D, &
1639 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1640 0 : checkflag=.true., rc=RC)
1641 0 : ASSERT_(RC==ESMF_SUCCESS)
1642 :
1643 0 : call HCO_ESMF_Get3DField(HCO_3DFld, hcoArray, my_IS, my_IE, my_JS, my_JE, 1, LM, flip=.true.)
1644 :
1645 0 : end subroutine HCO_Grid_CAM2HCO_3D
1646 : !EOC
1647 : !------------------------------------------------------------------------------
1648 : ! Harmonized Emissions Component (HEMCO) !
1649 : !------------------------------------------------------------------------------
1650 : !BOP
1651 : !
1652 : ! !IROUTINE: HCO_Grid_HCO2CAM_2D
1653 : !
1654 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_2D regrids a HEMCO lat-lon grid
1655 : ! field (i,j) to CAM physics mesh (i).
1656 : !\\
1657 : !\\
1658 : ! !INTERFACE:
1659 : !
1660 0 : subroutine HCO_Grid_HCO2CAM_2D(hcoArray, camArray)
1661 : !
1662 : ! !USES:
1663 : !
1664 0 : use ESMF, only: ESMF_FieldRegrid
1665 : use ESMF, only: ESMF_TERMORDER_SRCSEQ
1666 :
1667 : ! MPI Properties from CAM
1668 : use cam_logfile, only: iulog
1669 : use spmd_utils, only: masterproc
1670 : !
1671 : ! !INPUT PARAMETERS:
1672 : !
1673 : real(r8), intent(in) :: hcoArray(my_IS:my_IE,my_JS:my_JE)
1674 : !
1675 : ! !OUTPUT PARAMETERS:
1676 : !
1677 : real(r8), intent(inout) :: camArray(1:my_CE) ! Col and lev start on 1
1678 : !
1679 : ! !REMARKS:
1680 : ! (1) Used for 2-D fields, which are 1-D in CAM.
1681 : !
1682 : ! !REVISION HISTORY:
1683 : ! 31 Mar 2020 - H.P. Lin - Initial version
1684 : !EOP
1685 : !------------------------------------------------------------------------------
1686 : !BOC
1687 : !
1688 : ! !LOCAL VARIABLES:
1689 : !
1690 : character(len=*), parameter :: subname = 'HCO_Grid_HCO2CAM_2D'
1691 : integer :: RC
1692 :
1693 0 : call HCO_ESMF_Set2DHCO(HCO_2DFld, hcoArray, my_IS, my_IE, my_JS, my_JE)
1694 :
1695 : call ESMF_FieldRegrid(HCO_2DFld, CAM_2DFld, HCO2CAM_RouteHandle_2D, &
1696 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1697 0 : checkflag=.true., rc=RC)
1698 0 : ASSERT_(RC==ESMF_SUCCESS)
1699 :
1700 : ! HCO_ESMF_Get1DField(field_in, data_out, IS, IE)
1701 : ! (Physics "1D" fields on mesh are actually "2D" data)
1702 0 : call HCO_ESMF_Get1DField(CAM_2DFld, camArray, 1, my_CE)
1703 :
1704 0 : end subroutine HCO_Grid_HCO2CAM_2D
1705 : !EOC
1706 : !------------------------------------------------------------------------------
1707 : ! Harmonized Emissions Component (HEMCO) !
1708 : !------------------------------------------------------------------------------
1709 : !BOP
1710 : !
1711 : ! !IROUTINE: HCO_Grid_CAM2HCO_2D
1712 : !
1713 : ! !DESCRIPTION: Subroutine HCO\_Grid\_HCO2CAM\_2D regrids a HEMCO lat-lon grid
1714 : ! field (i,j,l) to CAM physics mesh (k,i).
1715 : !\\
1716 : !\\
1717 : ! !INTERFACE:
1718 : !
1719 0 : subroutine HCO_Grid_CAM2HCO_2D(camArray, hcoArray)
1720 : !
1721 : ! !USES:
1722 : !
1723 0 : use ESMF, only: ESMF_FieldRegrid, ESMF_TERMORDER_SRCSEQ
1724 :
1725 : ! MPI Properties from CAM
1726 : use cam_logfile, only: iulog
1727 : use spmd_utils, only: masterproc
1728 : !
1729 : ! !INPUT PARAMETERS:
1730 : !
1731 : real(r8), intent(in) :: camArray(1:my_CE) ! Col and lev start on 1
1732 : !
1733 : ! !OUTPUT PARAMETERS:
1734 : !
1735 : real(r8), intent(inout) :: hcoArray(my_IS:my_IE,my_JS:my_JE)
1736 : !
1737 : ! !REVISION HISTORY:
1738 : ! 31 Mar 2020 - H.P. Lin - Initial version
1739 : !EOP
1740 : !------------------------------------------------------------------------------
1741 : !BOC
1742 : !
1743 : ! !LOCAL VARIABLES:
1744 : !
1745 : character(len=*), parameter :: subname = 'HCO_Grid_CAM2HCO_2D'
1746 : integer :: RC
1747 :
1748 : integer :: J
1749 :
1750 : ! HCO_ESMF_Set2DCAM(field, data, CS, CE)
1751 0 : call HCO_ESMF_Set2DCAM(CAM_2DFld, camArray, 1, my_CE)
1752 :
1753 : call ESMF_FieldRegrid(CAM_2DFld, HCO_2DFld, CAM2HCO_RouteHandle_2D, &
1754 : termorderflag=ESMF_TERMORDER_SRCSEQ, &
1755 0 : checkflag=.true., rc=RC)
1756 0 : ASSERT_(RC==ESMF_SUCCESS)
1757 :
1758 0 : call HCO_ESMF_Get2DField(HCO_2DFld, hcoArray, my_IS, my_IE, my_JS, my_JE)
1759 :
1760 : ! Kludge for periodic point
1761 : ! It seems like the last point for the PET in the x-edge direction is messed up,
1762 : ! because it is the "periodic" point wrapping around the globe. This part is
1763 : ! not smooth and needs to be manually extrapolated.
1764 : !
1765 : ! This is a kludge as we want to revisit the regrid mechanism later.
1766 : ! For now fix it by copying the edge, not IDAVG
1767 0 : if(my_IE .eq. IM) then
1768 0 : do J = my_JS, my_JE
1769 0 : if(hcoArray(my_IE, J) .le. 0.000001_r8) then
1770 0 : hcoArray(my_IE, J) = hcoArray(my_IE-1, J)
1771 : endif
1772 : enddo
1773 : endif
1774 :
1775 0 : end subroutine HCO_Grid_CAM2HCO_2D
1776 : !EOC
1777 : !------------------------------------------------------------------------------
1778 : ! Harmonized Emissions Component (HEMCO) !
1779 : !------------------------------------------------------------------------------
1780 : !BOP
1781 : !
1782 : ! !IROUTINE: HCO_ESMF_Set2DHCO
1783 : !
1784 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set2DHCO sets values of a ESMF field on
1785 : ! the HEMCO lat-lon grid. (Internal use)
1786 : !\\
1787 : !\\
1788 : ! !INTERFACE:
1789 : !
1790 0 : subroutine HCO_ESMF_Set2DHCO(field, data, IS, IE, JS, JE)
1791 : !
1792 : ! !USES:
1793 : !
1794 0 : use ESMF, only: ESMF_FieldGet
1795 : !
1796 : ! !INPUT PARAMETERS:
1797 : !
1798 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1799 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
1800 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
1801 : real(r8), intent(in) :: data(IS:IE, JS:JE) ! Data to write in
1802 : ! !REMARKS:
1803 : ! Note there might need to be handling for periodic points?
1804 : !
1805 : ! !REVISION HISTORY:
1806 : ! 24 Feb 2020 - H.P. Lin - Initial version
1807 : !EOP
1808 : !------------------------------------------------------------------------------
1809 : !BOC
1810 : !
1811 : ! !LOCAL VARIABLES:
1812 : !
1813 : character(len=*), parameter :: subname = 'HCO_ESMF_Set2DHCO'
1814 : integer :: I, J
1815 : integer :: RC
1816 : integer :: lbnd(2), ubnd(2)
1817 0 : real(r8), pointer :: fptr(:,:)
1818 :
1819 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1820 : computationalLBound=lbnd, &
1821 0 : computationalUBound=ubnd, rc=RC)
1822 0 : ASSERT_(RC==ESMF_SUCCESS)
1823 0 : fptr(:,:) = 0.0_r8 ! Arbitrary missval
1824 0 : do J = lbnd(2), ubnd(2)
1825 0 : do I = lbnd(1), ubnd(1)
1826 0 : fptr(I,J) = data(I,J)
1827 : enddo
1828 : enddo
1829 :
1830 0 : end subroutine HCO_ESMF_Set2DHCO
1831 : !EOC
1832 : !------------------------------------------------------------------------------
1833 : ! Harmonized Emissions Component (HEMCO) !
1834 : !------------------------------------------------------------------------------
1835 : !BOP
1836 : !
1837 : ! !IROUTINE: HCO_ESMF_Set3DHCO
1838 : !
1839 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set3DHCO sets values of a ESMF field on
1840 : ! the HEMCO lat-lon grid. (Internal use)
1841 : !\\
1842 : !\\
1843 : ! !INTERFACE:
1844 : !
1845 0 : subroutine HCO_ESMF_Set3DHCO(field, data, IS, IE, JS, JE, KS, KE)
1846 : !
1847 : ! !USES:
1848 : !
1849 0 : use ESMF, only: ESMF_FieldGet
1850 : !
1851 : ! !INPUT PARAMETERS:
1852 : !
1853 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1854 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
1855 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
1856 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
1857 : real(r8), intent(in) :: data(IS:IE, JS:JE, KS:KE) ! Data to write in
1858 : ! !REMARKS:
1859 : ! Note there might need to be handling for periodic points?
1860 : !
1861 : ! !REVISION HISTORY:
1862 : ! 24 Feb 2020 - H.P. Lin - Initial version
1863 : !EOP
1864 : !------------------------------------------------------------------------------
1865 : !BOC
1866 : !
1867 : ! !LOCAL VARIABLES:
1868 : !
1869 : character(len=*), parameter :: subname = 'HCO_ESMF_Set3DHCO'
1870 : integer :: I, J, K
1871 : integer :: RC
1872 : integer :: lbnd(3), ubnd(3)
1873 0 : real(r8), pointer :: fptr(:,:,:)
1874 :
1875 0 : call ESMF_LogWrite("In HCO_ESMF_Set3DHCO before ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1876 :
1877 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1878 : computationalLBound=lbnd, &
1879 0 : computationalUBound=ubnd, rc=RC)
1880 0 : ASSERT_(RC==ESMF_SUCCESS)
1881 :
1882 0 : call ESMF_LogWrite("In HCO_ESMF_Set3DHCO after ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1883 0 : fptr(:,:,:) = 0.0_r8 ! Arbitrary missval
1884 :
1885 : ! 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)
1886 0 : do K = lbnd(3), ubnd(3)
1887 0 : do J = lbnd(2), ubnd(2)
1888 0 : do I = lbnd(1), ubnd(1)
1889 0 : fptr(I,J,K) = data(I,J,K)
1890 : enddo
1891 : enddo
1892 : enddo
1893 :
1894 0 : end subroutine HCO_ESMF_Set3DHCO
1895 : !EOC
1896 : !------------------------------------------------------------------------------
1897 : ! Harmonized Emissions Component (HEMCO) !
1898 : !------------------------------------------------------------------------------
1899 : !BOP
1900 : !
1901 : ! !IROUTINE: HCO_ESMF_Set2DCAM
1902 : !
1903 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set2DCAM sets values of a ESMF field on
1904 : ! the physics mesh. (Internal use)
1905 : !\\
1906 : !\\
1907 : ! !INTERFACE:
1908 : !
1909 0 : subroutine HCO_ESMF_Set2DCAM(field, data, CS, CE)
1910 : !
1911 : ! !USES:
1912 : !
1913 0 : use ESMF, only: ESMF_FieldGet
1914 : !
1915 : ! !INPUT PARAMETERS:
1916 : !
1917 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1918 : integer, intent(in) :: CS, CE ! Start and end indices (chunk)
1919 : real(r8), intent(in) :: data(CS:CE) ! Data to write in
1920 : ! !REMARKS:
1921 : ! Remember that "2-D" fields in physics chunk is actually 3-D (k,i), and 2-D
1922 : ! is actually just column-level data. The target is a mesh so in ESMF it is
1923 : ! the same.
1924 : !
1925 : ! !REVISION HISTORY:
1926 : ! 23 Feb 2020 - H.P. Lin - Initial version
1927 : !EOP
1928 : !------------------------------------------------------------------------------
1929 : !BOC
1930 : !
1931 : ! !LOCAL VARIABLES:
1932 : !
1933 : character(len=*), parameter :: subname = 'HCO_ESMF_Set2DCAM'
1934 : integer :: I
1935 : integer :: RC
1936 : integer :: lbnd(1), ubnd(1)
1937 0 : real(r8), pointer :: fptr(:)
1938 :
1939 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1940 : computationalLBound=lbnd, &
1941 0 : computationalUBound=ubnd, rc=RC)
1942 0 : ASSERT_(RC==ESMF_SUCCESS)
1943 0 : fptr(:) = 0.0_r8 ! Arbitrary missval
1944 0 : do I = lbnd(1), ubnd(1)
1945 0 : fptr(i) = data(i)
1946 : enddo
1947 :
1948 0 : end subroutine HCO_ESMF_Set2DCAM
1949 : !EOC
1950 : !------------------------------------------------------------------------------
1951 : ! Harmonized Emissions Component (HEMCO) !
1952 : !------------------------------------------------------------------------------
1953 : !BOP
1954 : !
1955 : ! !IROUTINE: HCO_ESMF_Set3DCAM
1956 : !
1957 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Set3DCAM sets values of a ESMF field on
1958 : ! the physics mesh. (Internal use)
1959 : !\\
1960 : !\\
1961 : ! !INTERFACE:
1962 : !
1963 0 : subroutine HCO_ESMF_Set3DCAM(field, data, KS, KE, CS, CE)
1964 : !
1965 : ! !USES:
1966 : !
1967 0 : use ESMF, only: ESMF_FieldGet
1968 : !
1969 : ! !INPUT PARAMETERS:
1970 : !
1971 : type(ESMF_Field), intent(in) :: field ! intent(in) because write to ptr
1972 : integer, intent(in) :: CS, CE ! Start and end indices (chunk)
1973 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
1974 : real(r8), intent(in) :: data(KS:KE, CS:CE) ! Data to write in
1975 : ! !REMARKS:
1976 : ! This is the actual 3-D data routine, pointers (i,k) on a mesh
1977 : !
1978 : ! !REVISION HISTORY:
1979 : ! 23 Feb 2020 - H.P. Lin - Initial version
1980 : !EOP
1981 : !------------------------------------------------------------------------------
1982 : !BOC
1983 : !
1984 : ! !LOCAL VARIABLES:
1985 : !
1986 : character(len=*), parameter :: subname = 'HCO_ESMF_Set3DCAM'
1987 : integer :: I, K
1988 : integer :: RC
1989 : integer :: lbnd(2), ubnd(2)
1990 0 : real(r8), pointer :: fptr(:,:)
1991 :
1992 : ! call ESMF_LogWrite("In HCO_ESMF_Set3DCAM before ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
1993 :
1994 : call ESMF_FieldGet(field, localDE=0, farrayPtr=fptr, &
1995 : computationalLBound=lbnd, &
1996 0 : computationalUBound=ubnd, rc=RC)
1997 0 : ASSERT_(RC==ESMF_SUCCESS)
1998 :
1999 :
2000 : ! call ESMF_LogWrite("In HCO_ESMF_Set3DCAM after ESMF_FieldGet", ESMF_LOGMSG_INFO, rc=RC)
2001 0 : fptr(:,:) = 0.0_r8 ! Arbitrary missval
2002 0 : do I = lbnd(2), ubnd(2)
2003 0 : do K = lbnd(1), ubnd(1)
2004 0 : fptr(k,i) = data(k,i)
2005 : enddo
2006 : enddo
2007 :
2008 0 : end subroutine HCO_ESMF_Set3DCAM
2009 : !EOC
2010 : !------------------------------------------------------------------------------
2011 : ! Harmonized Emissions Component (HEMCO) !
2012 : !------------------------------------------------------------------------------
2013 : !BOP
2014 : !
2015 : ! !IROUTINE: HCO_ESMF_Get1DField
2016 : !
2017 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get1DField gets a pointer to an 1-D ESMF
2018 : ! field.
2019 : !\\
2020 : !\\
2021 : ! !INTERFACE:
2022 : !
2023 0 : subroutine HCO_ESMF_Get1DField(field_in, data_out, IS, IE)
2024 : !
2025 : ! !USES:
2026 : !
2027 0 : use ESMF, only: ESMF_FieldGet
2028 : !
2029 : ! !INPUT PARAMETERS:
2030 : !
2031 : type(ESMF_Field), intent(in) :: field_in
2032 : integer, intent(in) :: IS, IE ! Start and end indices
2033 : !
2034 : ! !OUTPUT PARAMETERS:
2035 : !
2036 : real(r8), intent(out) :: data_out(IS:IE)
2037 : !
2038 : ! !REVISION HISTORY:
2039 : ! 23 Feb 2020 - H.P. Lin - Initial version
2040 : !EOP
2041 : !------------------------------------------------------------------------------
2042 : !BOC
2043 : !
2044 : ! !LOCAL VARIABLES:
2045 : !
2046 : character(len=*), parameter :: subname = 'HCO_ESMF_Get1DField'
2047 0 : real(r8), pointer :: fptr(:)
2048 : integer :: RC
2049 :
2050 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2051 0 : ASSERT_(RC==ESMF_SUCCESS)
2052 0 : data_out(:) = fptr(:)
2053 :
2054 0 : end subroutine HCO_ESMF_Get1DField
2055 : !EOC
2056 : !------------------------------------------------------------------------------
2057 : ! Harmonized Emissions Component (HEMCO) !
2058 : !------------------------------------------------------------------------------
2059 : !BOP
2060 : !
2061 : ! !IROUTINE: HCO_ESMF_Get2DField
2062 : !
2063 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get2DField gets a pointer to an 2-D ESMF
2064 : ! field.
2065 : !\\
2066 : !\\
2067 : ! !INTERFACE:
2068 : !
2069 0 : subroutine HCO_ESMF_Get2DField(field_in, data_out, IS, IE, JS, JE, flip)
2070 : !
2071 : ! !USES:
2072 : !
2073 0 : use ESMF, only: ESMF_FieldGet
2074 : !
2075 : ! !INPUT PARAMETERS:
2076 : !
2077 : type(ESMF_Field), intent(in) :: field_in
2078 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
2079 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
2080 : logical, optional, intent(in) :: flip ! Flip the FIRST dimension? Presence will flip
2081 : !
2082 : ! !OUTPUT PARAMETERS:
2083 : !
2084 : real(r8), intent(out) :: data_out(IS:IE, JS:JE)
2085 : !
2086 : ! !REMARKS:
2087 : ! If the flip argument is PROVIDED, then the FIRST dimension (IS:IE) will be flipped.
2088 : ! This is used when retrieving a CAM field regridded from HEMCO, when you have to flip
2089 : ! the vertical to fit correctly. The dimensions of a CAM 3d data is a 2d chunked array
2090 : ! with indices (k, i) so... the first dimension has to be flipped instead.
2091 : !
2092 : ! YES - THE FLIP ARGUMENT IS NOT READ. ITS MERE PRESENCE WILL TRIGGER .TRUE. FOR FLIPPING.
2093 : ! THIS IS INTENDED BEHAVIOR TO CIRCUMVENT AN INTEL COMPILER BUG.
2094 : !
2095 : ! !REVISION HISTORY:
2096 : ! 23 Feb 2020 - H.P. Lin - Initial version
2097 : ! 30 May 2020 - H.P. Lin - Add flip parameter
2098 : ! 03 Jan 2021 - H.P. Lin - I hate compiler bugs with PRESENT var
2099 : !EOP
2100 : !------------------------------------------------------------------------------
2101 : !BOC
2102 : !
2103 : ! !LOCAL VARIABLES:
2104 : !
2105 : character(len=*), parameter :: subname = 'HCO_ESMF_Get2DField'
2106 0 : real(r8), pointer :: fptr(:,:)
2107 : integer :: RC
2108 : logical :: flip1d = .false.
2109 :
2110 0 : if(present(flip)) then
2111 0 : flip1d = .true.
2112 : else
2113 0 : flip1d = .false.
2114 : endif
2115 :
2116 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2117 0 : ASSERT_(RC==ESMF_SUCCESS)
2118 :
2119 0 : if(flip1d) then
2120 0 : data_out(IS:IE:1, :) = fptr(IE:IS:-1, :)
2121 : else
2122 0 : data_out(:,:) = fptr(:,:)
2123 : endif
2124 :
2125 0 : end subroutine HCO_ESMF_Get2DField
2126 : !EOC
2127 : !------------------------------------------------------------------------------
2128 : ! Harmonized Emissions Component (HEMCO) !
2129 : !------------------------------------------------------------------------------
2130 : !BOP
2131 : !
2132 : ! !IROUTINE: HCO_ESMF_Get3DField
2133 : !
2134 : ! !DESCRIPTION: Subroutine HCO\_ESMF\_Get3DField gets a pointer to an 3-D ESMF
2135 : ! field.
2136 : !\\
2137 : !\\
2138 : ! !INTERFACE:
2139 : !
2140 0 : subroutine HCO_ESMF_Get3DField(field_in, data_out, IS, IE, JS, JE, KS, KE, flip)
2141 : !
2142 : ! !USES:
2143 : !
2144 0 : use ESMF, only: ESMF_FieldGet
2145 : !
2146 : ! !INPUT PARAMETERS:
2147 : !
2148 : type(ESMF_Field), intent(in) :: field_in
2149 : integer, intent(in) :: IS, IE ! Start and end indices (dim1)
2150 : integer, intent(in) :: JS, JE ! Start and end indices (dim2)
2151 : integer, intent(in) :: KS, KE ! Start and end indices (dim3)
2152 : logical, optional, intent(in) :: flip ! Flip the 3rd dimension?
2153 : !
2154 : ! !OUTPUT PARAMETERS:
2155 : !
2156 : real(r8), intent(out) :: data_out(IS:IE, JS:JE, KS:KE)
2157 : !
2158 : ! !REMARKS:
2159 : ! If the flip argument is set to .true., then the third dimension (KS:KE) will be flipped.
2160 : ! This is used when retrieving a HEMCO field regridded from CAM, when you have to flip
2161 : ! the vertical to fit correctly.
2162 : !
2163 : ! YES - THE FLIP ARGUMENT IS NOT READ. ITS MERE PRESENCE WILL TRIGGER .TRUE. FOR FLIPPING.
2164 : ! THIS IS INTENDED BEHAVIOR TO CIRCUMVENT AN INTEL COMPILER BUG.
2165 : !
2166 : ! !REVISION HISTORY:
2167 : ! 23 Feb 2020 - H.P. Lin - Initial version
2168 : ! 30 May 2020 - H.P. Lin - Add flip parameter
2169 : !EOP
2170 : !------------------------------------------------------------------------------
2171 : !BOC
2172 : !
2173 : ! !LOCAL VARIABLES:
2174 : !
2175 : character(len=*), parameter :: subname = 'HCO_ESMF_Get3DField'
2176 0 : real(r8), pointer :: fptr(:,:,:)
2177 : integer :: RC
2178 : logical :: flip3d = .false.
2179 :
2180 0 : if(present(flip)) then
2181 0 : flip3d = .true.
2182 : else
2183 0 : flip3d = .false.
2184 : endif
2185 :
2186 0 : call ESMF_FieldGet(field_in, localDE=0, farrayPtr=fptr, rc=RC)
2187 0 : ASSERT_(RC==ESMF_SUCCESS)
2188 :
2189 0 : if(flip3d) then
2190 0 : data_out(:,:,KS:KE:1) = fptr(:,:,KE:KS:-1)
2191 : else
2192 0 : data_out(:,:,:) = fptr(:,:,:)
2193 : endif
2194 :
2195 0 : end subroutine HCO_ESMF_Get3DField
2196 : !EOC
2197 0 : end module hco_esmf_grid
|
