24 module ensemble_manager_mod
30 use mpp_mod,
only : input_nml_file
31 use fms2_io_mod,
only : fms2_io_set_filename_appendix=>set_filename_appendix
37 integer,
parameter :: MAX_ENSEMBLE_SIZE = 100
40 integer,
allocatable,
dimension(:,:) :: ensemble_pelist
41 integer,
allocatable,
dimension(:,:) :: ensemble_pelist_ocean
42 integer,
allocatable,
dimension(:,:) :: ensemble_pelist_atmos
43 integer,
allocatable,
dimension(:,:) :: ensemble_pelist_land
44 integer,
allocatable,
dimension(:,:) :: ensemble_pelist_ice
45 integer,
allocatable,
dimension(:) :: ensemble_pelist_ocean_filter
46 integer,
allocatable,
dimension(:) :: ensemble_pelist_atmos_filter
47 integer,
allocatable,
dimension(:) :: ensemble_pelist_land_filter
48 integer,
allocatable,
dimension(:) :: ensemble_pelist_ice_filter
50 integer :: ensemble_size = 1
51 integer :: ensemble_id = 1
52 integer :: pe, total_npes_pm=0,ocean_npes_pm=0,atmos_npes_pm=0
53 integer :: land_npes_pm=0,ice_npes_pm=0
76 integer :: i, io_status, npes, ierr
78 namelist /ensemble_nml/ ensemble_size
80 read (input_nml_file, ensemble_nml, iostat=io_status)
83 if(ensemble_size < 1)
call mpp_error(fatal, &
84 'ensemble_manager_mod: ensemble_nml variable ensemble_size must be a positive integer')
85 if(ensemble_size > max_ensemble_size)
call mpp_error(fatal, &
86 'ensemble_manager_mod: ensemble_nml variable ensemble_size should be no larger than MAX_ENSEMBLE_SIZE, '// &
87 'change ensemble_size or increase MAX_ENSEMBLE_SIZE')
91 if (npes < ensemble_size)
then
92 call mpp_error(fatal,
'npes must be >= ensemble_size')
94 total_npes_pm = npes/ensemble_size
95 if (mod(npes, total_npes_pm) /= 0)
call mpp_error(fatal,
'ensemble_size must be divis by npes')
126 integer,
intent(inout) :: pelist(:,:)
127 character(len=*),
intent(in),
optional :: name
129 if (
size(pelist,1) < ensemble_size) &
130 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 1st index < ensemble_size')
132 if(
present(name))
then
135 if (
size(pelist,2) < ocean_npes_pm)&
136 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 2nd index < ocean_npes_pm')
138 pelist(1:ensemble_size,1:ocean_npes_pm) = &
139 ensemble_pelist_ocean(1:ensemble_size,1:ocean_npes_pm)
142 if (
size(pelist,2) < atmos_npes_pm)&
143 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 2nd index < atmos_npes_pm')
145 pelist(1:ensemble_size,1:atmos_npes_pm) = &
146 ensemble_pelist_atmos(1:ensemble_size,1:atmos_npes_pm)
149 if (
size(pelist,2) < land_npes_pm)&
150 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 2nd index < land_npes_pm')
152 pelist(1:ensemble_size,1:land_npes_pm) = &
153 ensemble_pelist_land(1:ensemble_size,1:land_npes_pm)
156 if (
size(pelist,2) < ice_npes_pm)&
157 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 2nd index < ice_npes_pm')
159 pelist(1:ensemble_size,1:ice_npes_pm) = &
160 ensemble_pelist_ice(1:ensemble_size,1:ice_npes_pm)
163 call mpp_error(fatal,
'get_ensemble_pelist: unknown argument name='//name)
166 if (
size(pelist,2) < total_npes_pm)&
167 call mpp_error(fatal,
'get_ensemble_pelist: size of pelist 2nd index < total_npes_pm')
169 pelist(1:ensemble_size,1:total_npes_pm) = &
170 ensemble_pelist(1:ensemble_size,1:total_npes_pm)
185 integer,
intent(inout) :: pelist(:)
186 character(len=*),
intent(in) :: name
190 if (
size(pelist) < ensemble_size * ocean_npes_pm)&
191 call mpp_error(fatal,
'get_ensemble_filter_pelist: size of pelist argument < ensemble_size * ocean_npes_pm')
193 pelist(1:ensemble_size*ocean_npes_pm) = &
194 ensemble_pelist_ocean_filter(1:ensemble_size*ocean_npes_pm)
197 if (
size(pelist) < ensemble_size * atmos_npes_pm)&
198 call mpp_error(fatal,
'get_ensemble_filter_pelist: size of pelist argument < ensemble_size * atmos_npes_pm')
200 pelist(1:ensemble_size*atmos_npes_pm) = &
201 ensemble_pelist_atmos_filter(1:ensemble_size*atmos_npes_pm)
204 if (
size(pelist) < ensemble_size * land_npes_pm)&
205 call mpp_error(fatal,
'get_ensemble_filter_pelist: size of pelist argument < ensemble_size * land_npes_pm')
207 pelist(1:ensemble_size*land_npes_pm) = &
208 ensemble_pelist_land_filter(1:ensemble_size*land_npes_pm)
211 if (
size(pelist) < ensemble_size * ice_npes_pm)&
212 call mpp_error(fatal,
'get_ensemble_filter_pelist: size of pelist argument < ensemble_size * ice_npes_pm')
214 pelist(1:ensemble_size*ice_npes_pm) = &
215 ensemble_pelist_ice_filter(1:ensemble_size*ice_npes_pm)
218 call mpp_error(fatal,
'get_ensemble_filter_pelist: unknown argument name='//name)
232 Atm_pelist, Ocean_pelist, Land_pelist, Ice_pelist)
233 logical,
intent(in) :: concurrent
234 integer,
intent(in) :: atmos_npes, ocean_npes
235 integer,
intent(in) :: land_npes, ice_npes
236 integer,
dimension(:),
intent(inout) :: atm_pelist, ocean_pelist
237 integer,
dimension(:),
intent(inout) :: land_pelist, ice_pelist
238 integer :: atmos_pe_start, atmos_pe_end, ocean_pe_start, ocean_pe_end
239 integer :: land_pe_start, land_pe_end, ice_pe_start, ice_pe_end
240 character(len=10) :: pelist_name, text
241 integer :: npes, n, m ,i
246 if(land_npes > atmos_npes)
call mpp_error(fatal,
'ensemble_manager_mod: land_npes > atmos_npes')
247 if(ice_npes > atmos_npes)
call mpp_error(fatal,
'ensemble_manager_mod: ice_npes > atmos_npes')
249 allocate(ensemble_pelist(ensemble_size,total_npes_pm))
250 allocate(ensemble_pelist_ocean(1:ensemble_size, 1:ocean_npes) )
251 allocate(ensemble_pelist_atmos(1:ensemble_size, 1:atmos_npes) )
252 allocate(ensemble_pelist_land(1:ensemble_size, 1:land_npes) )
253 allocate(ensemble_pelist_ice(1:ensemble_size, 1:ice_npes) )
259 if( concurrent .OR. atmos_npes+ocean_npes == npes )
then
260 ocean_pe_start = ensemble_size*atmos_npes
263 atmos_pe_end = atmos_pe_start + atmos_npes - 1
264 ocean_pe_end = ocean_pe_start + ocean_npes - 1
265 land_pe_end = land_pe_start + land_npes - 1
266 ice_pe_end = ice_pe_start + ice_npes - 1
267 ensemble_pelist_atmos(n, 1:atmos_npes) = (/(i,i=atmos_pe_start,atmos_pe_end)/)
268 ensemble_pelist_ocean(n, 1:ocean_npes) = (/(i,i=ocean_pe_start,ocean_pe_end)/)
269 ensemble_pelist_land(n, 1:land_npes) = (/(i,i=land_pe_start, land_pe_end)/)
270 ensemble_pelist_ice(n, 1:ice_npes) = (/(i,i=ice_pe_start, ice_pe_end)/)
271 ensemble_pelist(n, 1:atmos_npes) = ensemble_pelist_atmos(n, 1:atmos_npes)
272 if( concurrent .OR. atmos_npes+ocean_npes == npes ) &
273 ensemble_pelist(n, atmos_npes+1:npes) = ensemble_pelist_ocean(n, 1:ocean_npes)
274 if(any(ensemble_pelist(n,:) == pe)) ensemble_id = n
275 write(pelist_name,
'(a,i2.2)')
'_ens',n
277 atmos_pe_start = atmos_pe_end + 1
278 ocean_pe_start = ocean_pe_end + 1
279 land_pe_start = atmos_pe_start
280 ice_pe_start = atmos_pe_start
283 atm_pelist(:) = ensemble_pelist_atmos(ensemble_id,:)
284 ocean_pelist(:) = ensemble_pelist_ocean(ensemble_id,:)
285 land_pelist(:) = ensemble_pelist_land(ensemble_id,:)
286 ice_pelist(:) = ensemble_pelist_ice(ensemble_id,:)
337 write(pelist_name,
'(a,i2.2)')
'atm_ens',n
339 write(pelist_name,
'(a,i2.2)')
'ocn_ens',n
341 write(pelist_name,
'(a,i2.2)')
'lnd_ens',n
343 write(pelist_name,
'(a,i2.2)')
'ice_ens',n
347 write(pelist_name,
'(a,i2.2)')
'atm_ens',ensemble_id
349 write(pelist_name,
'(a,i2.2)')
'ocn_ens',ensemble_id
351 write(pelist_name,
'(a,i2.2)')
'lnd_ens',ensemble_id
353 write(pelist_name,
'(a,i2.2)')
'ice_ens',ensemble_id
357 ocean_npes_pm = ocean_npes
358 atmos_npes_pm = atmos_npes
359 land_npes_pm = land_npes
360 ice_npes_pm = ice_npes
363 allocate(ensemble_pelist_ocean_filter(ensemble_size * ocean_npes_pm))
364 allocate(ensemble_pelist_atmos_filter(ensemble_size * atmos_npes_pm))
365 allocate(ensemble_pelist_land_filter(ensemble_size * land_npes_pm))
366 allocate(ensemble_pelist_ice_filter(ensemble_size * ice_npes_pm))
369 i=(n-1)*ocean_npes_pm + m
370 ensemble_pelist_ocean_filter(i) = ensemble_pelist_ocean(n,m)
373 i=(n-1)*atmos_npes_pm + m
374 ensemble_pelist_atmos_filter(i) = ensemble_pelist_atmos(n,m)
377 i=(n-1)*land_npes_pm + m
378 ensemble_pelist_land_filter(i) = ensemble_pelist_land(n,m)
381 i=(n-1)*ice_npes_pm + m
382 ensemble_pelist_ice_filter(i) = ensemble_pelist_ice(n,m)
386 write(pelist_name,
'(a)')
'ocn_filter'
389 write(pelist_name,
'(a)')
'atm_filter'
392 write(pelist_name,
'(a)')
'lnd_filter'
395 write(pelist_name,
'(a)')
'ice_filter'
403 if (ensemble_size > 1)
then
404 write( text,
'(a,i2.2)' )
'ens_', ensemble_id
406 call fms2_io_set_filename_appendix(trim(text))
412 end module ensemble_manager_mod
integer function, public get_ensemble_id()
Getter function for ensemble_id.
integer function, dimension(6), public get_ensemble_size()
Returns ensemble size integer array.
subroutine, public ensemble_pelist_setup(concurrent, atmos_npes, ocean_npes, land_npes, ice_npes, Atm_pelist, Ocean_pelist, Land_pelist, Ice_pelist)
Sets up pe list for an ensemble.
subroutine, public ensemble_manager_init()
Initializes ensemble_manager_mod.
subroutine, public get_ensemble_pelist(pelist, name)
Gets pe list for current ensemble or a given ensemble component.
subroutine, public get_ensemble_filter_pelist(pelist, name)
Gets filter pelist for a given ensemble component.
integer function, public check_nml_error(IOSTAT, NML_NAME)
Checks the iostat argument that is returned after reading a namelist and determines if the error code...
integer function stdout()
This function returns the current standard fortran unit numbers for output.
subroutine mpp_declare_pelist(pelist, name, commID)
Declare a pelist.
integer function stdlog()
This function returns the current standard fortran unit numbers for log messages. Log messages,...
integer function mpp_npes()
Returns processor count for current pelist.
integer function mpp_pe()
Returns processor ID.