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