diff --git a/bld/build-namelist b/bld/build-namelist
index e670c6ca3c..6ac0d32c9d 100755
--- a/bld/build-namelist
+++ b/bld/build-namelist
@@ -705,6 +705,7 @@ if (!($simple_phys or $aqua_mode)) {
my $chem_nitrodep = chem_has_species($cfg, 'NO') and chem_has_species($cfg, 'NH3');
if ((!$chem_nitrodep) or ($chem =~ /geoschem/)) {
add_default($nl, 'stream_ndep_mesh_filename');
+ add_default($nl, 'stream_ndep_mapalgo');
add_default($nl, 'stream_ndep_data_filename', 'sim_year'=>$sim_year);
add_default($nl, 'stream_ndep_year_first', 'val'=>$sim_year_first);
add_default($nl, 'stream_ndep_year_last', 'val'=>$sim_year_last);
diff --git a/bld/namelist_files/namelist_defaults_cam.xml b/bld/namelist_files/namelist_defaults_cam.xml
index 38975992be..166d3e31d0 100644
--- a/bld/namelist_files/namelist_defaults_cam.xml
+++ b/bld/namelist_files/namelist_defaults_cam.xml
@@ -2117,6 +2117,7 @@
share/meshes/fv0.9x1.25_141008_polemod_ESMFmesh.nc
+bilinear
lnd/clm2/ndepdata/fndep_clm_hist_b.e21.BWHIST.f09_g17.CMIP6-historical-WACCM.ensmean_1849-2015_monthly_0.9x1.25_c180926.nc
atm/cam/chem/trop_strat_mam5_ts4_aero_spinup/ndep/fndep_clm_B1850C_MT4s.1850monthly.ne30_c251220.nc
diff --git a/bld/namelist_files/namelist_definition.xml b/bld/namelist_files/namelist_definition.xml
index 5ff967b2c4..861992aad0 100644
--- a/bld/namelist_files/namelist_definition.xml
+++ b/bld/namelist_files/namelist_definition.xml
@@ -8248,6 +8248,15 @@ Nitrogen deposition stream data filename.
Grid mesh file corresponding to stream_ndep_data_filename.
+
+Spatial mapping (regrid) method used to map the nitrogen deposition stream data from
+the stream grid to the model grid. Default 'bilinear'. Use 'redist' (pure index-based
+redistribution -- no interpolation weights or ESMF route handle) when the ndep data file
+is already on the model grid; this greatly reduces memory and setup cost at high
+resolution (e.g. ne1024pg2).
+
+
File containing MEGAN emissions factors.
diff --git a/src/chemistry/aerosol/soil_erod_mod.F90 b/src/chemistry/aerosol/soil_erod_mod.F90
index 33d6108761..483dce4114 100644
--- a/src/chemistry/aerosol/soil_erod_mod.F90
+++ b/src/chemistry/aerosol/soil_erod_mod.F90
@@ -24,26 +24,36 @@ subroutine soil_erod_init( dust_emis_fact, soil_erod_file )
use interpolate_data, only: lininterp_init, lininterp, lininterp_finish, interp_type
use ppgrid, only: begchunk, endchunk, pcols
use mo_constants, only: pi, d2r
- use pio, only: file_desc_t,pio_inq_dimid,pio_inq_dimlen,pio_get_var,pio_inq_varid, PIO_NOWRITE
use phys_grid, only: get_ncols_p, get_rlat_all_p, get_rlon_all_p
- use cam_pio_utils, only: cam_pio_openfile
use ioFileMod, only: getfil
+ use netcdf, only: nf90_open, nf90_nowrite, nf90_close, nf90_noerr, &
+ nf90_inq_dimid, nf90_inquire_dimension, &
+ nf90_inq_varid, nf90_get_var
+ use cam_shmem_mod, only: cam_shmem_alloc_r8_2d, cam_shmem_fence, &
+ cam_shmem_free, cam_shmem_is_leader, &
+ cam_shmem_leader_comm, cam_shmem_npes_per_node
+#ifdef SPMD
+ use mpishorthand, only: mpicom, mpiint, mpir8
+#endif
real(r8), intent(in) :: dust_emis_fact
character(len=*), intent(in) :: soil_erod_file
- real(r8), allocatable :: soil_erodibility_in(:,:) ! temporary input array
+ real(r8), pointer :: soil_erodibility_in(:,:) => null() ! node-shared input field
real(r8), allocatable :: dst_lons(:)
real(r8), allocatable :: dst_lats(:)
character(len=cl) :: infile
integer :: did, vid, nlat, nlon
- type(file_desc_t) :: ncid
+ integer :: ncid, iret
+ integer :: win_serod ! per-node shared-memory window for soil_erodibility_in
type(interp_type) :: lon_wgts, lat_wgts
real(r8) :: to_lats(pcols), to_lons(pcols)
integer :: c, ncols, ierr
real(r8), parameter :: zero=0._r8, twopi=2._r8*pi
+ win_serod = -1
+
soil_erod_fact = dust_emis_fact
! Summary to log file
@@ -57,27 +67,56 @@ subroutine soil_erod_init( dust_emis_fact, soil_erod_file )
! Get file name.
call getfil(soil_erod_file, infile, 0)
- call cam_pio_openfile (ncid, trim(infile), PIO_NOWRITE)
- ! Get input data resolution.
- ierr = pio_inq_dimid( ncid, 'lon', did )
- ierr = pio_inq_dimlen( ncid, did, nlon )
+ ! Read the source-grid soil erodibility on masterproc only (serial NetCDF) and
+ ! hold it in per-node MPI-3 shared memory: one physical copy per node instead
+ ! of one per MPI rank. At high resolution (e.g. ne1024) this removes
+ ! (ranks_per_node - 1) redundant nlon*nlat copies of the input field.
+ if (masterproc) then
+ iret = nf90_open(trim(infile), NF90_NOWRITE, ncid)
+ if (iret /= nf90_noerr) call endrun('soil_erod_init: failed to open '//trim(infile))
+ iret = nf90_inq_dimid( ncid, 'lon', did )
+ iret = nf90_inquire_dimension( ncid, did, len=nlon )
+ iret = nf90_inq_dimid( ncid, 'lat', did )
+ iret = nf90_inquire_dimension( ncid, did, len=nlat )
+ end if
- ierr = pio_inq_dimid( ncid, 'lat', did )
- ierr = pio_inq_dimlen( ncid, did, nlat )
+#ifdef SPMD
+ call mpibcast( nlon, 1, mpiint, 0, mpicom )
+ call mpibcast( nlat, 1, mpiint, 0, mpicom )
+#endif
allocate(dst_lons(nlon))
allocate(dst_lats(nlat))
- allocate(soil_erodibility_in(nlon,nlat))
+ call cam_shmem_alloc_r8_2d( soil_erodibility_in, win_serod, nlon, nlat )
+ if (masterproc) then
+ write(iulog,*) 'soil_erod_mod: soil_erodibility_in held in per-node shared memory; ', &
+ cam_shmem_npes_per_node()-1, ' redundant copies/node avoided'
+ end if
- ierr = pio_inq_varid( ncid, 'lon', vid )
- ierr = pio_get_var( ncid, vid, dst_lons )
+ ! Open the window epoch, fill on the node leader (masterproc), then publish.
+ call cam_shmem_fence( win_serod )
- ierr = pio_inq_varid( ncid, 'lat', vid )
- ierr = pio_get_var( ncid, vid, dst_lats )
+ if (masterproc) then
+ iret = nf90_inq_varid( ncid, 'lon', vid )
+ iret = nf90_get_var( ncid, vid, dst_lons )
+ iret = nf90_inq_varid( ncid, 'lat', vid )
+ iret = nf90_get_var( ncid, vid, dst_lats )
+ iret = nf90_inq_varid( ncid, 'mbl_bsn_fct_geo', vid )
+ iret = nf90_get_var( ncid, vid, soil_erodibility_in )
+ iret = nf90_close( ncid )
+ end if
- ierr = pio_inq_varid( ncid, 'mbl_bsn_fct_geo', vid )
- ierr = pio_get_var( ncid, vid, soil_erodibility_in )
+#ifdef SPMD
+ ! Source coordinates to every rank; the large field to every node leader. The
+ ! closing fence then publishes it to every rank on each node.
+ call mpibcast( dst_lons, nlon, mpir8, 0, mpicom )
+ call mpibcast( dst_lats, nlat, mpir8, 0, mpicom )
+ if (cam_shmem_is_leader()) then
+ call mpibcast( soil_erodibility_in, nlon*nlat, mpir8, 0, cam_shmem_leader_comm() )
+ end if
+#endif
+ call cam_shmem_fence( win_serod )
!-----------------------------------------------------------------------
! ... convert to radians and setup regridding
@@ -107,11 +146,8 @@ subroutine soil_erod_init( dust_emis_fact, soil_erod_file )
call lininterp_finish(lat_wgts)
call lininterp_finish(lon_wgts)
end do
- deallocate( soil_erodibility_in, stat=ierr )
- if( ierr /= 0 ) then
- write(iulog,*) 'soil_erod_init: failed to deallocate soil_erodibility_in, ierr = ',ierr
- call endrun('soil_erod_init: failed to deallocate soil_erodibility_in')
- end if
+ ! Release the node-shared input field (collective over the node communicator).
+ call cam_shmem_free( soil_erodibility_in, win_serod )
deallocate( dst_lats )
deallocate( dst_lons )
diff --git a/src/chemistry/mozart/chemistry.F90 b/src/chemistry/mozart/chemistry.F90
index 8181e0923e..5d475a2071 100644
--- a/src/chemistry/mozart/chemistry.F90
+++ b/src/chemistry/mozart/chemistry.F90
@@ -1371,6 +1371,7 @@ subroutine chem_final()
use species_sums_diags, only: species_sums_final
use mo_tuvx, only: tuvx_finalize, tuvx_active
use short_lived_species, only: short_lived_species_final
+ use mo_photo, only: photo_final
call mee_ion_final()
call rate_diags_final()
@@ -1379,6 +1380,7 @@ subroutine chem_final()
call tuvx_finalize()
end if
call short_lived_species_final()
+ call photo_final()
end subroutine chem_final
diff --git a/src/chemistry/mozart/mo_jlong.F90 b/src/chemistry/mozart/mo_jlong.F90
index dcc5bf77ae..70d993ffa9 100644
--- a/src/chemistry/mozart/mo_jlong.F90
+++ b/src/chemistry/mozart/mo_jlong.F90
@@ -13,6 +13,11 @@ module mo_jlong
use mpishorthand, only : mpicom,mpiint,mpir8, mpilog, mpir4
#endif
use spmd_utils, only : masterproc
+ use cam_shmem_mod, only : cam_shmem_alloc_r4_4d, cam_shmem_alloc_r4_5d, &
+ cam_shmem_alloc_r8_3d, cam_shmem_fence, &
+ cam_shmem_free, &
+ cam_shmem_is_leader, cam_shmem_leader_comm, &
+ cam_shmem_npes_per_node
implicit none
@@ -23,6 +28,7 @@ module mo_jlong
private
public :: jlong_init
+ public :: jlong_final
public :: jlong_timestep_init
public :: jlong
public :: numj
@@ -42,7 +48,10 @@ module mo_jlong
integer :: numsza ! number of zen angles in rsf
integer :: numalb ! number of albedos in rsf
integer :: numcolo3 ! number of o3 columns in rsf
- real(r4), allocatable :: xsqy(:,:,:,:)
+ ! Large read-only lookup tables held in per-node MPI shared memory
+ ! (one physical copy per node, mapped into every rank). These are
+ ! pointers (not allocatables) so they can alias a shared-memory window.
+ real(r4), pointer :: xsqy(:,:,:,:) => null()
real(r8), allocatable :: wc(:)
real(r8), allocatable :: we(:)
real(r8), allocatable :: wlintv(:)
@@ -50,9 +59,9 @@ module mo_jlong
real(r8), allocatable :: bde_o2_b(:)
real(r8), allocatable :: bde_o3_a(:)
real(r8), allocatable :: bde_o3_b(:)
- real(r8), allocatable :: xs_o2b(:,:,:)
- real(r8), allocatable :: xs_o3a(:,:,:)
- real(r8), allocatable :: xs_o3b(:,:,:)
+ real(r8), pointer :: xs_o2b(:,:,:) => null()
+ real(r8), pointer :: xs_o3a(:,:,:) => null()
+ real(r8), pointer :: xs_o3b(:,:,:) => null()
real(r8), allocatable :: p(:)
real(r8), allocatable :: del_p(:)
real(r8), allocatable :: prs(:)
@@ -64,7 +73,9 @@ module mo_jlong
real(r8), allocatable :: o3rat(:)
real(r8), allocatable :: del_o3rat(:)
real(r8), allocatable :: colo3(:)
- real(r4), allocatable :: rsf_tab(:,:,:,:,:)
+ real(r4), pointer :: rsf_tab(:,:,:,:,:) => null()
+ ! MPI shared-memory window handles for the tables above.
+ integer :: win_xsqy = -1, win_o2b = -1, win_o3a = -1, win_o3b = -1, win_rsf = -1
logical :: jlong_used = .false.
contains
@@ -113,9 +124,50 @@ subroutine jlong_init( xs_long_file, rsf_file, lng_indexer )
endif
jlong_used = .true.
-
+
end subroutine jlong_init
+ subroutine jlong_final
+!------------------------------------------------------------------------------
+! ... release the per-node shared-memory lookup tables and the small
+! per-rank work arrays. Collective over the node communicator
+! (cam_shmem_free -> MPI_Win_free) and called on every rank from the
+! chemistry finalize path, after the last photolysis timestep.
+!------------------------------------------------------------------------------
+ implicit none
+
+ if (.not. jlong_used) return ! jlong_used is identical on all ranks
+
+ call cam_shmem_free( xsqy, win_xsqy )
+ call cam_shmem_free( xs_o2b, win_o2b )
+ call cam_shmem_free( xs_o3a, win_o3a )
+ call cam_shmem_free( xs_o3b, win_o3b )
+ call cam_shmem_free( rsf_tab, win_rsf )
+
+ ! small per-rank allocatables (allocated in get_xsqy / get_rsf)
+ if( allocated(prs) ) deallocate(prs)
+ if( allocated(dprs) ) deallocate(dprs)
+ if( allocated(wc) ) deallocate(wc)
+ if( allocated(we) ) deallocate(we)
+ if( allocated(wlintv)) deallocate(wlintv)
+ if( allocated(etfphot)) deallocate(etfphot)
+ if( allocated(bde_o2_b)) deallocate(bde_o2_b)
+ if( allocated(bde_o3_a)) deallocate(bde_o3_a)
+ if( allocated(bde_o3_b)) deallocate(bde_o3_b)
+ if( allocated(p) ) deallocate(p)
+ if( allocated(del_p) ) deallocate(del_p)
+ if( allocated(sza) ) deallocate(sza)
+ if( allocated(del_sza)) deallocate(del_sza)
+ if( allocated(alb) ) deallocate(alb)
+ if( allocated(del_alb)) deallocate(del_alb)
+ if( allocated(o3rat) ) deallocate(o3rat)
+ if( allocated(del_o3rat)) deallocate(del_o3rat)
+ if( allocated(colo3) ) deallocate(colo3)
+
+ jlong_used = .false.
+
+ end subroutine jlong_final
+
subroutine get_xsqy( xs_long_file, lng_indexer )
!=============================================================================!
! PURPOSE: !
@@ -208,22 +260,37 @@ subroutine get_xsqy( xs_long_file, lng_indexer )
end if
end do
- !------------------------------------------------------------------------------
- ! ... allocate arrays
- !------------------------------------------------------------------------------
+ end if Masterproc_only
- allocate( xsqy(numj,nw,nt,np_xs),stat=iret )
- if( iret /= 0 ) then
- call alloc_err( iret, 'get_xsqy', 'xsqy', numj*nt*np_xs*nw )
- end if
- allocate( prs(np_xs),dprs(np_xs-1),stat=iret )
- if( iret /= 0 ) then
- call alloc_err( iret, 'get_xsqy', 'prs,dprs', np_xs )
- end if
- allocate( xs_o2b(nw,nt,np_xs),xs_o3a(nw,nt,np_xs),xs_o3b(nw,nt,np_xs),stat=iret )
- if( iret /= 0 ) then
- call alloc_err( iret, 'get_xsqy', 'xs_o2b ... xs_o3b', np_xs )
- end if
+#ifdef SPMD
+ call mpibcast( numj, 1, mpiint, 0, mpicom )
+ call mpibcast( nt, 1, mpiint, 0, mpicom )
+ call mpibcast( nw, 1, mpiint, 0, mpicom )
+ call mpibcast( np_xs, 1, mpiint, 0, mpicom )
+#endif
+
+ !------------------------------------------------------------------------------
+ ! ... allocate the large read-only tables in per-node shared memory
+ ! (one copy per node); prs/dprs are small and kept per-rank.
+ !------------------------------------------------------------------------------
+ call cam_shmem_alloc_r4_4d( xsqy, win_xsqy, numj, nw, nt, np_xs )
+ call cam_shmem_alloc_r8_3d( xs_o2b, win_o2b, nw, nt, np_xs )
+ call cam_shmem_alloc_r8_3d( xs_o3a, win_o3a, nw, nt, np_xs )
+ call cam_shmem_alloc_r8_3d( xs_o3b, win_o3b, nw, nt, np_xs )
+ allocate( prs(np_xs), dprs(np_xs-1), stat=iret )
+ if( iret /= 0 ) call alloc_err( iret, 'get_xsqy', 'prs,dprs', np_xs )
+ if( masterproc ) then
+ write(iulog,*) 'mo_jlong: xsqy/xs_o2b/o3a/o3b held in per-node shared memory; ', &
+ cam_shmem_npes_per_node()-1, ' redundant copies/node avoided'
+ end if
+
+ ! Open the window epoch; node leaders then fill the shared tables below.
+ call cam_shmem_fence( win_xsqy )
+ call cam_shmem_fence( win_o2b )
+ call cam_shmem_fence( win_o3a )
+ call cam_shmem_fence( win_o3b )
+
+ Masterproc_read : if( masterproc ) then
!------------------------------------------------------------------------------
! ... read cross sections
!------------------------------------------------------------------------------
@@ -266,43 +333,25 @@ subroutine get_xsqy( xs_long_file, lng_indexer )
iret = nf90_inq_varid( ncid, 'pressure', varid )
iret = nf90_get_var( ncid, varid, prs )
iret = nf90_close( ncid )
- end if Masterproc_only
+ end if Masterproc_read
#ifdef SPMD
-! call mpibarrier( mpicom )
- call mpibcast( numj, 1, mpiint, 0, mpicom )
- call mpibcast( nt, 1, mpiint, 0, mpicom )
- call mpibcast( nw, 1, mpiint, 0, mpicom )
- call mpibcast( np_xs, 1, mpiint, 0, mpicom )
- call mpibcast( lng_indexer, phtcnt, mpiint, 0, mpicom )
-#endif
- if( .not. masterproc ) then
- !------------------------------------------------------------------------------
- ! ... allocate arrays
- !------------------------------------------------------------------------------
- allocate( xsqy(numj,nw,nt,np_xs),stat=iret )
- if( iret /= nf90_noerr) then
- write(iulog,*) 'get_xsqy : failed to allocate xsqy ; error = ',iret
- call endrun
- end if
- allocate( prs(np_xs),dprs(np_xs-1),stat=iret )
- if( iret /= nf90_noerr) then
- write(iulog,*) 'get_xsqy : failed to allocate prs,dprs ; error = ',iret
- call endrun
- end if
- allocate( xs_o2b(nw,nt,np_xs),xs_o3a(nw,nt,np_xs),xs_o3b(nw,nt,np_xs),stat=iret )
- if( iret /= 0 ) then
- call alloc_err( iret, 'get_xsqy', 'xs_o2b ... xs_o3b', np_xs )
- end if
+ ! Distribute the tables to the leader of every other node (leaders only),
+ ! then publish to every rank on each node via the closing window fence.
+ if( cam_shmem_is_leader() ) then
+ call mpibcast( xsqy, numj*nt*np_xs*nw, mpir4, 0, cam_shmem_leader_comm() )
+ call mpibcast( xs_o2b, nw*nt*np_xs, mpir8, 0, cam_shmem_leader_comm() )
+ call mpibcast( xs_o3a, nw*nt*np_xs, mpir8, 0, cam_shmem_leader_comm() )
+ call mpibcast( xs_o3b, nw*nt*np_xs, mpir8, 0, cam_shmem_leader_comm() )
end if
-#ifdef SPMD
-! call mpibarrier( mpicom )
- call mpibcast( prs, np_xs, mpir8, 0, mpicom )
- call mpibcast( xsqy, numj*nt*np_xs*nw, mpir4, 0, mpicom )
- call mpibcast( xs_o2b, nw*nt*np_xs, mpir8, 0, mpicom )
- call mpibcast( xs_o3a, nw*nt*np_xs, mpir8, 0, mpicom )
- call mpibcast( xs_o3b, nw*nt*np_xs, mpir8, 0, mpicom )
+ call mpibcast( lng_indexer, phtcnt, mpiint, 0, mpicom )
+ call mpibcast( prs, np_xs, mpir8, 0, mpicom )
#endif
+ ! Publish the filled tables to every rank on each node.
+ call cam_shmem_fence( win_xsqy )
+ call cam_shmem_fence( win_o2b )
+ call cam_shmem_fence( win_o3a )
+ call cam_shmem_fence( win_o3b )
dprs(:np_xs-1) = 1._r8/(prs(1:np_xs-1) - prs(2:np_xs))
end subroutine get_xsqy
@@ -415,11 +464,14 @@ subroutine get_rsf(rsf_file)
if( iret /= 0 ) then
call alloc_err( iret, 'get_rsf', 'colo3', nump )
end if
- allocate( rsf_tab(nw,nump,numsza,numcolo3,numalb),stat=iret )
- if( iret /= 0 ) then
- write(iulog,*) 'get_rsf : dimensions = ',nw,nump,numsza,numcolo3,numalb
- call alloc_err( iret, 'get_rsf', 'rsf_tab', numalb*numcolo3*numsza*nump )
+ ! rsf_tab is large and read-only; hold one copy per node in shared memory.
+ call cam_shmem_alloc_r4_5d( rsf_tab, win_rsf, nw, nump, numsza, numcolo3, numalb )
+ if( masterproc ) then
+ write(iulog,*) 'mo_jlong: rsf_tab held in per-node shared memory; ', &
+ cam_shmem_npes_per_node()-1, ' redundant copies/node avoided'
end if
+ ! Open the window epoch before the node leaders fill rsf_tab below.
+ call cam_shmem_fence( win_rsf )
Masterproc_only2 : if( masterproc ) then
!------------------------------------------------------------------------------
@@ -470,10 +522,15 @@ subroutine get_rsf(rsf_file)
call mpibcast( alb, numalb, mpir8, 0, mpicom )
call mpibcast( o3rat, numcolo3, mpir8, 0, mpicom )
call mpibcast( colo3, nump, mpir8, 0, mpicom )
- do w = 1,nw
- call mpibcast( rsf_tab(w,:,:,:,:), numalb*numcolo3*numsza*nump, mpir4, 0, mpicom )
- enddo
+ ! Distribute rsf_tab to node leaders only (per wavelength slice), then publish.
+ if( cam_shmem_is_leader() ) then
+ do w = 1,nw
+ call mpibcast( rsf_tab(w,:,:,:,:), numalb*numcolo3*numsza*nump, mpir4, 0, cam_shmem_leader_comm() )
+ enddo
+ end if
#endif
+ ! Publish the filled rsf_tab to every rank on each node.
+ call cam_shmem_fence( win_rsf )
#ifdef USE_BDE
if (masterproc) write(iulog,*) 'Jlong using bdes'
bde_o2_b(:) = max( 0._r8, hc*(wc_o2_b - wc(:))/(wc_o2_b*wc(:)) )
diff --git a/src/chemistry/mozart/mo_photo.F90 b/src/chemistry/mozart/mo_photo.F90
index 01ca12c06b..b032a6b7ab 100644
--- a/src/chemistry/mozart/mo_photo.F90
+++ b/src/chemistry/mozart/mo_photo.F90
@@ -18,7 +18,7 @@ module mo_photo
private
- public :: photo_inti, table_photo
+ public :: photo_inti, photo_final, table_photo
public :: set_ub_col
public :: setcol
public :: photo_timestep_init
@@ -534,6 +534,18 @@ subroutine photo_inti( xs_coef_file, xs_short_file, xs_long_file, rsf_file, &
end subroutine photo_inti
+ subroutine photo_final
+ !----------------------------------------------------------------------
+ ! ... release the photolysis lookup tables held in per-node shared
+ ! memory. Mirrors photo_inti -> jlong_init; called on every rank
+ ! from chem_final, before MPI shutdown.
+ !----------------------------------------------------------------------
+ use mo_jlong, only : jlong_final
+
+ call jlong_final
+
+ end subroutine photo_final
+
subroutine table_photo( photos, pmid, pdel, temper, zmid, zint, &
col_dens, zen_angle, srf_alb, lwc, clouds, &
esfact, vmr, invariants, ncol, lchnk, pbuf )
diff --git a/src/chemistry/utils/tracer_data.F90 b/src/chemistry/utils/tracer_data.F90
index ce6843eafa..310e03d6b5 100644
--- a/src/chemistry/utils/tracer_data.F90
+++ b/src/chemistry/utils/tracer_data.F90
@@ -26,6 +26,8 @@ module tracer_data
pio_get_var, pio_get_att, pio_nowrite, pio_inq_dimlen, &
pio_inq_vardimid, pio_inq_dimlen, pio_closefile, &
pio_inquire_variable
+ use cam_shmem_mod, only : cam_shmem_alloc_r8_3d, cam_shmem_fence, &
+ cam_shmem_free, cam_shmem_is_leader, cam_shmem_npes_per_node
implicit none
@@ -133,6 +135,20 @@ module tracer_data
logical :: top_bndry = .false.
logical :: top_layer = .false.
logical :: stepTime = .false. ! Do not interpolate in time, but use stepwise times
+ ! Bookkeeping for in-memory carry-forward of the overlapping time slice.
+ ! When the bracketing window advances by exactly one record, the new lower
+ ! record equals the previous upper record already resident in input(2), so
+ ! it can be copied forward instead of re-read from disk.
+ integer :: loaded_recno(2) = -1 ! file record numbers resident in input(1:2)
+ integer :: loaded_fh(2) = -1 ! PIO file handle (curr vs next) owning each
+ ! Per-node MPI shared-memory buffer holding ONE global source field
+ ! (nlon x nlat x nlev) in flight during a read. The node leader reads the
+ ! field once (serial netCDF) and every rank on the node interpolates from
+ ! this single shared copy, instead of all ranks redundantly reading and
+ ! storing the full field. Allocated once in trcdata_init for lat-lon
+ ! (non-unstructured, non-zonal) files; reused across fields and updates.
+ real(r8), pointer, dimension(:,:,:) :: shr_src => null()
+ integer :: win_shr_src = -1
endtype trfile
integer, public, parameter :: MAXTRCRS = 100
@@ -393,6 +409,24 @@ subroutine trcdata_init( specifier, filename, filelist, datapath, flds, file, &
! Hackish workaround is to make a copy...
lev_dimid = old_dimid
+ ! Allocate the per-node shared source-field buffer used by read_2d_trc /
+ ! read_3d_trc. Only lat-lon files (which go through the redundant
+ ! "every rank reads the whole field then interpolates" path) need it; files
+ ! already on the unstructured physics grid (read_physgrid_*) and zonal-average
+ ! files (read_za_trc) do not. One physical copy per node; reused across all
+ ! fields of this file and all subsequent updates (freed at MPI_Finalize).
+ if ( .not.file%unstructured .and. .not.file%zonal_ave .and. &
+ .not.associated(file%shr_src) ) then
+ call cam_shmem_alloc_r8_3d( file%shr_src, file%win_shr_src, &
+ file%nlon, file%nlat, max(file%nlev,1) )
+ if (masterproc) then
+ write(iulog,*) 'trcdata_init: ', trim(file%curr_filename), &
+ ' source field in per-node shared memory (', &
+ file%nlon,'x',file%nlat,'x',max(file%nlev,1),'); ', &
+ cam_shmem_npes_per_node()-1,' redundant copies/node avoided'
+ end if
+ end if
+
if (file%has_ps) then
allocate( file%hyam(file%nlev), file%hybm(file%nlev), stat=astat )
@@ -1285,6 +1319,7 @@ subroutine read_next_trcdata( flds, file )
integer :: strt3(3) ! array of starting indices
type(file_desc_t) :: fids(4)
logical :: times_found
+ logical :: carry ! .true. when input(1) can be carried from the resident input(2)
integer :: cur_yr, cur_mon, cur_day, cur_sec, yr1, yr2, mon, date, sec
real(r8) :: series1_time, series2_time
@@ -1367,8 +1402,28 @@ subroutine read_next_trcdata( flds, file )
! Set up hyperslab corners
!
+ ! If the bracketing window advanced by exactly one record, the new lower
+ ! record (recnos(1)) equals the previous upper record already resident in
+ ! input(2)/ps_in(2). Carry it forward in memory instead of re-reading it
+ ! from disk (this halves the per-update forcing I/O, time-axis broadcast and
+ ! interpolation work). Only the 2-record path is eligible.
+ carry = ( file%interp_recs == 2 .and. file%initialized .and. &
+ recnos(1) == file%loaded_recno(2) .and. fids(1)%fh == file%loaded_fh(2) )
+
do i=1,file%interp_recs
+ if ( carry .and. i == 1 ) then
+ ! input(2) -> input(1) (and ps_in(2) -> ps_in(1)) is bit-identical to
+ ! the value the old code re-read for this record.
+ do f = 1,nflds
+ flds(f)%input(1)%data(:,:,:) = flds(f)%input(2)%data(:,:,:)
+ end do
+ if ( file%has_ps ) then
+ file%ps_in(1)%data(:,:) = file%ps_in(2)%data(:,:)
+ end if
+ cycle
+ end if
+
strt4(:) = 1
strt3(:) = 1
@@ -1442,6 +1497,17 @@ subroutine read_next_trcdata( flds, file )
enddo
+ ! Remember which records are now resident so the next update can carry the
+ ! overlapping slice forward. Only the 2-record path supports carry; other
+ ! paths (stepTime=1 record, fill_in_months=4 records) invalidate it.
+ if ( file%interp_recs == 2 ) then
+ file%loaded_recno(1:2) = recnos(1:2)
+ file%loaded_fh(1:2) = (/ fids(1)%fh, fids(2)%fh /)
+ else
+ file%loaded_recno(:) = -1
+ file%loaded_fh(:) = -1
+ end if
+
end subroutine read_next_trcdata
!------------------------------------------------------------------------
@@ -1455,6 +1521,9 @@ subroutine read_2d_trc( fid, vid, loc_arr, strt, cnt, file, order )
use dycore, only: dycore_is
use polar_avg, only: polar_average
use horizontal_interpolate, only : xy_interp
+ use ioFileMod, only: getfil
+ use netcdf, only: nf90_open, nf90_close, nf90_get_var, &
+ nf90_inq_varid, nf90_nowrite, nf90_noerr
implicit none
type(file_desc_t), intent(in) :: fid
@@ -1472,29 +1541,78 @@ subroutine read_2d_trc( fid, vid, loc_arr, strt, cnt, file, order )
type(interp_type) :: lon_wgts, lat_wgts
integer :: lons(pcols), lats(pcols)
real(r8) :: file_lats(file%nlat)
+ logical :: use_shmem
+ character(len=shr_kind_cl) :: srcname, srcpath, locfn, varname
+ integer :: nfid, nfvid, iret
nullify(wrk2d_in)
- allocate( wrk2d(cnt(1),cnt(2)), stat=ierr )
- if( ierr /= 0 ) then
- write(iulog,*) 'read_2d_trc: wrk2d allocation error = ',ierr
- call endrun
- end if
- if(order(1)/=1 .or. order(2)/=2 .or. cnt(1)/=file%nlon .or. cnt(2)/=file%nlat) then
- allocate( wrk2d_in(file%nlon, file%nlat), stat=ierr )
+ ! Per-node shared read: leader reads the whole field once (serial netCDF)
+ ! into the shared buffer; every rank interpolates from that single copy.
+ ! Falls back to the per-rank PIO read when the buffer is absent.
+ use_shmem = associated(file%shr_src)
+
+ if (use_shmem) then
+ iret = pio_inquire_variable(fid, vid, name=varname)
+ if (fid%fh == file%curr_fileid%fh) then
+ srcname = file%curr_filename
+ else
+ srcname = file%next_filename
+ end if
+ if (len_trim(file%pathname) == 0) then
+ srcpath = trim(srcname)
+ else
+ srcpath = trim(file%pathname) // '/' // trim(srcname)
+ end if
+
+ call cam_shmem_fence(file%win_shr_src)
+ if (cam_shmem_is_leader()) then
+ call getfil( srcpath, locfn, 0 )
+ iret = nf90_open( trim(locfn), nf90_nowrite, nfid )
+ if (iret /= nf90_noerr) call endrun('read_2d_trc: nf90_open failed: '//trim(locfn))
+ iret = nf90_inq_varid( nfid, trim(varname), nfvid )
+ if (iret /= nf90_noerr) call endrun('read_2d_trc: nf90_inq_varid failed: '//trim(varname))
+ allocate( wrk2d(cnt(1),cnt(2)), stat=ierr )
+ if( ierr /= 0 ) then
+ write(iulog,*) 'read_2d_trc: wrk2d allocation error = ',ierr
+ call endrun
+ end if
+ iret = nf90_get_var( nfid, nfvid, wrk2d, start=strt, count=cnt )
+ if (iret /= nf90_noerr) call endrun('read_2d_trc: nf90_get_var failed: '//trim(varname))
+ iret = nf90_close( nfid )
+ if(order(1)/=1 .or. order(2)/=2 .or. cnt(1)/=file%nlon .or. cnt(2)/=file%nlat) then
+ file%shr_src(:,:,1) = reshape( wrk2d(:,:),(/file%nlon,file%nlat/), order=order )
+ else
+ file%shr_src(:,:,1) = wrk2d(:,:)
+ end if
+ deallocate(wrk2d)
+ end if
+ call cam_shmem_fence(file%win_shr_src)
+ wrk2d_in => file%shr_src(:,:,1)
+
+ else
+ allocate( wrk2d(cnt(1),cnt(2)), stat=ierr )
if( ierr /= 0 ) then
- write(iulog,*) 'read_2d_trc: wrk2d_in allocation error = ',ierr
+ write(iulog,*) 'read_2d_trc: wrk2d allocation error = ',ierr
call endrun
end if
- end if
+ if(order(1)/=1 .or. order(2)/=2 .or. cnt(1)/=file%nlon .or. cnt(2)/=file%nlat) then
+ allocate( wrk2d_in(file%nlon, file%nlat), stat=ierr )
+ if( ierr /= 0 ) then
+ write(iulog,*) 'read_2d_trc: wrk2d_in allocation error = ',ierr
+ call endrun
+ end if
+ end if
- ierr = pio_get_var( fid, vid, strt, cnt, wrk2d )
- if(associated(wrk2d_in)) then
- wrk2d_in = reshape( wrk2d(:,:),(/file%nlon,file%nlat/), order=order )
- deallocate(wrk2d)
- else
- wrk2d_in => wrk2d
+
+ ierr = pio_get_var( fid, vid, strt, cnt, wrk2d )
+ if(associated(wrk2d_in)) then
+ wrk2d_in = reshape( wrk2d(:,:),(/file%nlon,file%nlat/), order=order )
+ deallocate(wrk2d)
+ else
+ wrk2d_in => wrk2d
+ end if
end if
! PGI 13.9 bug workaround.
@@ -1559,10 +1677,12 @@ subroutine read_2d_trc( fid, vid, loc_arr, strt, cnt, file, order )
end if
- if(allocated(wrk2d)) then
- deallocate(wrk2d)
- else
- deallocate(wrk2d_in)
+ if (.not. use_shmem) then
+ if(allocated(wrk2d)) then
+ deallocate(wrk2d)
+ else
+ deallocate(wrk2d_in)
+ end if
end if
if(dycore_is('LR')) call polar_average(loc_arr)
end subroutine read_2d_trc
@@ -1706,6 +1826,9 @@ subroutine read_3d_trc( fid, vid, loc_arr, strt, cnt, file, order)
use dycore, only : dycore_is
use polar_avg, only : polar_average
use horizontal_interpolate, only : xy_interp
+ use ioFileMod, only : getfil
+ use netcdf, only : nf90_open, nf90_close, nf90_get_var, &
+ nf90_inq_varid, nf90_nowrite, nf90_noerr
implicit none
@@ -1726,28 +1849,79 @@ subroutine read_3d_trc( fid, vid, loc_arr, strt, cnt, file, order)
real(r8) :: to_lons(pcols), to_lats(pcols)
real(r8), parameter :: zero=0_r8, twopi=2_r8*pi
type(interp_type) :: lon_wgts, lat_wgts
+ logical :: use_shmem
+ character(len=shr_kind_cl) :: srcname, srcpath, locfn, varname
+ integer :: nfid, nfvid, iret
loc_arr(:,:,:) = 0._r8
nullify(wrk3d_in)
- allocate(wrk3d(cnt(1),cnt(2),cnt(3)), stat=ierr)
- if( ierr /= 0 ) then
- write(iulog,*) 'read_3d_trc: wrk3d allocation error = ',ierr
- call endrun
- end if
- ierr = pio_get_var( fid, vid, strt, cnt, wrk3d )
+ ! When the per-node shared source buffer is present (lat-lon files), the
+ ! node leader reads the whole field once (serial netCDF) into the shared
+ ! buffer and every rank interpolates from that single copy. Otherwise fall
+ ! back to the original path where every rank reads the field via PIO.
+ use_shmem = associated(file%shr_src)
+
+ if (use_shmem) then
+ iret = pio_inquire_variable(fid, vid, name=varname)
+ if (fid%fh == file%curr_fileid%fh) then
+ srcname = file%curr_filename
+ else
+ srcname = file%next_filename
+ end if
+ if (len_trim(file%pathname) == 0) then
+ srcpath = trim(srcname)
+ else
+ srcpath = trim(file%pathname) // '/' // trim(srcname)
+ end if
- if(order(1)/=1 .or. order(2)/=2 .or. order(3)/=3 .or. &
- cnt(1)/=file%nlon.or.cnt(2)/=file%nlat.or.cnt(3)/=file%nlev) then
- allocate(wrk3d_in(file%nlon,file%nlat,file%nlev),stat=ierr)
+ call cam_shmem_fence(file%win_shr_src) ! open the window epoch
+ if (cam_shmem_is_leader()) then
+ call getfil( srcpath, locfn, 0 )
+ iret = nf90_open( trim(locfn), nf90_nowrite, nfid )
+ if (iret /= nf90_noerr) call endrun('read_3d_trc: nf90_open failed: '//trim(locfn))
+ iret = nf90_inq_varid( nfid, trim(varname), nfvid )
+ if (iret /= nf90_noerr) call endrun('read_3d_trc: nf90_inq_varid failed: '//trim(varname))
+ allocate(wrk3d(cnt(1),cnt(2),cnt(3)), stat=ierr)
+ if( ierr /= 0 ) then
+ write(iulog,*) 'read_3d_trc: wrk3d allocation error = ',ierr
+ call endrun
+ end if
+ iret = nf90_get_var( nfid, nfvid, wrk3d, start=strt, count=cnt )
+ if (iret /= nf90_noerr) call endrun('read_3d_trc: nf90_get_var failed: '//trim(varname))
+ iret = nf90_close( nfid )
+ if(order(1)/=1 .or. order(2)/=2 .or. order(3)/=3 .or. &
+ cnt(1)/=file%nlon.or.cnt(2)/=file%nlat.or.cnt(3)/=file%nlev) then
+ file%shr_src(:,:,:) = reshape( wrk3d(:,:,:),(/file%nlon,file%nlat,file%nlev/), order=order )
+ else
+ file%shr_src(:,:,:) = wrk3d(:,:,:)
+ end if
+ deallocate(wrk3d)
+ end if
+ call cam_shmem_fence(file%win_shr_src) ! publish to all ranks on node
+ wrk3d_in => file%shr_src
+
+ else
+ allocate(wrk3d(cnt(1),cnt(2),cnt(3)), stat=ierr)
if( ierr /= 0 ) then
write(iulog,*) 'read_3d_trc: wrk3d allocation error = ',ierr
call endrun
end if
- wrk3d_in = reshape( wrk3d(:,:,:),(/file%nlon,file%nlat,file%nlev/), order=order )
- deallocate(wrk3d)
- else
- wrk3d_in => wrk3d
+
+ ierr = pio_get_var( fid, vid, strt, cnt, wrk3d )
+
+ if(order(1)/=1 .or. order(2)/=2 .or. order(3)/=3 .or. &
+ cnt(1)/=file%nlon.or.cnt(2)/=file%nlat.or.cnt(3)/=file%nlev) then
+ allocate(wrk3d_in(file%nlon,file%nlat,file%nlev),stat=ierr)
+ if( ierr /= 0 ) then
+ write(iulog,*) 'read_3d_trc: wrk3d allocation error = ',ierr
+ call endrun
+ end if
+ wrk3d_in = reshape( wrk3d(:,:,:),(/file%nlon,file%nlat,file%nlev/), order=order )
+ deallocate(wrk3d)
+ else
+ wrk3d_in => wrk3d
+ end if
end if
! If weighting by latitude, then perform horizontal interpolation by using weight_x, weight_y
@@ -1796,15 +1970,17 @@ subroutine read_3d_trc( fid, vid, loc_arr, strt, cnt, file, order)
end do
endif
- if(allocated(wrk3d)) then
- deallocate( wrk3d, stat=astat )
- else
- deallocate( wrk3d_in, stat=astat )
+ if (.not. use_shmem) then
+ if(allocated(wrk3d)) then
+ deallocate( wrk3d, stat=astat )
+ else
+ deallocate( wrk3d_in, stat=astat )
+ end if
+ if( astat/= 0 ) then
+ write(iulog,*) 'read_3d_trc: failed to deallocate wrk3d array; error = ',astat
+ call endrun
+ endif
end if
- if( astat/= 0 ) then
- write(iulog,*) 'read_3d_trc: failed to deallocate wrk3d array; error = ',astat
- call endrun
- endif
if(dycore_is('LR')) call polar_average(file%nlev, loc_arr)
end subroutine read_3d_trc
diff --git a/src/control/cam_comp.F90 b/src/control/cam_comp.F90
index a040762067..ead822a993 100644
--- a/src/control/cam_comp.F90
+++ b/src/control/cam_comp.F90
@@ -28,6 +28,7 @@ module cam_comp
use perf_mod
use cam_logfile, only: iulog
use cam_abortutils, only: endrun
+use cam_shmem_mod, only: cam_shmem_init
implicit none
private
@@ -70,7 +71,7 @@ subroutine cam_init( &
!
!-----------------------------------------------------------------------
- use cam_initfiles, only: cam_initfiles_open
+ use cam_initfiles, only: cam_initfiles_open, cam_initfiles_close
use dyn_grid, only: dyn_grid_init
use phys_grid, only: phys_grid_init
use physpkg, only: phys_register, phys_init
@@ -132,6 +133,12 @@ subroutine cam_init( &
call init_pio_subsystem()
+ ! Build the cam_shmem per-node/leader communicators now, at this early all-ranks
+ ! collective point, so the first MPI_Comm_split_type(SHARED) is not paid cold mid-
+ ! physprop (rad_cnst_init). Needed for the theta-l dycore, which (unlike
+ ! theta-l_kokkos/COMPOSE) does no early SHARED comm-split to warm it.
+ call cam_shmem_init()
+
! Initializations using data passed from coupler.
call cam_ctrl_init( &
caseid_in=caseid, &
@@ -195,6 +202,13 @@ subroutine cam_init( &
call phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
+ ! Initial-condition and topo files are fully read by the end of phys_init
+ ! (which is their last use); close them now to free PIO/NetCDF file handles
+ ! for the rest of the run rather than holding them open until cam_final.
+ if (initial_run_in) then
+ call cam_initfiles_close()
+ end if
+
call stepon_init(dyn_in, dyn_out)
call offline_driver_init()
@@ -414,10 +428,8 @@ subroutine cam_final( cam_out, cam_in )
!-----------------------------------------------------------------------
use stepon, only: stepon_final
use physpkg, only: phys_final
- use cam_initfiles, only: cam_initfiles_close
use camsrfexch, only: atm2hub_deallocate, hub2atm_deallocate
use ionosphere_interface, only: ionosphere_final
- use cam_control_mod, only: initial_run
!
! Arguments
@@ -433,9 +445,8 @@ subroutine cam_final( cam_out, cam_in )
call stepon_final(dyn_in, dyn_out)
call ionosphere_final()
- if (initial_run) then
- call cam_initfiles_close()
- end if
+ ! Note: the initial-condition and topo files (fh_ini/fh_topo) are now closed at
+ ! the end of cam_init, immediately after their last use in phys_init.
call hub2atm_deallocate(cam_in)
call atm2hub_deallocate(cam_out)
diff --git a/src/cpl/nuopc/atm_stream_ndep.F90 b/src/cpl/nuopc/atm_stream_ndep.F90
index f54509b269..e8a063353b 100644
--- a/src/cpl/nuopc/atm_stream_ndep.F90
+++ b/src/cpl/nuopc/atm_stream_ndep.F90
@@ -43,6 +43,7 @@ module atm_stream_ndep
integer :: stream_ndep_year_first ! first year in stream to use
integer :: stream_ndep_year_last ! last year in stream to use
integer :: stream_ndep_year_align ! align stream_year_firstndep with
+ character(len=CS) :: stream_ndep_mapalgo ! stream->model mapping: bilinear|redist|nn|consf|consd|none
!==============================================================================
contains
@@ -66,6 +67,7 @@ subroutine stream_ndep_readnl(nlfile)
namelist /ndep_stream_nl/ &
stream_ndep_data_filename, &
stream_ndep_mesh_filename, &
+ stream_ndep_mapalgo, &
stream_ndep_year_first, &
stream_ndep_year_last, &
stream_ndep_year_align
@@ -76,6 +78,15 @@ subroutine stream_ndep_readnl(nlfile)
stream_ndep_year_first = 1 ! first year in stream to use
stream_ndep_year_last = 1 ! last year in stream to use
stream_ndep_year_align = 1 ! align stream_ndep_year_first with this model year
+ ! Spatial mapping from the stream grid to the model grid. Initialized to an
+ ! invalid sentinel rather than a usable value so that an unset namelist variable
+ ! fails loudly instead of running on a hidden compiled-in default; build-namelist
+ ! supplies the value (default 'bilinear'). Use 'bilinear' when the ndep file is
+ ! on a different (coarser) grid than the model. Use 'redist' (pure index-based
+ ! redistribution -- no interpolation weights, no ESMF route handle) when the ndep
+ ! file is already on the model grid; this avoids building an ESMF regrid over the
+ ! full grid and is far cheaper in memory at high resolution (e.g. ne1024pg2).
+ stream_ndep_mapalgo = 'UNSET'
! For now variable list in stream data file is hard-wired
stream_varlist_ndep = (/'NDEP_NHx_month', 'NDEP_NOy_month'/)
@@ -105,6 +116,8 @@ subroutine stream_ndep_readnl(nlfile)
if (ierr /= 0) call endrun(trim(subname)//": FATAL: mpi_bcast: stream_ndep_year_last")
call mpi_bcast(stream_ndep_year_align, 1, mpi_integer, 0, mpicom, ierr)
if (ierr /= 0) call endrun(trim(subname)//": FATAL: mpi_bcast: stream_ndep_year_align")
+ call mpi_bcast(stream_ndep_mapalgo, len(stream_ndep_mapalgo), mpi_character, 0, mpicom, ierr)
+ if (ierr /= 0) call endrun(trim(subname)//": FATAL: mpi_bcast: stream_ndep_mapalgo")
ndep_stream_active = len_trim(stream_ndep_data_filename)>0 .and. stream_ndep_data_filename/='UNSET'
@@ -118,11 +131,21 @@ subroutine stream_ndep_readnl(nlfile)
return
endif
+ ! The ndep stream is active, so stream_ndep_mapalgo must have been supplied by
+ ! build-namelist (add_default, alongside stream_ndep_data_filename). Still holding
+ ! the sentinel means the namelist is inconsistent with an active ndep stream; fail
+ ! here rather than silently mapping with an unintended algorithm.
+ if (trim(stream_ndep_mapalgo) == 'UNSET') then
+ call endrun(trim(subname)//": FATAL: ndep stream is active but stream_ndep_mapalgo"// &
+ " is not set"//errMsg(sourcefile, __LINE__))
+ end if
+
if (masterproc) then
write(iulog,'(a)' ) ' '
write(iulog,'(a,i8)') 'stream ndep settings:'
write(iulog,'(a,a)' ) ' stream_ndep_data_filename = ',trim(stream_ndep_data_filename)
write(iulog,'(a,a)' ) ' stream_ndep_mesh_filename = ',trim(stream_ndep_mesh_filename)
+ write(iulog,'(a,a)' ) ' stream_ndep_mapalgo = ',trim(stream_ndep_mapalgo)
write(iulog,'(a,a,a)') ' stream_varlist_ndep = ',trim(stream_varlist_ndep(1)), trim(stream_varlist_ndep(2))
write(iulog,'(a,i8)') ' stream_ndep_year_first = ',stream_ndep_year_first
write(iulog,'(a,i8)') ' stream_ndep_year_last = ',stream_ndep_year_last
@@ -152,28 +175,35 @@ subroutine stream_ndep_init(model_mesh, model_clock, rc)
! Read in units
call stream_ndep_check_units(stream_ndep_data_filename)
- ! Initialize the cdeps data type sdat_ndep
- call shr_strdata_init_from_inline(sdat_ndep, &
- my_task = iam, &
- logunit = iulog, &
- compname = 'ATM', &
- model_clock = model_clock, &
- model_mesh = model_mesh, &
- stream_meshfile = trim(stream_ndep_mesh_filename), &
- stream_filenames = (/trim(stream_ndep_data_filename)/), &
- stream_yearFirst = stream_ndep_year_first, &
- stream_yearLast = stream_ndep_year_last, &
- stream_yearAlign = stream_ndep_year_align, &
- stream_fldlistFile = stream_varlist_ndep, &
- stream_fldListModel = stream_varlist_ndep, &
- stream_lev_dimname = 'null', &
- stream_mapalgo = 'bilinear', &
- stream_offset = 0, &
- stream_taxmode = 'cycle', &
- stream_dtlimit = 1.0e30_r8, &
- stream_tintalgo = 'linear', &
- stream_name = 'Nitrogen deposition data ', &
- rc = rc)
+ ! Initialize the cdeps data type sdat_ndep.
+ ! When the ndep file is already on the model grid, set mapalgo='redist' in the
+ ! namelist: CDEPS then reuses the already-built model mesh for the stream rather
+ ! than constructing a duplicate full ESMF mesh from stream_ndep_mesh_filename --
+ ! at ne1024pg2 that duplicate mesh (1.6 GB file read + mesh build) is a large
+ ! init memory/time cost. For other mapalgos (e.g. 'bilinear' with a coarser file)
+ ! the stream is on a different grid, so CDEPS still creates the stream mesh from
+ ! its own file.
+ call shr_strdata_init_from_inline(sdat_ndep, &
+ my_task = iam, &
+ logunit = iulog, &
+ compname = 'ATM', &
+ model_clock = model_clock, &
+ model_mesh = model_mesh, &
+ stream_meshfile = trim(stream_ndep_mesh_filename), &
+ stream_filenames = (/trim(stream_ndep_data_filename)/), &
+ stream_yearFirst = stream_ndep_year_first, &
+ stream_yearLast = stream_ndep_year_last, &
+ stream_yearAlign = stream_ndep_year_align, &
+ stream_fldlistFile = stream_varlist_ndep, &
+ stream_fldListModel = stream_varlist_ndep, &
+ stream_lev_dimname = 'null', &
+ stream_mapalgo = trim(stream_ndep_mapalgo), &
+ stream_offset = 0, &
+ stream_taxmode = 'cycle', &
+ stream_dtlimit = 1.0e30_r8, &
+ stream_tintalgo = 'linear', &
+ stream_name = 'Nitrogen deposition data ', &
+ rc = rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) then
call ESMF_Finalize(endflag=ESMF_END_ABORT)
end if
diff --git a/src/physics/cam/phys_prop.F90 b/src/physics/cam/phys_prop.F90
index 8175f46380..feb369a922 100644
--- a/src/physics/cam/phys_prop.F90
+++ b/src/physics/cam/phys_prop.F90
@@ -20,6 +20,10 @@ module phys_prop
use cam_logfile, only: iulog
use cam_abortutils, only: endrun
+use cam_shmem_mod, only: cam_shmem_alloc_r8_2d, cam_shmem_alloc_r8_4d, &
+ cam_shmem_alloc_r8_5d, cam_shmem_free, &
+ cam_shmem_fence, cam_shmem_is_leader, &
+ cam_shmem_npes_per_node
implicit none
private
@@ -31,6 +35,7 @@ module phys_prop
physprop_accum_unique_files, &! Make a list of the unique set of files that contain properties
! This is an initialization step that must be done before calling physprop_init
physprop_init, &! Initialization -- read the input datasets
+ physprop_final, &! Finalization -- free per-node shared-memory tables
physprop_get_id, &! Return ID used to access the property data from the input files
physprop_get ! Return data for specified ID
@@ -120,12 +125,24 @@ module phys_prop
real(r8) :: rhcrystal ! crystalization relative humidity for mode
real(r8) :: rhdeliques ! deliquescence relative humidity for mode
+ ! MPI-3 shared-memory window handles for the large tables above that are held
+ ! one-copy-per-node (see cam_shmem_mod). -1 means "not shared" (table either
+ ! unused for this entry's optics method, or kept as an ordinary allocate).
+ integer :: win_extpsw = -1, win_abspsw = -1, win_asmpsw = -1, win_absplw = -1 ! modal
+ integer :: win_sw_hygro_ext = -1, win_sw_hygro_ssa = -1 ! hygro(scopic)
+ integer :: win_sw_hygro_asm = -1, win_lw_hygro_abs = -1
+ integer :: win_sw_hygro_ext_wtp = -1, win_sw_hygro_ssa_wtp = -1 ! hygroscopic_wtp
+ integer :: win_sw_hygro_asm_wtp = -1, win_lw_hygro_abs_wtp = -1
+ integer :: win_cs_ext = -1, win_cs_ssa = -1, win_cs_asm = -1, win_cs_abs = -1 ! coreshell
+ integer :: win_r_sw_ext = -1, win_r_sw_scat = -1 ! volcanic_radius
+ integer :: win_r_sw_ascat = -1, win_r_lw_abs = -1
+
endtype physprop_type
! This module stores data in an array of physprop_type structures. The way this data
! is accessed outside the module is via a physprop ID, which is an index into the array.
integer :: numphysprops = 0 ! an incremental total across ALL clim and diag constituents
-type (physprop_type), pointer :: physprop(:)
+type (physprop_type), pointer :: physprop(:) => null()
! Temporary storage location for filenames in namelist, and construction of dynamic index
! to properties. The unique filenames specified in the namelist are the identifiers of
@@ -223,6 +240,7 @@ subroutine physprop_init()
! nulls which aren't dealt with by trim()
integer :: ierr ! error codes from mpi
+ integer :: npernode ! ranks sharing a node (for shared-memory report)
!------------------------------------------------------------------------------------
@@ -291,10 +309,63 @@ subroutine physprop_init()
call pio_closefile(nc_id)
end do
+
+ ! report shared-memory savings (collective: every rank calls the helper so the
+ ! lazy node-communicator setup is not entered on masterproc alone)
+ npernode = cam_shmem_npes_per_node()
+ if (masterproc) then
+ write(iulog,*) 'phys_prop: large aerosol-optics tables held in per-node shared memory; ', &
+ npernode-1, ' redundant copies/node avoided'
+ end if
+
end subroutine physprop_init
!================================================================================================
+subroutine physprop_final
+
+ ! Release the per-node shared-memory optics tables (collective over the node
+ ! communicator via MPI_Win_free) and the physprop array spine. Called on every
+ ! rank from phys_final, after radiation's last use of the tables. Each
+ ! cam_shmem_free is a no-op for entries whose optics method did not share that
+ ! table (win == -1 on every rank), so the free sequence stays collective-safe.
+
+ integer :: i
+
+ if (.not. associated(physprop)) return
+
+ do i = 1, numphysprops
+ call cam_shmem_free(physprop(i)%extpsw, physprop(i)%win_extpsw)
+ call cam_shmem_free(physprop(i)%abspsw, physprop(i)%win_abspsw)
+ call cam_shmem_free(physprop(i)%asmpsw, physprop(i)%win_asmpsw)
+ call cam_shmem_free(physprop(i)%absplw, physprop(i)%win_absplw)
+ call cam_shmem_free(physprop(i)%sw_hygro_ext, physprop(i)%win_sw_hygro_ext)
+ call cam_shmem_free(physprop(i)%sw_hygro_ssa, physprop(i)%win_sw_hygro_ssa)
+ call cam_shmem_free(physprop(i)%sw_hygro_asm, physprop(i)%win_sw_hygro_asm)
+ call cam_shmem_free(physprop(i)%lw_hygro_abs, physprop(i)%win_lw_hygro_abs)
+ call cam_shmem_free(physprop(i)%sw_hygro_ext_wtp, physprop(i)%win_sw_hygro_ext_wtp)
+ call cam_shmem_free(physprop(i)%sw_hygro_ssa_wtp, physprop(i)%win_sw_hygro_ssa_wtp)
+ call cam_shmem_free(physprop(i)%sw_hygro_asm_wtp, physprop(i)%win_sw_hygro_asm_wtp)
+ call cam_shmem_free(physprop(i)%lw_hygro_abs_wtp, physprop(i)%win_lw_hygro_abs_wtp)
+ call cam_shmem_free(physprop(i)%sw_hygro_coreshell_ext, physprop(i)%win_cs_ext)
+ call cam_shmem_free(physprop(i)%sw_hygro_coreshell_ssa, physprop(i)%win_cs_ssa)
+ call cam_shmem_free(physprop(i)%sw_hygro_coreshell_asm, physprop(i)%win_cs_asm)
+ call cam_shmem_free(physprop(i)%lw_hygro_coreshell_abs, physprop(i)%win_cs_abs)
+ call cam_shmem_free(physprop(i)%r_sw_ext, physprop(i)%win_r_sw_ext)
+ call cam_shmem_free(physprop(i)%r_sw_scat, physprop(i)%win_r_sw_scat)
+ call cam_shmem_free(physprop(i)%r_sw_ascat, physprop(i)%win_r_sw_ascat)
+ call cam_shmem_free(physprop(i)%r_lw_abs, physprop(i)%win_r_lw_abs)
+ end do
+
+ deallocate(physprop)
+ nullify(physprop)
+ if (allocated(uniquefilenames)) deallocate(uniquefilenames)
+ numphysprops = 0
+
+end subroutine physprop_final
+
+!================================================================================================
+
integer function physprop_get_id(filename)
! Look for filename in the global list of unique filenames (module data uniquefilenames).
@@ -588,10 +659,14 @@ subroutine hygro_optics_init(phys_prop, nc_id)
real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
!------------------------------------------------------------------------------------
- allocate(phys_prop%sw_hygro_ext(nrh,nswbands))
- allocate(phys_prop%sw_hygro_ssa(nrh,nswbands))
- allocate(phys_prop%sw_hygro_asm(nrh,nswbands))
+ ! sw hygroscopic-growth tables held one copy per node in shared memory
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ext, phys_prop%win_sw_hygro_ext, nrh, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ssa, phys_prop%win_sw_hygro_ssa, nrh, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_asm, phys_prop%win_sw_hygro_asm, nrh, nswbands)
allocate(phys_prop%lw_abs(nlwbands))
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm)
ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_idx_id)
@@ -637,21 +712,28 @@ subroutine hygro_optics_init(phys_prop, nc_id)
ierr = pio_get_var(nc_id, lw_ext_id, phys_prop%lw_abs)
- ! interpolate onto cam's rh mesh
- do kbnd = 1,nswbands
- do krh = 1, nrh
- rh = 1.0_r8 / nrh * (krh - 1)
- phys_prop%sw_hygro_ext(krh,kbnd) = &
- exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
- * fsw_ext(1, kbnd)
- phys_prop%sw_hygro_ssa(krh,kbnd) = &
- lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
- * fsw_ssa(1, kbnd)
- phys_prop%sw_hygro_asm(krh,kbnd) = &
- lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
- * fsw_asm(1, kbnd)
+ ! interpolate onto cam's rh mesh (leader fills the node-shared copy)
+ if (cam_shmem_is_leader()) then
+ do kbnd = 1,nswbands
+ do krh = 1, nrh
+ rh = 1.0_r8 / nrh * (krh - 1)
+ phys_prop%sw_hygro_ext(krh,kbnd) = &
+ exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
+ * fsw_ext(1, kbnd)
+ phys_prop%sw_hygro_ssa(krh,kbnd) = &
+ lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
+ * fsw_ssa(1, kbnd)
+ phys_prop%sw_hygro_asm(krh,kbnd) = &
+ lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
+ * fsw_asm(1, kbnd)
+ enddo
enddo
- enddo
+ end if
+
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm)
deallocate (fsw_ext, fsw_asm, fsw_ssa, frh)
@@ -765,16 +847,22 @@ subroutine volcanic_radius_optics_init(phys_prop, nc_id)
integer :: sw_ext_id, sw_scat_id, sw_ascat_id, lw_abs_id
integer :: swbands, nbnd, n_mu_samples
integer :: ierr ! error flag
+ real(r8), allocatable :: tmp2(:,:) ! temp PIO read buffer (shared tables)
!------------------------------------------------------------------------------------
ierr = pio_inq_dimid(nc_id, 'mu_samples', mu_did)
ierr = pio_inq_dimlen(nc_id, mu_did, n_mu_samples)
- allocate (phys_prop%r_sw_ext(nswbands,n_mu_samples))
- allocate (phys_prop%r_sw_scat(nswbands,n_mu_samples))
- allocate (phys_prop%r_sw_ascat(nswbands,n_mu_samples))
- allocate (phys_prop%r_lw_abs(nlwbands,n_mu_samples))
+ ! radius-dependent tables held one copy per node in shared memory
+ call cam_shmem_alloc_r8_2d(phys_prop%r_sw_ext, phys_prop%win_r_sw_ext, nswbands, n_mu_samples)
+ call cam_shmem_alloc_r8_2d(phys_prop%r_sw_scat, phys_prop%win_r_sw_scat, nswbands, n_mu_samples)
+ call cam_shmem_alloc_r8_2d(phys_prop%r_sw_ascat, phys_prop%win_r_sw_ascat, nswbands, n_mu_samples)
+ call cam_shmem_alloc_r8_2d(phys_prop%r_lw_abs, phys_prop%win_r_lw_abs, nlwbands, n_mu_samples)
allocate (phys_prop%mu(n_mu_samples))
+ call cam_shmem_fence(phys_prop%win_r_sw_ext)
+ call cam_shmem_fence(phys_prop%win_r_sw_scat)
+ call cam_shmem_fence(phys_prop%win_r_sw_ascat)
+ call cam_shmem_fence(phys_prop%win_r_lw_abs)
ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
@@ -797,12 +885,27 @@ subroutine volcanic_radius_optics_init(phys_prop, nc_id)
ierr = pio_inq_varid(nc_id, 'babs_lw', lw_abs_id)
ierr = pio_inq_varid(nc_id, 'mu_samples', mu_id)
- ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%r_sw_ext)
- ierr = pio_get_var(nc_id, sw_scat_id, phys_prop%r_sw_scat)
- ierr = pio_get_var(nc_id, sw_ascat_id, phys_prop%r_sw_ascat)
- ierr = pio_get_var(nc_id, lw_abs_id, phys_prop%r_lw_abs)
+ ! all ranks read (PIO is collective); leader copies into the shared windows
+ allocate(tmp2(nswbands,n_mu_samples))
+ ierr = pio_get_var(nc_id, sw_ext_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%r_sw_ext = tmp2
+ ierr = pio_get_var(nc_id, sw_scat_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%r_sw_scat = tmp2
+ ierr = pio_get_var(nc_id, sw_ascat_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%r_sw_ascat = tmp2
+ deallocate(tmp2)
+ allocate(tmp2(nlwbands,n_mu_samples))
+ ierr = pio_get_var(nc_id, lw_abs_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%r_lw_abs = tmp2
+ deallocate(tmp2)
ierr = pio_get_var(nc_id, mu_id, phys_prop%mu)
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(phys_prop%win_r_sw_ext)
+ call cam_shmem_fence(phys_prop%win_r_sw_scat)
+ call cam_shmem_fence(phys_prop%win_r_sw_ascat)
+ call cam_shmem_fence(phys_prop%win_r_lw_abs)
+
! read bulk aero props
call bulk_props_init(phys_prop, nc_id)
@@ -888,10 +991,15 @@ subroutine hygroscopic_optics_init(phys_prop, nc_id)
character(len=*), parameter :: sub = 'hygroscopic_optics_init'
!------------------------------------------------------------------------------------
- allocate(phys_prop%sw_hygro_ext(nrh,nswbands))
- allocate(phys_prop%sw_hygro_ssa(nrh,nswbands))
- allocate(phys_prop%sw_hygro_asm(nrh,nswbands))
- allocate(phys_prop%lw_hygro_abs(nrh,nlwbands))
+ ! hygroscopic-growth tables held one copy per node in shared memory
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ext, phys_prop%win_sw_hygro_ext, nrh, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ssa, phys_prop%win_sw_hygro_ssa, nrh, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_asm, phys_prop%win_sw_hygro_asm, nrh, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%lw_hygro_abs, phys_prop%win_lw_hygro_abs, nrh, nlwbands)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm)
+ call cam_shmem_fence(phys_prop%win_lw_hygro_abs)
ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_idx_id)
ierr = pio_inq_dimlen(nc_id, rh_idx_id, nfilerh)
@@ -925,29 +1033,37 @@ subroutine hygroscopic_optics_init(phys_prop, nc_id)
ierr = pio_get_var(nc_id, sw_asm_id, fsw_asm)
ierr = pio_get_var(nc_id, lw_ext_id, flw_abs)
- ! interpolate onto cam's rh mesh
- do kbnd = 1,nswbands
- do krh = 1, nrh
- rh = 1.0_r8 / nrh * (krh - 1)
- phys_prop%sw_hygro_ext(krh,kbnd) = &
- exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
- * fsw_ext(1, kbnd)
- phys_prop%sw_hygro_ssa(krh,kbnd) = &
- lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
- * fsw_ssa(1, kbnd)
- phys_prop%sw_hygro_asm(krh,kbnd) = &
- lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
- * fsw_asm(1, kbnd)
+ ! interpolate onto cam's rh mesh (leader fills the node-shared copy)
+ if (cam_shmem_is_leader()) then
+ do kbnd = 1,nswbands
+ do krh = 1, nrh
+ rh = 1.0_r8 / nrh * (krh - 1)
+ phys_prop%sw_hygro_ext(krh,kbnd) = &
+ exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
+ * fsw_ext(1, kbnd)
+ phys_prop%sw_hygro_ssa(krh,kbnd) = &
+ lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
+ * fsw_ssa(1, kbnd)
+ phys_prop%sw_hygro_asm(krh,kbnd) = &
+ lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
+ * fsw_asm(1, kbnd)
+ enddo
enddo
- enddo
- do kbnd = 1,nlwbands
- do krh = 1, nrh
- rh = 1.0_r8 / nrh * (krh - 1)
- phys_prop%lw_hygro_abs(krh,kbnd) = &
- exp_interpol( frh, flw_abs(:,kbnd) / flw_abs(1,kbnd), rh ) &
- * flw_abs(1, kbnd)
+ do kbnd = 1,nlwbands
+ do krh = 1, nrh
+ rh = 1.0_r8 / nrh * (krh - 1)
+ phys_prop%lw_hygro_abs(krh,kbnd) = &
+ exp_interpol( frh, flw_abs(:,kbnd) / flw_abs(1,kbnd), rh ) &
+ * flw_abs(1, kbnd)
+ enddo
enddo
- enddo
+ end if
+
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm)
+ call cam_shmem_fence(phys_prop%win_lw_hygro_abs)
deallocate (fsw_ext, fsw_asm, fsw_ssa, flw_abs, frh)
@@ -1141,32 +1257,38 @@ subroutine modal_optics_init(props, ncid)
ierr = pio_inq_dimid(ncid, 'refindex_im', did)
ierr = pio_inq_dimlen(ncid, did, props%prefi)
- ! Allocate arrays
+ ! Allocate arrays. The four large specific-optics tables are held one copy
+ ! per node in shared memory; the small refractive-index tables stay per-rank.
+ call cam_shmem_alloc_r8_4d(props%extpsw, props%win_extpsw, props%ncoef,props%prefr,props%prefi,nswbands)
+ call cam_shmem_alloc_r8_4d(props%abspsw, props%win_abspsw, props%ncoef,props%prefr,props%prefi,nswbands)
+ call cam_shmem_alloc_r8_4d(props%asmpsw, props%win_asmpsw, props%ncoef,props%prefr,props%prefi,nswbands)
+ call cam_shmem_alloc_r8_4d(props%absplw, props%win_absplw, props%ncoef,props%prefr,props%prefi,nlwbands)
allocate( &
- props%extpsw(props%ncoef,props%prefr,props%prefi,nswbands), &
- props%abspsw(props%ncoef,props%prefr,props%prefi,nswbands), &
- props%asmpsw(props%ncoef,props%prefr,props%prefi,nswbands), &
- props%absplw(props%ncoef,props%prefr,props%prefi,nlwbands), &
props%refrtabsw(props%prefr,nswbands), &
props%refitabsw(props%prefi,nswbands), &
props%refrtablw(props%prefr,nlwbands), &
props%refitablw(props%prefi,nlwbands) )
+ ! open the window epoch; the node leader fills the shared tables below
+ call cam_shmem_fence(props%win_extpsw)
+ call cam_shmem_fence(props%win_abspsw)
+ call cam_shmem_fence(props%win_asmpsw)
+ call cam_shmem_fence(props%win_absplw)
! allocate temp to remove the mode dimension from the sw variables
allocate(rval(props%ncoef,props%prefr,props%prefi,1,nswbands))
ierr = pio_inq_varid(ncid, 'extpsw', vid)
ierr = pio_get_var(ncid, vid, rval)
- props%extpsw = rval(:,:,:,1,:)
+ if (cam_shmem_is_leader()) props%extpsw = rval(:,:,:,1,:)
ierr = pio_inq_varid(ncid, 'abspsw', vid)
ierr = pio_get_var(ncid, vid, rval)
- props%abspsw = rval(:,:,:,1,:)
+ if (cam_shmem_is_leader()) props%abspsw = rval(:,:,:,1,:)
ierr = pio_inq_varid(ncid, 'asmpsw', vid)
ierr = pio_get_var(ncid, vid, rval)
- props%asmpsw = rval(:,:,:,1,:)
+ if (cam_shmem_is_leader()) props%asmpsw = rval(:,:,:,1,:)
deallocate(rval)
@@ -1175,10 +1297,16 @@ subroutine modal_optics_init(props, ncid)
ierr = pio_inq_varid(ncid, 'absplw', vid)
ierr = pio_get_var(ncid, vid, rval)
- props%absplw = rval(:,:,:,1,:)
+ if (cam_shmem_is_leader()) props%absplw = rval(:,:,:,1,:)
deallocate(rval)
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(props%win_extpsw)
+ call cam_shmem_fence(props%win_abspsw)
+ call cam_shmem_fence(props%win_asmpsw)
+ call cam_shmem_fence(props%win_absplw)
+
ierr = pio_inq_varid(ncid, 'refindex_real_sw', vid)
ierr = pio_get_var(ncid, vid, props%refrtabsw)
@@ -1550,6 +1678,7 @@ subroutine hygroscopic_coreshell_optics_init(phys_prop, nc_id)
integer :: nrh ! number of rh values in file
integer :: nfrac ! number of core/shell ratio values in file
integer :: nbcdust,nkap
+ real(r8), allocatable :: tmp5(:,:,:,:,:) ! temp PIO read buffer (shared tables)
real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
character(len=*), parameter :: sub = 'hygroscopic_coreshell_optics_init'
@@ -1582,18 +1711,23 @@ subroutine hygroscopic_coreshell_optics_init(phys_prop, nc_id)
ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_id)
ierr = pio_inq_dimlen(nc_id, rh_id, phys_prop%nrelh)
- allocate(phys_prop%sw_hygro_coreshell_ext(phys_prop%nrelh,nswbands, &
- phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
- allocate(phys_prop%sw_hygro_coreshell_ssa(phys_prop%nrelh,nswbands, &
- phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
- allocate(phys_prop%sw_hygro_coreshell_asm(phys_prop%nrelh,nswbands, &
- phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
- allocate(phys_prop%lw_hygro_coreshell_abs(phys_prop%nrelh,nlwbands, &
- phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
+ ! large coreshell tables held one copy per node in shared memory
+ call cam_shmem_alloc_r8_5d(phys_prop%sw_hygro_coreshell_ext, phys_prop%win_cs_ext, &
+ phys_prop%nrelh,nswbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap)
+ call cam_shmem_alloc_r8_5d(phys_prop%sw_hygro_coreshell_ssa, phys_prop%win_cs_ssa, &
+ phys_prop%nrelh,nswbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap)
+ call cam_shmem_alloc_r8_5d(phys_prop%sw_hygro_coreshell_asm, phys_prop%win_cs_asm, &
+ phys_prop%nrelh,nswbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap)
+ call cam_shmem_alloc_r8_5d(phys_prop%lw_hygro_coreshell_abs, phys_prop%win_cs_abs, &
+ phys_prop%nrelh,nlwbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap)
allocate(phys_prop%corefrac(phys_prop%nfrac))
allocate(phys_prop%bcdust(phys_prop%nbcdust))
allocate(phys_prop%kap(phys_prop%nkap))
allocate(phys_prop%relh(phys_prop%nrelh))
+ call cam_shmem_fence(phys_prop%win_cs_ext)
+ call cam_shmem_fence(phys_prop%win_cs_ssa)
+ call cam_shmem_fence(phys_prop%win_cs_asm)
+ call cam_shmem_fence(phys_prop%win_cs_abs)
ierr = pio_inq_varid(nc_id, 'rh', rh_id)
ierr = pio_inq_varid(nc_id, 'coreshellratio', coreshell_id) ! modified by Pengfei for coreshell
@@ -1605,15 +1739,30 @@ subroutine hygroscopic_coreshell_optics_init(phys_prop, nc_id)
ierr = pio_inq_varid(nc_id, 'asm_sw_coreshell', sw_asm_id)
ierr = pio_inq_varid(nc_id, 'abs_lw_coreshell', lw_abs_id)
- ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%sw_hygro_coreshell_ext)
- ierr = pio_get_var(nc_id, sw_ssa_id, phys_prop%sw_hygro_coreshell_ssa)
- ierr = pio_get_var(nc_id, sw_asm_id, phys_prop%sw_hygro_coreshell_asm)
- ierr = pio_get_var(nc_id, lw_abs_id, phys_prop%lw_hygro_coreshell_abs)
+ ! all ranks read (PIO is collective); leader copies into the shared windows
+ allocate(tmp5(phys_prop%nrelh,nswbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
+ ierr = pio_get_var(nc_id, sw_ext_id, tmp5)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_coreshell_ext = tmp5
+ ierr = pio_get_var(nc_id, sw_ssa_id, tmp5)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_coreshell_ssa = tmp5
+ ierr = pio_get_var(nc_id, sw_asm_id, tmp5)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_coreshell_asm = tmp5
+ deallocate(tmp5)
+ allocate(tmp5(phys_prop%nrelh,nlwbands,phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
+ ierr = pio_get_var(nc_id, lw_abs_id, tmp5)
+ if (cam_shmem_is_leader()) phys_prop%lw_hygro_coreshell_abs = tmp5
+ deallocate(tmp5)
ierr = pio_get_var(nc_id, kap_id, phys_prop%kap)
ierr = pio_get_var(nc_id, rh_id, phys_prop%relh)
ierr = pio_get_var(nc_id, dstbc_id, phys_prop%bcdust)
ierr = pio_get_var(nc_id, coreshell_id, phys_prop%corefrac)
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(phys_prop%win_cs_ext)
+ call cam_shmem_fence(phys_prop%win_cs_ssa)
+ call cam_shmem_fence(phys_prop%win_cs_asm)
+ call cam_shmem_fence(phys_prop%win_cs_abs)
+
! read refractive index data if available
call refindex_aer_init(phys_prop, nc_id)
@@ -1634,6 +1783,7 @@ subroutine hygroscopic_wtp_optics_init(phys_prop, nc_id)
integer :: lw_band_id, sw_band_id, did
integer :: sw_ext_wtp_id, sw_ssa_wtp_id, sw_asm_wtp_id, lw_ext_wtp_id, wtp_id
integer :: nbnd, swbands
+ real(r8), allocatable :: tmp2(:,:) ! temp PIO read buffer (shared tables)
real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
character(len=*), parameter :: sub = 'hygroscopic_wtp_optics_init'
@@ -1645,11 +1795,16 @@ subroutine hygroscopic_wtp_optics_init(phys_prop, nc_id)
ierr = pio_inq_dimlen(nc_id, did, phys_prop%nwtp)
- allocate(phys_prop%sw_hygro_ext_wtp(phys_prop%nwtp,nswbands))
- allocate(phys_prop%sw_hygro_ssa_wtp(phys_prop%nwtp,nswbands))
- allocate(phys_prop%sw_hygro_asm_wtp(phys_prop%nwtp,nswbands))
- allocate(phys_prop%lw_hygro_abs_wtp(phys_prop%nwtp,nlwbands))
+ ! large weight-percent tables held one copy per node in shared memory
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ext_wtp, phys_prop%win_sw_hygro_ext_wtp, phys_prop%nwtp, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_ssa_wtp, phys_prop%win_sw_hygro_ssa_wtp, phys_prop%nwtp, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%sw_hygro_asm_wtp, phys_prop%win_sw_hygro_asm_wtp, phys_prop%nwtp, nswbands)
+ call cam_shmem_alloc_r8_2d(phys_prop%lw_hygro_abs_wtp, phys_prop%win_lw_hygro_abs_wtp, phys_prop%nwtp, nlwbands)
allocate(phys_prop%wgtpct(phys_prop%nwtp))
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext_wtp)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa_wtp)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm_wtp)
+ call cam_shmem_fence(phys_prop%win_lw_hygro_abs_wtp)
ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
@@ -1667,12 +1822,27 @@ subroutine hygroscopic_wtp_optics_init(phys_prop, nc_id)
ierr = pio_inq_varid(nc_id, 'abs_lw_wtp', lw_ext_wtp_id)
ierr = pio_inq_varid(nc_id, 'wgtpct', wtp_id)
- ierr = pio_get_var(nc_id, sw_ext_wtp_id, phys_prop%sw_hygro_ext_wtp)
- ierr = pio_get_var(nc_id, sw_ssa_wtp_id, phys_prop%sw_hygro_ssa_wtp)
- ierr = pio_get_var(nc_id, sw_asm_wtp_id, phys_prop%sw_hygro_asm_wtp)
- ierr = pio_get_var(nc_id, lw_ext_wtp_id, phys_prop%lw_hygro_abs_wtp)
+ ! all ranks read (PIO is collective); leader copies into the shared windows
+ allocate(tmp2(phys_prop%nwtp,nswbands))
+ ierr = pio_get_var(nc_id, sw_ext_wtp_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_ext_wtp = tmp2
+ ierr = pio_get_var(nc_id, sw_ssa_wtp_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_ssa_wtp = tmp2
+ ierr = pio_get_var(nc_id, sw_asm_wtp_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%sw_hygro_asm_wtp = tmp2
+ deallocate(tmp2)
+ allocate(tmp2(phys_prop%nwtp,nlwbands))
+ ierr = pio_get_var(nc_id, lw_ext_wtp_id, tmp2)
+ if (cam_shmem_is_leader()) phys_prop%lw_hygro_abs_wtp = tmp2
+ deallocate(tmp2)
ierr = pio_get_var(nc_id, wtp_id, phys_prop%wgtpct)
+ ! publish the filled tables to every rank on the node
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ext_wtp)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_ssa_wtp)
+ call cam_shmem_fence(phys_prop%win_sw_hygro_asm_wtp)
+ call cam_shmem_fence(phys_prop%win_lw_hygro_abs_wtp)
+
! read refractive index data if available
call refindex_aer_init(phys_prop, nc_id)
diff --git a/src/physics/cam/physpkg.F90 b/src/physics/cam/physpkg.F90
index 28a9714c9c..7a51cf61f1 100644
--- a/src/physics/cam/physpkg.F90
+++ b/src/physics/cam/physpkg.F90
@@ -1333,6 +1333,7 @@ subroutine phys_final( phys_state, phys_tend, pbuf2d )
use phys_grid_ctem, only : phys_grid_ctem_final
use nudging, only: Nudge_Model, nudging_final
use ctem_diags_mod, only: ctem_diags_final
+ use phys_prop, only : physprop_final
!-----------------------------------------------------------------------
!
@@ -1352,6 +1353,7 @@ subroutine phys_final( phys_state, phys_tend, pbuf2d )
deallocate(phys_state)
deallocate(phys_tend)
call chem_final
+ call physprop_final
call carma_final
call wv_sat_final
call microp_aero_final()
diff --git a/src/physics/cam7/physpkg.F90 b/src/physics/cam7/physpkg.F90
index d00643c8d0..be6ab67f94 100644
--- a/src/physics/cam7/physpkg.F90
+++ b/src/physics/cam7/physpkg.F90
@@ -1330,6 +1330,7 @@ subroutine phys_final( phys_state, phys_tend, pbuf2d )
use nudging, only: Nudge_Model, nudging_final
use hemco_interface, only: HCOI_Chunk_Final
use ctem_diags_mod, only: ctem_diags_final
+ use phys_prop, only: physprop_final
!-----------------------------------------------------------------------
!
@@ -1349,6 +1350,7 @@ subroutine phys_final( phys_state, phys_tend, pbuf2d )
deallocate(phys_state)
deallocate(phys_tend)
call chem_final
+ call physprop_final
call carma_final
call wv_sat_final
call microp_aero_final()
diff --git a/src/utils/cam_shmem_mod.F90 b/src/utils/cam_shmem_mod.F90
new file mode 100644
index 0000000000..609e954cb1
--- /dev/null
+++ b/src/utils/cam_shmem_mod.F90
@@ -0,0 +1,452 @@
+module cam_shmem_mod
+!-------------------------------------------------------------------------------
+! Per-node MPI-3 shared-memory allocation helper for large, read-only lookup
+! tables that are identical on every MPI rank (e.g. photolysis cross-section /
+! radiative-source-function tables, aerosol-optics property tables). Instead of
+! every rank holding its own copy, one physical copy is allocated per
+! shared-memory node and mapped (read-only) into every rank on that node. This
+! frees (ranks_per_node - 1) copies per node.
+!
+! Usage pattern (collective over the global communicator):
+! 1. call cam_shmem_alloc_rX_Nd(ptr, win, dims...) ! all ranks
+! 2. call cam_shmem_fence(win) ! open epoch
+! 3. node leaders fill the table (read from file / broadcast among leaders)
+! 4. call cam_shmem_fence(win) ! publish; all ranks may now read
+! After step 4 the data is static and may be read for the rest of the run with no
+! further synchronization. Free it with cam_shmem_free(ptr, win) at finalize
+! (collective over the node communicator); otherwise MPI_Finalize reclaims it.
+!
+! This module uses the F90 'mpi' module (not mpif.h). MPI-3.0 makes the TYPE(C_PTR)
+! overloads of MPI_WIN_ALLOCATE_SHARED / MPI_WIN_SHARED_QUERY optional in the 'mpi'
+! module, and some implementations (cray-mpich) do not expose those two routines
+! there at all, so explicit interfaces for them are declared below instead.
+!-------------------------------------------------------------------------------
+
+ use shr_kind_mod, only: r4 => shr_kind_r4, r8 => shr_kind_r8
+ use cam_abortutils, only: endrun
+
+#ifdef SPMD
+ use spmd_utils, only: mpicom
+ use mpi, only: MPI_ADDRESS_KIND, MPI_COMM_NULL, MPI_COMM_TYPE_SHARED, &
+ MPI_INFO_NULL, MPI_SUCCESS, MPI_UNDEFINED, MPI_WIN_NULL, &
+ mpi_comm_rank, mpi_comm_size, mpi_comm_split, &
+ mpi_comm_split_type, mpi_win_fence, mpi_win_free
+ use, intrinsic :: iso_c_binding, only: c_ptr, c_f_pointer
+#endif
+
+ implicit none
+ private
+
+ public :: cam_shmem_alloc_r4_4d ! allocate node-shared real(r4) rank-4 table
+ public :: cam_shmem_alloc_r4_5d ! allocate node-shared real(r4) rank-5 table
+ public :: cam_shmem_alloc_r8_2d ! allocate node-shared real(r8) rank-2 table
+ public :: cam_shmem_alloc_r8_3d ! allocate node-shared real(r8) rank-3 table
+ public :: cam_shmem_alloc_r8_4d ! allocate node-shared real(r8) rank-4 table
+ public :: cam_shmem_alloc_r8_5d ! allocate node-shared real(r8) rank-5 table
+ public :: cam_shmem_free ! free a node-shared table (MPI_Win_free)
+ public :: cam_shmem_fence ! synchronize a window (publish writes)
+ public :: cam_shmem_is_leader ! .true. on the leader (rank 0) of this node
+ public :: cam_shmem_leader_comm ! communicator containing only node leaders
+ public :: cam_shmem_npes_per_node ! number of ranks sharing this node
+ public :: cam_shmem_init ! force one-time node-comm setup at an early collective point
+
+ ! Generic finalizer: free whatever node-shared table the pointer aliases. Safe
+ ! to call blindly on an unallocated table (win == -1 / null pointer -> no-op).
+ interface cam_shmem_free
+ module procedure cam_shmem_free_r4_4d
+ module procedure cam_shmem_free_r4_5d
+ module procedure cam_shmem_free_r8_2d
+ module procedure cam_shmem_free_r8_3d
+ module procedure cam_shmem_free_r8_4d
+ module procedure cam_shmem_free_r8_5d
+ end interface cam_shmem_free
+
+ integer, parameter :: bytes_r4 = 4
+ integer, parameter :: bytes_r8 = 8
+
+#ifdef SPMD
+ ! MPI_WIN_ALLOCATE_SHARED and MPI_WIN_SHARED_QUERY take a TYPE(C_PTR) baseptr in
+ ! MPI-3, but that overload is not exposed by every implementation's F90 'mpi'
+ ! module (cray-mpich does not expose these two routines at all), so they cannot be
+ ! named in the 'only' list above. Declare explicit interfaces here rather than
+ ! relying on implicit ones, so the argument types are still checked. Both forms
+ ! resolve to the same mpi_win_allocate_shared_ / mpi_win_shared_query_ bindings,
+ ! so this is equivalent on implementations that do expose them.
+ interface
+ subroutine mpi_win_allocate_shared(nbytes, disp_unit, info, comm, baseptr, win, ierror)
+ import :: MPI_ADDRESS_KIND, c_ptr
+ integer(kind=MPI_ADDRESS_KIND), intent(in) :: nbytes
+ integer, intent(in) :: disp_unit
+ integer, intent(in) :: info
+ integer, intent(in) :: comm
+ type(c_ptr), intent(out) :: baseptr
+ integer, intent(out) :: win
+ integer, intent(out) :: ierror
+ end subroutine mpi_win_allocate_shared
+
+ subroutine mpi_win_shared_query(win, rank, segsize, disp_unit, baseptr, ierror)
+ import :: MPI_ADDRESS_KIND, c_ptr
+ integer, intent(in) :: win
+ integer, intent(in) :: rank
+ integer(kind=MPI_ADDRESS_KIND), intent(out) :: segsize
+ integer, intent(out) :: disp_unit
+ type(c_ptr), intent(out) :: baseptr
+ integer, intent(out) :: ierror
+ end subroutine mpi_win_shared_query
+ end interface
+
+ logical, save :: initialized = .false.
+ integer, save :: node_comm = MPI_COMM_NULL ! ranks sharing a node
+ integer, save :: leader_comm = MPI_COMM_NULL ! one rank per node (the leaders)
+ integer, save :: node_rank = 0
+ integer, save :: node_size = 1
+ logical, save :: is_leader = .true.
+#else
+ ! Non-SPMD build: single task, nothing is shared.
+ integer, save :: node_size = 1
+ logical, save :: is_leader = .true.
+#endif
+
+contains
+
+!===============================================================================
+
+#ifdef SPMD
+ subroutine init_comms()
+ ! Lazily build the node-local and node-leader communicators. Collective
+ ! over mpicom; safe to call from every shared-memory allocation request.
+ integer :: ierr, color
+
+ if (initialized) return
+
+ call mpi_comm_split_type(mpicom, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, &
+ node_comm, ierr)
+ call mpi_comm_rank(node_comm, node_rank, ierr)
+ call mpi_comm_size(node_comm, node_size, ierr)
+ is_leader = (node_rank == 0)
+
+ ! Communicator of node leaders only (used to distribute data between nodes).
+ ! masterproc (global rank 0) is a node leader and is rank 0 here.
+ if (is_leader) then
+ color = 0
+ else
+ color = MPI_UNDEFINED
+ end if
+ call mpi_comm_split(mpicom, color, 0, leader_comm, ierr)
+
+ initialized = .true.
+ end subroutine init_comms
+
+!===============================================================================
+
+ subroutine shmem_alloc_bytes(nbytes, disp_unit, win, baseptr)
+ ! Allocate a node-shared window of nbytes (only the node leader requests
+ ! storage; peers map the leader's segment) and return the base C pointer of
+ ! the leader's contiguous segment. Type/shape-agnostic core shared by all
+ ! of the typed cam_shmem_alloc_* wrappers.
+ integer(kind=MPI_ADDRESS_KIND), intent(in) :: nbytes
+ integer, intent(in) :: disp_unit
+ integer, intent(out) :: win
+ type(c_ptr), intent(out) :: baseptr
+
+ integer(kind=MPI_ADDRESS_KIND) :: winsize, qsize
+ integer :: ierr, qdisp
+
+ call init_comms()
+
+ if (is_leader) then
+ winsize = nbytes
+ else
+ winsize = 0_MPI_ADDRESS_KIND
+ end if
+
+ call mpi_win_allocate_shared(winsize, disp_unit, MPI_INFO_NULL, node_comm, &
+ baseptr, win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_allocate_shared failed')
+
+ ! Non-leaders learn the address of the leader's (rank 0) contiguous segment.
+ if (.not. is_leader) then
+ call mpi_win_shared_query(win, 0, qsize, qdisp, baseptr, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_shared_query failed')
+ end if
+ end subroutine shmem_alloc_bytes
+#endif
+
+!===============================================================================
+
+ subroutine cam_shmem_init()
+ ! Force the one-time node-local / node-leader communicator setup NOW, at a
+ ! controlled early all-ranks collective point, instead of lazily on the first
+ ! cam_shmem_alloc_* during physics init. At large rank counts the first
+ ! MPI_Comm_split_type(MPI_COMM_TYPE_SHARED) can be expensive; doing it early
+ ! (rather than mid-physprop) keeps it predictable. Idempotent; non-SPMD no-op.
+#ifdef SPMD
+ call init_comms()
+#endif
+ end subroutine cam_shmem_init
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r4_4d(ptr, win, n1, n2, n3, n4)
+ real(r4), pointer, intent(out) :: ptr(:,:,:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2, n3, n4
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(n3,MPI_ADDRESS_KIND)*int(n4,MPI_ADDRESS_KIND) &
+ *int(bytes_r4,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r4, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2, n3, n4])
+#else
+ win = -1
+ allocate(ptr(n1,n2,n3,n4))
+#endif
+ end subroutine cam_shmem_alloc_r4_4d
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r4_5d(ptr, win, n1, n2, n3, n4, n5)
+ real(r4), pointer, intent(out) :: ptr(:,:,:,:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2, n3, n4, n5
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(n3,MPI_ADDRESS_KIND)*int(n4,MPI_ADDRESS_KIND) &
+ *int(n5,MPI_ADDRESS_KIND)*int(bytes_r4,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r4, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2, n3, n4, n5])
+#else
+ win = -1
+ allocate(ptr(n1,n2,n3,n4,n5))
+#endif
+ end subroutine cam_shmem_alloc_r4_5d
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r8_2d(ptr, win, n1, n2)
+ real(r8), pointer, intent(out) :: ptr(:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(bytes_r8,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r8, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2])
+#else
+ win = -1
+ allocate(ptr(n1,n2))
+#endif
+ end subroutine cam_shmem_alloc_r8_2d
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r8_3d(ptr, win, n1, n2, n3)
+ real(r8), pointer, intent(out) :: ptr(:,:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2, n3
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(n3,MPI_ADDRESS_KIND)*int(bytes_r8,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r8, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2, n3])
+#else
+ win = -1
+ allocate(ptr(n1,n2,n3))
+#endif
+ end subroutine cam_shmem_alloc_r8_3d
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r8_4d(ptr, win, n1, n2, n3, n4)
+ real(r8), pointer, intent(out) :: ptr(:,:,:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2, n3, n4
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(n3,MPI_ADDRESS_KIND)*int(n4,MPI_ADDRESS_KIND) &
+ *int(bytes_r8,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r8, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2, n3, n4])
+#else
+ win = -1
+ allocate(ptr(n1,n2,n3,n4))
+#endif
+ end subroutine cam_shmem_alloc_r8_4d
+
+!===============================================================================
+
+ subroutine cam_shmem_alloc_r8_5d(ptr, win, n1, n2, n3, n4, n5)
+ real(r8), pointer, intent(out) :: ptr(:,:,:,:,:)
+ integer, intent(out) :: win
+ integer, intent(in) :: n1, n2, n3, n4, n5
+#ifdef SPMD
+ integer(kind=MPI_ADDRESS_KIND) :: nbytes
+ type(c_ptr) :: baseptr
+ nbytes = int(n1,MPI_ADDRESS_KIND)*int(n2,MPI_ADDRESS_KIND) &
+ *int(n3,MPI_ADDRESS_KIND)*int(n4,MPI_ADDRESS_KIND) &
+ *int(n5,MPI_ADDRESS_KIND)*int(bytes_r8,MPI_ADDRESS_KIND)
+ call shmem_alloc_bytes(nbytes, bytes_r8, win, baseptr)
+ call c_f_pointer(baseptr, ptr, [n1, n2, n3, n4, n5])
+#else
+ win = -1
+ allocate(ptr(n1,n2,n3,n4,n5))
+#endif
+ end subroutine cam_shmem_alloc_r8_5d
+
+!===============================================================================
+
+ subroutine cam_shmem_fence(win)
+ ! Collective over the node communicator; synchronizes the window so that
+ ! stores by the leader become visible to all ranks on the node.
+ integer, intent(in) :: win
+#ifdef SPMD
+ integer :: ierr
+ call mpi_win_fence(0, win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_fence failed')
+#endif
+ end subroutine cam_shmem_fence
+
+!===============================================================================
+! cam_shmem_free generic procedures. Each frees the node-shared window (in SPMD)
+! and disassociates the pointer, or deallocates it (non-SPMD). Collective over
+! the node communicator in SPMD; a no-op when win == -1 (table never shared).
+!===============================================================================
+
+ subroutine cam_shmem_free_r4_4d(ptr, win)
+ real(r4), pointer :: ptr(:,:,:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r4_4d
+
+ subroutine cam_shmem_free_r4_5d(ptr, win)
+ real(r4), pointer :: ptr(:,:,:,:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r4_5d
+
+ subroutine cam_shmem_free_r8_2d(ptr, win)
+ real(r8), pointer :: ptr(:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r8_2d
+
+ subroutine cam_shmem_free_r8_3d(ptr, win)
+ real(r8), pointer :: ptr(:,:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r8_3d
+
+ subroutine cam_shmem_free_r8_4d(ptr, win)
+ real(r8), pointer :: ptr(:,:,:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r8_4d
+
+ subroutine cam_shmem_free_r8_5d(ptr, win)
+ real(r8), pointer :: ptr(:,:,:,:,:)
+ integer, intent(inout) :: win
+#ifdef SPMD
+ integer :: ierr
+ if (win /= -1 .and. win /= MPI_WIN_NULL) then
+ call mpi_win_free(win, ierr)
+ if (ierr /= MPI_SUCCESS) call endrun('cam_shmem_mod: MPI_Win_free failed')
+ end if
+ if (associated(ptr)) nullify(ptr)
+#else
+ if (associated(ptr)) deallocate(ptr)
+#endif
+ win = -1
+ end subroutine cam_shmem_free_r8_5d
+
+!===============================================================================
+
+ logical function cam_shmem_is_leader()
+#ifdef SPMD
+ call init_comms()
+#endif
+ cam_shmem_is_leader = is_leader
+ end function cam_shmem_is_leader
+
+!===============================================================================
+
+ integer function cam_shmem_leader_comm()
+#ifdef SPMD
+ call init_comms()
+ cam_shmem_leader_comm = leader_comm
+#else
+ cam_shmem_leader_comm = -1
+#endif
+ end function cam_shmem_leader_comm
+
+!===============================================================================
+
+ integer function cam_shmem_npes_per_node()
+#ifdef SPMD
+ call init_comms()
+#endif
+ cam_shmem_npes_per_node = node_size
+ end function cam_shmem_npes_per_node
+
+!===============================================================================
+
+end module cam_shmem_mod