diff --git a/.gitignore b/.gitignore index f306e88a3c..2db1c24b31 100644 --- a/.gitignore +++ b/.gitignore @@ -203,6 +203,7 @@ test_quad_reg_interp test_table_read test_ran_unif test_kde_dist +test_window # Directories to NOT IGNORE ... same as executable names # as far as I know, these must be listed after the executables diff --git a/CHANGELOG.rst b/CHANGELOG.rst index b466545457..02bc0ddbc6 100644 --- a/CHANGELOG.rst +++ b/CHANGELOG.rst @@ -22,17 +22,43 @@ individual files. The changes are now listed with the most recent at the top. -**October 22 2024 :: Bug-fixes: WRF and GOES. Tag 11.8.2** +**November 12 2024 :: MPAS bug-fixes. Tag v11.8.5** + +- Fixed 2m and 10m fields not being updated - set istatus for VERTISHEIGHT + and VERTISLEVEL for convert_vert_distrib_state. +- Fixed vertical location in convert_vertical_obs to use zGridFace. + +**November 8 2024 :: POP initial ensemble available from GDEX. Tag v11.8.4** + +Documentation update: + + - POP initial ensemble available from GDEX. + +Bug-fixes: + + - Removed unnecessary loops around calendar types in time_manager_mod. + - Removed unused routine from normal_distribution_mod which case giving compilation warnings. + - Replaced broadcast_minmax calls with all_reduce_min_max. + + +**November 7 2024 :: MPI window memory reduction. Tag v11.8.3** + +- Removes unnecessary copy of state into mpi window. +- Removes cray pointer version of the mpi window. +- | Fortran-testanything included in developer tests. + | *From dennisdjensen: see developer_tests/contrib/fortran-testanything/LICENSE.txt* + +**October 22 2024 :: Bug-fixes: WRF and GOES. Tag v11.8.2** - Force THM to be the WRF-DART temperature variable - Remove offset on GOES observation converter -**September 27 2024 :: MOM6 mask bug-fix. Tag 11.8.1** +**September 27 2024 :: MOM6 mask bug-fix. Tag v11.8.1** - Fix for MOM6 CESM3 workhorse 2/3 degree grid TL319_t232 to mask missing geolon|lat|u|v|t values -**September 10 2024 :: MARBL_column. Tag 11.8.0** +**September 10 2024 :: MARBL_column. Tag v11.8.0** - Interface for MARBL_column for DART: @@ -48,7 +74,7 @@ Bugfix: - fix for IO for NetCDF files when only some variables have the unlimited dimension -**August 29 2024 :: Bug fixes for shortest_time_between_assimilations and get_close_init. Tag 11.7.1** +**August 29 2024 :: Bug fixes for shortest_time_between_assimilations and get_close_init. Tag v11.7.1** Bug fixes: @@ -64,7 +90,7 @@ Doc fixes: - GitHub template for reporting documentation issues -**August 26 2024 :: KQCEF. Tag 11.7.0** +**August 26 2024 :: KQCEF. Tag v11.7.0** - Adds a Quantile-Conserving Ensemble Filter Based on Kernel-Density Estimation to DART. - New distribution module kde_distribution_mod. diff --git a/assimilation_code/modules/assimilation/normal_distribution_mod.f90 b/assimilation_code/modules/assimilation/normal_distribution_mod.f90 index 6b0656c62d..b3a7c33290 100644 --- a/assimilation_code/modules/assimilation/normal_distribution_mod.f90 +++ b/assimilation_code/modules/assimilation/normal_distribution_mod.f90 @@ -493,51 +493,6 @@ subroutine set_normal_params_from_ens(ens, num, p) end subroutine set_normal_params_from_ens -!------------------------------------------------------------------------ -subroutine inv_cdf_quadrature_like(quantiles, ens, likelihood, ens_size, cdf, p, x_out) - -interface - function cdf(x, p) - use types_mod, only : r8 - use distribution_params_mod, only : distribution_params_type - real(r8) :: cdf - real(r8), intent(in) :: x - type(distribution_params_type), intent(in) :: p - end function -end interface - -integer, intent(in) :: ens_size -real(r8), intent(in) :: quantiles(ens_size) -real(r8), intent(in) :: ens(ens_size) -real(r8), intent(in) :: likelihood(ens_size) -type(distribution_params_type), intent(in) :: p -real(r8), intent(out) :: x_out(ens_size) - -integer :: i -real(r8) :: quad_like(ens_size + 1), q_ens(ens_size + 1) - -! Assume that the quantiles and the corresponding ens are sorted - -! Get the likelihood for each of the ens_size + 1 intervals -do i = 2, ens_size - quad_like(i) = (likelihood(i - 1) + likelihood(i)) / 2.0_r8 -end do -quad_like(1) = likelihood(1) -quad_like(ens_size + 1) = likelihood(ens_size) - -! Compute the quantiles at the ensemble boundaries for the posterior -q_ens(1) = quad_like(1) * quantiles(1) -do i = 2, ens_size - q_ens(i) = q_ens(i - 1) + quad_like(i) * (quantiles(i) - quantiles(i - 1)) -end do -q_ens(ens_size + 1) = q_ens(ens_size) + & - quad_like(ens_size + 1) * (1.0_r8 - quantiles(ens_size)) - -! Normalize so that this is a posterior cdf -q_ens = q_ens / q_ens(ens_size + 1) - -end subroutine inv_cdf_quadrature_like - !------------------------------------------------------------------------ end module normal_distribution_mod diff --git a/assimilation_code/modules/utilities/cray_win_mod.f90 b/assimilation_code/modules/utilities/cray_win_mod.f90 deleted file mode 100644 index fa04c1284b..0000000000 --- a/assimilation_code/modules/utilities/cray_win_mod.f90 +++ /dev/null @@ -1,219 +0,0 @@ -! DART software - Copyright UCAR. This open source software is provided -! by UCAR, "as is", without charge, subject to all terms of use at -! http://www.image.ucar.edu/DAReS/DART/DART_download - -!> Contains the window information for the state. Two windows: -!> One for all copies, one for the mean. -!> Not sure whether we should just have one window to avoid multiple synchronizations. - -module window_mod - -!> \defgroup window window_mod -!> @{ -use mpi_utilities_mod, only : datasize, my_task_id -use types_mod, only : r8 -use ensemble_manager_mod, only : ensemble_type, map_pe_to_task, get_var_owner_index, & - copies_in_window, init_ensemble_manager, & - get_allow_transpose, end_ensemble_manager, & - set_num_extra_copies, all_copies_to_all_vars, & - all_vars_to_all_copies - -use mpi - -implicit none - -private -public :: create_mean_window, create_state_window, free_mean_window, & - free_state_window, data_count, mean_win, state_win, current_win, & - mean_ens_handle, NO_WINDOW, MEAN_WINDOW, STATE_WINDOW - -! mpi window handles -integer :: state_win !! window for the forward operator -integer :: mean_win !! window for the mean -integer :: current_win !! keep track of current window, start out assuming an invalid window -!>@todo the number of copies in the window is sloppy. You need to make this better. - -! parameters for keeping track of which window is open -integer, parameter :: NO_WINDOW = -1 -integer, parameter :: MEAN_WINDOW = 0 -integer, parameter :: STATE_WINDOW = 2 - -integer :: data_count !! number of copies in the window -integer(KIND=MPI_ADDRESS_KIND) window_size -logical :: use_distributed_mean = .false. ! initialize to false - -real(r8) :: duplicate_state(*) ! duplicate array for cray pointer fwd -pointer(a, duplicate_state) - -real(r8) :: duplicate_mean(*) ! duplicate array for cray pointer vert convert -pointer(b, duplicate_mean) -type(ensemble_type) :: mean_ens_handle - -contains - -!------------------------------------------------------------- -!> Create the window for the ensemble complete state vector -!> I think you have to pass it the state ensemble handle -subroutine create_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handle) - -type(ensemble_type), intent(inout) :: state_ens_handle !< state ensemble handle -type(ensemble_type), intent(inout), optional :: fwd_op_ens_handle -type(ensemble_type), intent(inout), optional :: qc_ens_handle - -integer :: ii, jj, count, ierr -integer :: bytesize !< size in bytes of each element in the window -integer :: my_num_vars - -! Find out how many copies to put in the window -! copies_in_window is not necessarily equal to ens_handle%num_copies -data_count = copies_in_window(state_ens_handle) - -if (get_allow_transpose(state_ens_handle)) then - call all_copies_to_all_vars(state_ens_handle) - if (present(fwd_op_ens_handle)) then - call all_copies_to_all_vars(fwd_op_ens_handle) - endif - if (present(qc_ens_handle)) then - call all_copies_to_all_vars(qc_ens_handle) - endif -else - ! find how many variables I have - my_num_vars = state_ens_handle%my_num_vars - - ! allocate some RDMA accessible memory - ! using MPI_ALLOC_MEM because the MPI standard allows vendors to require MPI_ALLOC_MEM for remote memory access - call mpi_type_size(datasize, bytesize, ierr) - window_size = my_num_vars*data_count*bytesize - a = malloc(my_num_vars*data_count) - call MPI_ALLOC_MEM(window_size, MPI_INFO_NULL, a, ierr) - - count = 1 - ! create a duplicate copies array for remote memory access - ! Doing this because you cannot use a cray pointer with an allocatable array - ! Can't do array assignment with a cray pointer, so you need to loop - do ii = 1, my_num_vars - do jj = 1, data_count - duplicate_state(count) = state_ens_handle%copies(jj,ii) - count = count + 1 - enddo - enddo - - ! Expose local memory to RMA operation by other processes in a communicator. - call mpi_win_create(duplicate_state, window_size, bytesize, MPI_INFO_NULL, mpi_comm_world, state_win, ierr) - -endif - -! Set the current window to the state window -current_win = STATE_WINDOW - -data_count = copies_in_window(state_ens_handle) - -end subroutine create_state_window - -!------------------------------------------------------------- -!> Create the window for the ensemble complete state vector -!> I think you have to pass it the state ensemble handle -subroutine create_mean_window(state_ens_handle, mean_copy, distribute_mean) - -type(ensemble_type), intent(in) :: state_ens_handle -integer, intent(in) :: mean_copy -logical, intent(in) :: distribute_mean - -integer :: ii, ierr -integer :: bytesize -integer :: my_num_vars - -! find out how many variables I have -my_num_vars = state_ens_handle%my_num_vars - -! create an ensemble handle of just the mean copy. -use_distributed_mean = distribute_mean - -if (use_distributed_mean) then - - call init_ensemble_manager(mean_ens_handle, 1, state_ens_handle%num_vars) ! distributed ensemble - call set_num_extra_copies(mean_ens_handle, 0) - mean_ens_handle%copies(1,:) = state_ens_handle%copies(mean_copy, :) - - ! find out how many variables I have - my_num_vars = mean_ens_handle%my_num_vars - call mpi_type_size(datasize, bytesize, ierr) - window_size = my_num_vars*bytesize - ! allocate some RDMA accessible memory - ! using MPI_ALLOC_MEM because the MPI standard allows vendors to require MPI_ALLOC_MEM for remote memory access - ! Have a look at MPI-3, I think this removes cray pointers. - b = malloc(my_num_vars) - call MPI_ALLOC_MEM(window_size, MPI_INFO_NULL, b, ierr) - - do ii = 1, my_num_vars - duplicate_mean(ii) = mean_ens_handle%copies(1, ii) - enddo - - ! Expose local memory to RMA operation by other processes in a communicator. - call mpi_win_create(duplicate_mean, window_size, bytesize, MPI_INFO_NULL, mpi_comm_world, mean_win, ierr) - -else - - call init_ensemble_manager(mean_ens_handle, 1, state_ens_handle%num_vars, transpose_type_in = 3) - call set_num_extra_copies(mean_ens_handle, 0) - mean_ens_handle%copies(1,:) = state_ens_handle%copies(mean_copy, :) - call all_copies_to_all_vars(mean_ens_handle) ! this is a transpose-duplicate - -endif - -! Set the current window to the state window -current_win = MEAN_WINDOW - -data_count = copies_in_window(mean_ens_handle) ! One - -end subroutine create_mean_window - -!------------------------------------------------------------- -!> End epoch of state access. -!> Need to transpose qc and fwd operator back to copy complete -subroutine free_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handle) - -type(ensemble_type), intent(inout) :: state_ens_handle -type(ensemble_type), intent(inout), optional :: fwd_op_ens_handle -type(ensemble_type), intent(inout), optional :: qc_ens_handle - -integer :: ierr - -if(get_allow_transpose(state_ens_handle)) then ! the forward operators were done var complete - !transpose back if allowing transposes - if (present(fwd_op_ens_handle)) & - call all_vars_to_all_copies(fwd_op_ens_handle) - if (present(qc_ens_handle)) & - call all_vars_to_all_copies(qc_ens_handle) -else - ! close mpi window - call mpi_win_free(state_win, ierr) - call MPI_FREE_MEM(duplicate_state, ierr) -endif - -current_win = NO_WINDOW - -end subroutine free_state_window - -!--------------------------------------------------------- -!> Free the mpi window -subroutine free_mean_window() - -integer :: ierr - -if(get_allow_transpose(mean_ens_handle)) then - call end_ensemble_manager(mean_ens_handle) -else - call mpi_win_free(mean_win, ierr) - call MPI_FREE_MEM(duplicate_mean, ierr) - call end_ensemble_manager(mean_ens_handle) -endif - -current_win = NO_WINDOW - -end subroutine free_mean_window - -!--------------------------------------------------------- -!> @} -end module window_mod - diff --git a/assimilation_code/modules/utilities/distributed_state_mod.f90 b/assimilation_code/modules/utilities/distributed_state_mod.f90 index 340afa1510..9038bc3f7d 100644 --- a/assimilation_code/modules/utilities/distributed_state_mod.f90 +++ b/assimilation_code/modules/utilities/distributed_state_mod.f90 @@ -110,7 +110,7 @@ subroutine get_fwd(x, my_index, state_ens_handle) !x = get_local_state(element_index) x = state_ens_handle%copies(1:data_count, element_index) else - call get_from_fwd(owner_of_state, state_win, element_index, data_count, x) + call get_from_fwd(owner_of_state, state_win, element_index, state_ens_handle%num_copies, data_count, x) endif endif diff --git a/assimilation_code/modules/utilities/mpi_utilities_mod.f90 b/assimilation_code/modules/utilities/mpi_utilities_mod.f90 index 90797018e8..5702541ed2 100644 --- a/assimilation_code/modules/utilities/mpi_utilities_mod.f90 +++ b/assimilation_code/modules/utilities/mpi_utilities_mod.f90 @@ -12,11 +12,6 @@ !> and allows programs to swap in the null version to compile the !> same source files into a serial program. !> -!> The names of these routines were intentionally picked to be -!> more descriptive to someone who doesn't the MPI interfaces. -!> e.g. MPI_AllReduce() may not immediately tell a user what -!> it does, but broadcast_minmax() is hopefully more helpful. -!> !> If you add any routines or change any arguments in this file !> you must make the same changes in the null version. These two !> modules have the same module name and must have identical @@ -138,9 +133,8 @@ module mpi_utilities_mod task_sync, array_broadcast, send_to, receive_from, iam_task0, & broadcast_send, broadcast_recv, shell_execute, sleep_seconds, & sum_across_tasks, get_dart_mpi_comm, datasize, send_minmax_to, & - get_from_fwd, get_from_mean, broadcast_minmax, broadcast_flag, & - start_mpi_timer, read_mpi_timer, send_sum_to, get_global_max, & - all_reduce_min_max ! deprecated, replace by broadcast_minmax + get_from_fwd, get_from_mean, broadcast_flag, start_mpi_timer, & + read_mpi_timer, send_sum_to, get_global_max, all_reduce_min_max character(len=*), parameter :: source = 'mpi_utilities_mod.f90' @@ -1467,26 +1461,11 @@ end subroutine send_minmax_to !----------------------------------------------------------------------------- -!> cover routine which is deprecated. when all user code replaces this -!> with broadcast_minmax(), remove this. - -subroutine all_reduce_min_max(min_var, max_var, num_elements) - -integer, intent(in) :: num_elements -real(r8), intent(inout) :: min_var(num_elements) -real(r8), intent(inout) :: max_var(num_elements) - -call broadcast_minmax(min_var, max_var, num_elements) - -end subroutine all_reduce_min_max - -!----------------------------------------------------------------------------- - !> Find min and max of each element of an array, put the result on every task. !> Overwrites arrays min_var, max_var with the minimum and maximum for each !> element across all tasks. -subroutine broadcast_minmax(min_var, max_var, num_elements) +subroutine all_reduce_min_max(min_var, max_var, num_elements) integer, intent(in) :: num_elements real(r8), intent(inout) :: min_var(num_elements) @@ -1496,13 +1475,13 @@ subroutine broadcast_minmax(min_var, max_var, num_elements) if ( .not. module_initialized ) then write(errstring, *) 'initialize_mpi_utilities() must be called first' - call error_handler(E_ERR,'broadcast_minmax', errstring, source) + call error_handler(E_ERR,'all_reduce_min_max', errstring, source) endif call mpi_allreduce(MPI_IN_PLACE, min_var, num_elements, datasize, MPI_MIN, get_dart_mpi_comm(), errcode) call mpi_allreduce(MPI_IN_PLACE, max_var, num_elements, datasize, MPI_MAX, get_dart_mpi_comm(), errcode) -end subroutine broadcast_minmax +end subroutine all_reduce_min_max !----------------------------------------------------------------------------- !> Broadcast logical @@ -1968,13 +1947,14 @@ end subroutine get_from_mean !----------------------------------------------------------------------------- -subroutine get_from_fwd(owner, window, mindex, num_rows, x) +subroutine get_from_fwd(owner, window, mindex, rows_in_window, num_rows, x) integer, intent(in) :: owner ! task in the window that owns the memory integer, intent(in) :: window ! window object integer, intent(in) :: mindex ! index in the tasks memory -integer, intent(in) :: num_rows ! number of rows in the window -real(r8), intent(out) :: x(:) ! result +integer, intent(in) :: rows_in_window ! number of rows in the window +integer, intent(in) :: num_rows ! number of rows to get from the window +real(r8), intent(out) :: x(num_rows) ! result integer(KIND=MPI_ADDRESS_KIND) :: target_disp integer :: errcode @@ -1983,7 +1963,7 @@ subroutine get_from_fwd(owner, window, mindex, num_rows, x) ! to have occured until the call to mpi_win_unlock. ! => Don't do anything with x in between mpi_get and mpi_win_lock -target_disp = (mindex - 1)*num_rows +target_disp = (mindex - 1)*rows_in_window call mpi_win_lock(MPI_LOCK_SHARED, owner, 0, window, errcode) call mpi_get(x, num_rows, datasize, owner, target_disp, num_rows, datasize, window, errcode) call mpi_win_unlock(owner, window, errcode) diff --git a/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 b/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 index dee8c618d1..5f0ac14204 100644 --- a/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 +++ b/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 @@ -12,11 +12,6 @@ !> and allows programs to swap in the null version to compile the !> same source files into a serial program. !> -!> The names of these routines were intentionally picked to be -!> more descriptive to someone who doesn't the MPI interfaces. -!> e.g. MPI_AllReduce() may not immediately tell a user what -!> it does, but broadcast_minmax() is hopefully more helpful. -!> !> If you add any routines or change any arguments in this file !> you must make the same changes in the null version. These two !> modules have the same module name and must have identical @@ -138,9 +133,8 @@ module mpi_utilities_mod task_sync, array_broadcast, send_to, receive_from, iam_task0, & broadcast_send, broadcast_recv, shell_execute, sleep_seconds, & sum_across_tasks, get_dart_mpi_comm, datasize, send_minmax_to, & - get_from_fwd, get_from_mean, broadcast_minmax, broadcast_flag, & - start_mpi_timer, read_mpi_timer, send_sum_to, get_global_max, & - all_reduce_min_max ! deprecated, replace by broadcast_minmax + get_from_fwd, get_from_mean, broadcast_flag, start_mpi_timer, & + read_mpi_timer, send_sum_to, get_global_max, all_reduce_min_max character(len=*), parameter :: source = 'mpi_utilities_mod.f90' @@ -1467,26 +1461,11 @@ end subroutine send_minmax_to !----------------------------------------------------------------------------- -!> cover routine which is deprecated. when all user code replaces this -!> with broadcast_minmax(), remove this. - -subroutine all_reduce_min_max(min_var, max_var, num_elements) - -integer, intent(in) :: num_elements -real(r8), intent(inout) :: min_var(num_elements) -real(r8), intent(inout) :: max_var(num_elements) - -call broadcast_minmax(min_var, max_var, num_elements) - -end subroutine all_reduce_min_max - -!----------------------------------------------------------------------------- - !> Find min and max of each element of an array, put the result on every task. !> Overwrites arrays min_var, max_var with the minimum and maximum for each !> element across all tasks. -subroutine broadcast_minmax(min_var, max_var, num_elements) +subroutine all_reduce_min_max(min_var, max_var, num_elements) integer, intent(in) :: num_elements real(r8), intent(inout) :: min_var(num_elements) @@ -1496,13 +1475,13 @@ subroutine broadcast_minmax(min_var, max_var, num_elements) if ( .not. module_initialized ) then write(errstring, *) 'initialize_mpi_utilities() must be called first' - call error_handler(E_ERR,'broadcast_minmax', errstring, source) + call error_handler(E_ERR,'all_reduce_min_max', errstring, source) endif call mpi_allreduce(MPI_IN_PLACE, min_var, num_elements, datasize, MPI_MIN, get_dart_mpi_comm(), errcode) call mpi_allreduce(MPI_IN_PLACE, max_var, num_elements, datasize, MPI_MAX, get_dart_mpi_comm(), errcode) -end subroutine broadcast_minmax +end subroutine all_reduce_min_max !----------------------------------------------------------------------------- !> Broadcast logical @@ -1969,12 +1948,13 @@ end subroutine get_from_mean !----------------------------------------------------------------------------- -subroutine get_from_fwd(owner, window, mindex, num_rows, x) +subroutine get_from_fwd(owner, window, mindex, rows_in_window, num_rows, x) integer, intent(in) :: owner ! task in the window that owns the memory type(MPI_Win), intent(in) :: window ! window object integer, intent(in) :: mindex ! index in the tasks memory -integer, intent(in) :: num_rows ! number of rows in the window +integer, intent(in) :: rows_in_window ! number of rows in the window +integer, intent(in) :: num_rows ! number of rows to get from the window real(r8), intent(out) :: x(num_rows) ! result integer(KIND=MPI_ADDRESS_KIND) :: target_disp @@ -1984,7 +1964,7 @@ subroutine get_from_fwd(owner, window, mindex, num_rows, x) ! to have occured until the call to mpi_win_unlock. ! => Don't do anything with x in between mpi_get and mpi_win_lock -target_disp = (mindex - 1)*num_rows +target_disp = (mindex - 1)*rows_in_window call mpi_win_lock(MPI_LOCK_SHARED, owner, 0, window, errcode) call mpi_get(x, num_rows, datasize, owner, target_disp, num_rows, datasize, window, errcode) call mpi_win_unlock(owner, window, errcode) diff --git a/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 b/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 index cc729c23aa..4c6e781d49 100644 --- a/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 +++ b/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 @@ -109,9 +109,8 @@ module mpi_utilities_mod task_sync, array_broadcast, send_to, receive_from, iam_task0, & broadcast_send, broadcast_recv, shell_execute, sleep_seconds, & sum_across_tasks, get_dart_mpi_comm, datasize, send_minmax_to, & - get_from_fwd, get_from_mean, broadcast_minmax, broadcast_flag, & - start_mpi_timer, read_mpi_timer, send_sum_to, get_global_max, & - all_reduce_min_max ! deprecated, replace by broadcast_minmax + get_from_fwd, get_from_mean, broadcast_flag, start_mpi_timer, & + read_mpi_timer, send_sum_to, get_global_max, all_reduce_min_max character(len=*), parameter :: source = 'null_mpi_utilities_mod.f90' @@ -432,32 +431,17 @@ end subroutine send_minmax_to !----------------------------------------------------------------------------- -!> cover routine which is deprecated. when all user code replaces this -!> with broadcast_minmax(), remove this. - -subroutine all_reduce_min_max(min_var, max_var, num_elements) - -integer, intent(in) :: num_elements -real(r8), intent(inout) :: min_var(num_elements) -real(r8), intent(inout) :: max_var(num_elements) - -call broadcast_minmax(min_var, max_var, num_elements) - -end subroutine all_reduce_min_max - -!----------------------------------------------------------------------------- - !> Find min and max of each element of an array across tasks, put the result on every task. !> For this null_mpi_version min_var and max_var are unchanged because there is !> only 1 task. -subroutine broadcast_minmax(min_var, max_var, num_elements) +subroutine all_reduce_min_max(min_var, max_var, num_elements) integer, intent(in) :: num_elements real(r8), intent(inout) :: min_var(num_elements) real(r8), intent(inout) :: max_var(num_elements) -end subroutine broadcast_minmax +end subroutine all_reduce_min_max !----------------------------------------------------------------------------- @@ -634,13 +618,14 @@ end subroutine get_from_mean !----------------------------------------------------------------------------- -subroutine get_from_fwd(owner, window, mindex, num_rows, x) +subroutine get_from_fwd(owner, window, mindex, rows_in_window, num_rows, x) integer, intent(in) :: owner ! task in the window that owns the memory integer, intent(in) :: window ! window object integer, intent(in) :: mindex ! index in the tasks memory -integer, intent(in) :: num_rows ! number of rows in the window -real(r8), intent(out) :: x(num_rows) ! result +integer, intent(in) :: rows_in_window ! number of rows in the window +integer, intent(in) :: num_rows ! number of rows to get from the window +real(r8), intent(out) :: x(num_rows) ! result call error_handler(E_ERR,'get_from_fwd', 'cannot be used in serial mode', source) diff --git a/assimilation_code/modules/utilities/time_manager_mod.f90 b/assimilation_code/modules/utilities/time_manager_mod.f90 index db35f4fd33..cdc88517aa 100644 --- a/assimilation_code/modules/utilities/time_manager_mod.f90 +++ b/assimilation_code/modules/utilities/time_manager_mod.f90 @@ -665,7 +665,8 @@ end function repeat_alarm !========================================================================= subroutine set_calendar_type_integer(mytype) - +!------------------------------------------------------------------------ +! ! Selects calendar for default mapping from time to date - if you know ! the magic integer for the calendar of interest. @@ -684,7 +685,8 @@ end subroutine set_calendar_type_integer subroutine set_calendar_type_string(calstring) - +!------------------------------------------------------------------------ +! ! Selects calendar for default mapping from time to date - given a string. character(len=*), intent(in) :: calstring @@ -693,8 +695,6 @@ subroutine set_calendar_type_string(calstring) character(len=len(calstring)) :: str1 character(len=max_calendar_string_length) :: cstring -logical :: found_calendar = .false. -integer :: i if ( .not. module_initialized ) call time_manager_init @@ -714,49 +714,25 @@ subroutine set_calendar_type_string(calstring) ! We must check for the gregorian_mars calendar before ! the gregorian calendar for similar reasons. -WhichCalendar : do i = 0, max_type - - if ( cstring == 'NO_CALENDAR' ) then - calendar_type = NO_CALENDAR - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'NO CALENDAR' ) then ! allow this as a synonym - calendar_type = NO_CALENDAR - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'NONE' ) then ! also allow this - calendar_type = NO_CALENDAR - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'THIRTY_DAY_MONTHS' ) then - calendar_type = THIRTY_DAY_MONTHS - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'JULIAN' ) then - calendar_type = JULIAN - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'NOLEAP' ) then - calendar_type = NOLEAP - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'GREGORIAN_MARS' ) then - calendar_type = GREGORIAN_MARS - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'SOLAR_MARS' ) then - calendar_type = SOLAR_MARS - found_calendar = .true. - exit WhichCalendar - elseif ( cstring == 'GREGORIAN' ) then - calendar_type = GREGORIAN - found_calendar = .true. - exit WhichCalendar - endif - -enddo WhichCalendar - -if( .not. found_calendar ) then +if ( cstring == 'NO_CALENDAR' ) then + calendar_type = NO_CALENDAR +elseif ( cstring == 'NO CALENDAR' ) then ! allow this as a synonym + calendar_type = NO_CALENDAR +elseif ( cstring == 'NONE' ) then ! also allow this + calendar_type = NO_CALENDAR +elseif ( cstring == 'THIRTY_DAY_MONTHS' ) then + calendar_type = THIRTY_DAY_MONTHS +elseif ( cstring == 'JULIAN' ) then + calendar_type = JULIAN +elseif ( cstring == 'NOLEAP' ) then + calendar_type = NOLEAP +elseif ( cstring == 'GREGORIAN_MARS' ) then + calendar_type = GREGORIAN_MARS +elseif ( cstring == 'SOLAR_MARS' ) then + calendar_type = SOLAR_MARS +elseif ( cstring == 'GREGORIAN' ) then + calendar_type = GREGORIAN +else write(errstring,*)'Unknown calendar ',calstring call error_handler(E_ERR,'set_calendar_type_string',errstring,source) endif @@ -785,23 +761,19 @@ subroutine get_calendar_string(mystring) ! ! Returns default calendar type for mapping from time to date. -character(len=*), intent(OUT) :: mystring - -integer :: i +character(len=*), intent(out) :: mystring if ( .not. module_initialized ) call time_manager_init -mystring = ' ' +mystring = '' -do i = 0,max_type - if (calendar_type == JULIAN) mystring = 'JULIAN' - if (calendar_type == NOLEAP) mystring = 'NOLEAP' - if (calendar_type == GREGORIAN) mystring = 'GREGORIAN' - if (calendar_type == NO_CALENDAR) mystring = 'NO_CALENDAR' - if (calendar_type == GREGORIAN_MARS) mystring = 'GREGORIAN_MARS' - if (calendar_type == SOLAR_MARS) mystring = 'SOLAR_MARS' - if (calendar_type == THIRTY_DAY_MONTHS) mystring = 'THIRTY_DAY_MONTHS' -enddo +if (calendar_type == JULIAN) mystring = 'JULIAN' +if (calendar_type == NOLEAP) mystring = 'NOLEAP' +if (calendar_type == GREGORIAN) mystring = 'GREGORIAN' +if (calendar_type == NO_CALENDAR) mystring = 'NO_CALENDAR' +if (calendar_type == GREGORIAN_MARS) mystring = 'GREGORIAN_MARS' +if (calendar_type == SOLAR_MARS) mystring = 'SOLAR_MARS' +if (calendar_type == THIRTY_DAY_MONTHS) mystring = 'THIRTY_DAY_MONTHS' if (len_trim(mystring) < 3) then write(errstring,*)'unknown calendar type ', calendar_type diff --git a/assimilation_code/modules/utilities/no_cray_win_mod.f90 b/assimilation_code/modules/utilities/win_mod.f90 similarity index 87% rename from assimilation_code/modules/utilities/no_cray_win_mod.f90 rename to assimilation_code/modules/utilities/win_mod.f90 index 6090621f83..e965261006 100644 --- a/assimilation_code/modules/utilities/no_cray_win_mod.f90 +++ b/assimilation_code/modules/utilities/win_mod.f90 @@ -2,7 +2,6 @@ ! by UCAR, "as is", without charge, subject to all terms of use at ! http://www.image.ucar.edu/DAReS/DART/DART_download -!> Window without cray pointer. Should you point the window at contigous memory? module window_mod !> \defgroup window window_mod @@ -30,21 +29,16 @@ module window_mod integer :: current_win !< keep track of current window, start out assuming an invalid window ! parameters for keeping track of which window is open -!>@todo should this be in the window_mod? you will have to change in both cray -!> and non cray versions integer, parameter :: NO_WINDOW = -1 integer, parameter :: MEAN_WINDOW = 0 integer, parameter :: STATE_WINDOW = 2 -integer :: data_count !! number of copies in the window +integer :: data_count !! number of copies required integer(KIND=MPI_ADDRESS_KIND) :: window_size logical :: use_distributed_mean = .false. ! initialize to false ! Global memory to stick the mpi window to. ! Need a simply contiguous piece of memory to pass to mpi_win_create -! Openmpi 1.10.0 will not compile with ifort 16 if -! you create a window with a 2d array. -real(r8), allocatable :: contiguous_fwd(:) real(r8), allocatable :: mean_1d(:) type(ensemble_type) :: mean_ens_handle @@ -65,8 +59,7 @@ subroutine create_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handl integer :: bytesize !< size in bytes of each element in the window integer :: my_num_vars !< my number of vars -! Find out how many copies to put in the window -! copies_in_window is not necessarily equal to ens_handle%num_copies +! Find out how many copies to get, maybe different to state_ens_handle%num_copies data_count = copies_in_window(state_ens_handle) if (get_allow_transpose(state_ens_handle)) then @@ -82,20 +75,15 @@ subroutine create_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handl my_num_vars = state_ens_handle%my_num_vars call mpi_type_size(datasize, bytesize, ierr) - window_size = my_num_vars*data_count*bytesize - - allocate(contiguous_fwd(data_count*my_num_vars)) - contiguous_fwd = reshape(state_ens_handle%copies(1:data_count, :), (/my_num_vars*data_count/)) + window_size = my_num_vars*state_ens_handle%num_copies*bytesize ! Expose local memory to RMA operation by other processes in a communicator. - call mpi_win_create(contiguous_fwd, window_size, bytesize, MPI_INFO_NULL, get_dart_mpi_comm(), state_win, ierr) + call mpi_win_create(state_ens_handle%copies, window_size, bytesize, MPI_INFO_NULL, get_dart_mpi_comm(), state_win, ierr) endif ! Set the current window to the state window current_win = STATE_WINDOW -data_count = copies_in_window(state_ens_handle) - end subroutine create_state_window !------------------------------------------------------------- @@ -163,7 +151,6 @@ subroutine free_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handle) else ! close mpi window call mpi_win_free(state_win, ierr) - deallocate(contiguous_fwd) endif current_win = NO_WINDOW diff --git a/assimilation_code/modules/utilities/no_cray_winf08_mod.f90 b/assimilation_code/modules/utilities/winf08_mod.f90 similarity index 87% rename from assimilation_code/modules/utilities/no_cray_winf08_mod.f90 rename to assimilation_code/modules/utilities/winf08_mod.f90 index 7ecfd90830..54f14de3db 100644 --- a/assimilation_code/modules/utilities/no_cray_winf08_mod.f90 +++ b/assimilation_code/modules/utilities/winf08_mod.f90 @@ -2,7 +2,6 @@ ! by UCAR, "as is", without charge, subject to all terms of use at ! http://www.image.ucar.edu/DAReS/DART/DART_download -!> Window without cray pointer. Should you point the window at contigous memory? module window_mod !> \defgroup window window_mod @@ -30,21 +29,16 @@ module window_mod integer :: current_win !< keep track of current window, start out assuming an invalid window ! parameters for keeping track of which window is open -!>@todo should this be in the window_mod? you will have to change in both cray -!> and non cray versions integer, parameter :: NO_WINDOW = -1 integer, parameter :: MEAN_WINDOW = 0 integer, parameter :: STATE_WINDOW = 2 -integer :: data_count !! number of copies in the window +integer :: data_count !! number of copies required integer(KIND=MPI_ADDRESS_KIND) :: window_size logical :: use_distributed_mean = .false. ! initialize to false ! Global memory to stick the mpi window to. ! Need a simply contiguous piece of memory to pass to mpi_win_create -! Openmpi 1.10.0 will not compile with ifort 16 if -! you create a window with a 2d array. -real(r8), allocatable :: contiguous_fwd(:) real(r8), allocatable :: mean_1d(:) type(ensemble_type) :: mean_ens_handle @@ -65,8 +59,7 @@ subroutine create_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handl integer :: bytesize !< size in bytes of each element in the window integer :: my_num_vars !< my number of vars -! Find out how many copies to put in the window -! copies_in_window is not necessarily equal to ens_handle%num_copies +! Find out how many copies to get, maybe different to state_ens_handle%num_copies data_count = copies_in_window(state_ens_handle) if (get_allow_transpose(state_ens_handle)) then @@ -82,20 +75,15 @@ subroutine create_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handl my_num_vars = state_ens_handle%my_num_vars call mpi_type_size(datasize, bytesize, ierr) - window_size = my_num_vars*data_count*bytesize - - allocate(contiguous_fwd(data_count*my_num_vars)) - contiguous_fwd = reshape(state_ens_handle%copies(1:data_count, :), (/my_num_vars*data_count/)) + window_size = my_num_vars*state_ens_handle%num_copies*bytesize ! Expose local memory to RMA operation by other processes in a communicator. - call mpi_win_create(contiguous_fwd, window_size, bytesize, MPI_INFO_NULL, get_dart_mpi_comm(), state_win, ierr) + call mpi_win_create(state_ens_handle%copies, window_size, bytesize, MPI_INFO_NULL, get_dart_mpi_comm(), state_win, ierr) endif ! Set the current window to the state window current_win = STATE_WINDOW -data_count = copies_in_window(state_ens_handle) - end subroutine create_state_window !------------------------------------------------------------- @@ -163,7 +151,6 @@ subroutine free_state_window(state_ens_handle, fwd_op_ens_handle, qc_ens_handle) else ! close mpi window call mpi_win_free(state_win, ierr) - deallocate(contiguous_fwd) endif current_win = NO_WINDOW diff --git a/build_templates/buildconvfunctions.sh b/build_templates/buildconvfunctions.sh index 8b194c7c90..91de56d7e5 100644 --- a/build_templates/buildconvfunctions.sh +++ b/build_templates/buildconvfunctions.sh @@ -25,7 +25,6 @@ # # The GSI obs converter needs mpi # mpisrc="null_mpi" -# windowsrc="" # m="" #------------------------- set -e @@ -34,7 +33,6 @@ source "$DART"/build_templates/buildpreprocess.sh # Defaults mpisrc="null_mpi" -windowsrc="" m="" LIBRARIES="" EXTRA="" @@ -123,28 +121,22 @@ local mpi="$DART"/assimilation_code/modules/utilities/mpi_utilities_mod.f90 local mpif08="$DART"/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 local nullmpi="$DART"/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 local nullwin="$DART"/assimilation_code/modules/utilities/null_win_mod.f90 -local craywin="$DART"/assimilation_code/modules/utilities/cray_win_mod.f90 -local nocraywin="$DART"/assimilation_code/modules/utilities/no_cray_win_mod.f90 -local no_cray_winf08="$DART"/assimilation_code/modules/utilities/no_cray_winf08_mod.f90 +local win="$DART"/assimilation_code/modules/utilities/win_mod.f90 +local winf08="$DART"/assimilation_code/modules/utilities/winf08_mod.f90 if [ "$mpisrc" == "mpi" ]; then core=${core//$nullmpi/} core=${core//$nullwin/} core=${core//$mpif08/} - core=${core//$no_cray_winf08/} - if [ "$windowsrc" == "craywin" ]; then - core=${core//$nocraywin/} - else #nocraywin - core=${core//$craywin/} - fi + core=${core//$winf08/} + else #nompi core=${core//$mpi/} core=${core//$mpif08/} - core=${core//$nocraywin/} - core=${core//$no_cray_winf08/} - core=${core//$craywin/} + core=${core//$win/} + core=${core//$winf08/} fi convsrc="${core} ${conv} ${obserrsrc} ${modelsrc} ${misc} ${loc}" diff --git a/build_templates/buildfunctions.sh b/build_templates/buildfunctions.sh index f8188a8cad..1440922da4 100644 --- a/build_templates/buildfunctions.sh +++ b/build_templates/buildfunctions.sh @@ -105,7 +105,6 @@ fi # Default to build with mpi (non f08 version) mpisrc=mpi -windowsrc=no_cray_win m="-w" # mkmf wrapper arg # if the first argument is help, nompi, mpi, mpif08, clean @@ -116,20 +115,17 @@ case $1 in nompi) mpisrc="null_mpi" - windowsrc="" m="" shift 1 ;; mpi) mpisrc="mpi" - windowsrc="no_cray_win" shift 1 ;; mpif08) mpisrc="mpif08" - windowsrc="no_cray_winf08" shift 1 ;; @@ -168,37 +164,29 @@ local mpi="$DART"/assimilation_code/modules/utilities/mpi_utilities_mod.f90 local mpif08="$DART"/assimilation_code/modules/utilities/mpif08_utilities_mod.f90 local nullmpi="$DART"/assimilation_code/modules/utilities/null_mpi_utilities_mod.f90 local nullwin="$DART"/assimilation_code/modules/utilities/null_win_mod.f90 -local craywin="$DART"/assimilation_code/modules/utilities/cray_win_mod.f90 -local nocraywin="$DART"/assimilation_code/modules/utilities/no_cray_win_mod.f90 -local no_cray_winf08="$DART"/assimilation_code/modules/utilities/no_cray_winf08_mod.f90 +local win="$DART"/assimilation_code/modules/utilities/win_mod.f90 +local winf08="$DART"/assimilation_code/modules/utilities/winf08_mod.f90 if [ "$mpisrc" == "mpi" ]; then core=${core//$nullmpi/} core=${core//$nullwin/} core=${core//$mpif08/} - core=${core//$no_cray_winf08/} - if [ "$windowsrc" == "craywin" ]; then - core=${core//$nocraywin/} - else #nocraywin - core=${core//$craywin/} - fi + core=${core//$winf08/} elif [ "$mpisrc" == "mpif08" ]; then core=${core//$nullmpi/} core=${core//$nullwin/} core=${core//$mpi/} - core=${core//$craywin/} - core=${core//$nocraywin/} + core=${core//$win/} else #nompi core=${core//$mpi/} core=${core//$mpif08/} - core=${core//$nocraywin/} - core=${core//$no_cray_winf08/} - core=${core//$craywin/} + core=${core//$win/} + core=${core//$winf08/} fi dartsrc="${core} ${modelsrc} ${loc} ${misc}" @@ -234,8 +222,10 @@ done function dartbuild() { local program +local devlibs if [ $dev_test -eq 0 ]; then + devlibs="" #look in $program directory for {main}.f90 if [ $1 == "obs_diag" ]; then program=$DART/assimilation_code/programs/obs_diag/$LOCATION @@ -247,11 +237,13 @@ if [ $dev_test -eq 0 ]; then else # For developer tests {main}.f90 is in developer_tests program=$DART/developer_tests/$TEST/$1.f90 + devlibs=$DART/developer_tests/contrib/fortran-testanything fi $DART/build_templates/mkmf -x -a $DART $m -p $1 \ $dartsrc \ $EXTRA \ + $devlibs \ $program } @@ -291,7 +283,6 @@ if [ ! -z "$single_prog" ] ; then # build a single program elif [[ " ${serial_programs[*]} " =~ " ${single_prog} " ]]; then echo "building serial dart program " $single_prog mpisrc="null_mpi" - windowsrc="" m="" findsrc dartbuild $single_prog @@ -304,7 +295,6 @@ if [ ! -z "$single_prog" ] ; then # build a single program elif [[ " ${model_serial_programs[*]} " =~ " ${single_prog} " ]];then echo "building model program" $single_prog mpisrc="null_mpi" - windowsrc="" m="" findsrc modelbuild $single_prog @@ -340,7 +330,6 @@ done [ $mpisrc == "mpi" ] && \rm -f *.o *.mod mpisrc="null_mpi" -windowsrc="" m="" # Serial programs diff --git a/conf.py b/conf.py index 902f7875d9..8365db092b 100644 --- a/conf.py +++ b/conf.py @@ -21,7 +21,7 @@ author = 'Data Assimilation Research Section' # The full version, including alpha/beta/rc tags -release = '11.8.2' +release = '11.8.5' root_doc = 'index' # -- General configuration --------------------------------------------------- diff --git a/developer_tests/contrib/fortran-testanything/LICENSE.txt b/developer_tests/contrib/fortran-testanything/LICENSE.txt new file mode 100644 index 0000000000..aaf4092f79 --- /dev/null +++ b/developer_tests/contrib/fortran-testanything/LICENSE.txt @@ -0,0 +1,14 @@ +Copyright 2015 Dennis Decker Jensen + +Permission to use, copy, modify, and distribute this software for any +purpose with or without fee is hereby granted, provided that the above +copyright notice and this permission notice appear in all copies. + +THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES +WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF +MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR +ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES +WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN +ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF +OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. + diff --git a/developer_tests/contrib/fortran-testanything/is_i.inc b/developer_tests/contrib/fortran-testanything/is_i.inc new file mode 100644 index 0000000000..7f98a09089 --- /dev/null +++ b/developer_tests/contrib/fortran-testanything/is_i.inc @@ -0,0 +1,24 @@ +! Template parameter: wp (working precision) +! Template free identifiers: testline, tests +subroutine is(got, expected, msg) + integer(kind=wp), intent(in) :: got, expected + character(len=*), intent(in), optional :: msg + character(len=:), allocatable :: testmsg, idmsg + character(len=120) gotmsg, expectedmsg + logical good + + if (present(msg)) then + allocate(character(len=len_trim(msg)+20) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,A,A)') 'Failed test: "', trim(msg), '"' + testmsg = trim(msg) + else + allocate(character(len=30) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,I0)') 'Failed test no. ', tests + 1 + testmsg = "" + end if + write (unit=gotmsg, fmt='(A,I0)') ' got: ', got + write (unit=expectedmsg, fmt='(A,I0)') 'expected: ', expected + + good = got == expected + call testline(good, testmsg, idmsg, gotmsg, expectedmsg) +end diff --git a/developer_tests/contrib/fortran-testanything/is_r.inc b/developer_tests/contrib/fortran-testanything/is_r.inc new file mode 100644 index 0000000000..98599716b4 --- /dev/null +++ b/developer_tests/contrib/fortran-testanything/is_r.inc @@ -0,0 +1,83 @@ +! Template parameter: wp (working precision) +! Template free identifiers: testline, tests +subroutine isabs(got, expected, eps, msg) + real(kind=wp), intent(in) :: got, expected + character(len=*), intent(in), optional :: msg + real(kind=wp), intent(in), optional :: eps + character(len=:), allocatable :: testmsg, idmsg + character(len=120) gotmsg, expectedmsg + real(kind=wp) tolerance + logical good + + if (present(msg)) then + allocate(character(len=len_trim(msg)+20) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,A,A)') 'Failed test: "', trim(msg), '"' + testmsg = trim(msg) + else + allocate(character(len=30) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,I0)') 'Failed test no. ', tests + 1 + testmsg = "" + end if + write (unit=gotmsg, fmt='(A,G0)') ' got: ', got + write (unit=expectedmsg, fmt='(A,G0)') 'expected: ', expected + + if (present(eps)) then + tolerance = eps + else + tolerance = epsilon(got) + end if + ! eps = 0.5e-10_wp + ! Absolute accuracy within the 10 least significant digits + good = abs(got - expected) < tolerance + call testline(good, testmsg, idmsg, gotmsg, expectedmsg) +end + +subroutine isrel(got, expected, eps, msg) + real(kind=wp), intent(in) :: got, expected + character(len=*), intent(in), optional :: msg + real(kind=wp), intent(in), optional :: eps + real(kind=wp) tolerance + + ! eps = (abs(a) + abs(b)) * 0.5e-10_wp + ! Relative accuracy within the 10 most significant digits + tolerance = (abs(got) + abs(expected)) + if (present(eps)) then + tolerance = tolerance * eps + else + tolerance = tolerance * epsilon(got) + end if + call isabs(got, expected, tolerance, msg) +end + +subroutine isnear(got, expected, eps, msg) + real(kind=wp), intent(in) :: got, expected + character(len=*), intent(in), optional :: msg + real(kind=wp), intent(in), optional :: eps + character(len=:), allocatable :: testmsg, idmsg + character(len=120) gotmsg, expectedmsg + real(kind=wp) tolerance + logical good + + if (present(msg)) then + allocate(character(len=len_trim(msg)+20) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,A,A)') 'Failed test: "', trim(msg), '"' + testmsg = trim(msg) + else + allocate(character(len=30) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,I0)') 'Failed test no. ', tests + 1 + testmsg = "" + end if + write (unit=gotmsg, fmt='(A,G0)') ' got: ', got + write (unit=expectedmsg, fmt='(A,G0)') 'expected: ', expected + + if (present(eps)) then + tolerance = eps + else + tolerance = epsilon(got) ! minimun eps for which 1 + eps /= 1 + end if + ! Relative accuracy around 1.0_wp + ! Semantics of isnear means using <=, and not <, c.f. epsilon(got) + good = abs(got / expected - 1.0_wp) <= tolerance + call testline(good, testmsg, idmsg, gotmsg, expectedmsg) +end + diff --git a/developer_tests/contrib/fortran-testanything/test.f90 b/developer_tests/contrib/fortran-testanything/test.f90 new file mode 100644 index 0000000000..5b565779df --- /dev/null +++ b/developer_tests/contrib/fortran-testanything/test.f90 @@ -0,0 +1,373 @@ +! Copyright 2015 Dennis Decker Jensen +! See and +! Tectonics: gfortran -g -Wall -Wextra -std=f2008ts -c test.f08 + +module test_base + use, intrinsic :: iso_fortran_env, only: output_unit, error_unit + implicit none + + ! Kept as variables instead of aliases, + ! so that test output or diagnostic output can be redirected + integer :: test_unit = output_unit, diag_unit = error_unit + + integer :: tests = 0, todos = 0 + character(len=120) :: todomsg = "" + + interface todo + module procedure todo_i, todo_s, todo_s_i, todo + end interface + +contains + + subroutine diag(msg) + character(len=*), intent(in) :: msg + write (diag_unit, '("# ",A)') trim(msg) ! only trailing spaces + end + + subroutine note(msg) + character(len=*), intent(in) :: msg + write (test_unit, '("# ",A)') trim(msg) + end + + subroutine testline(ok, msg, idmsg, gotmsg, expectedmsg) + logical, intent(in) :: ok + character(len=*), intent(in) :: msg, idmsg, gotmsg, expectedmsg + + tests = tests + 1 + if (.not. ok) call out("not ") + write (test_unit, '("ok ",I0)', advance="NO") tests + + if (msg /= "" .or. todos > 0) call out(" - ") + + if (msg /= "") call out(trim(msg)) + + if (todos > 0) then + todos = todos - 1 + if (msg /= "") call out(" ") + call out("# TODO") + if (todomsg .ne. "") then + call out(": ") + call out(trim(todomsg)) + end if + end if + if (todos == 0) todomsg = "" + + write (test_unit, *) "" + + if (.not. ok) then + ! 3 spaces prepended = 4 spaces indentation after # on diag + if (idmsg /= "") call diag(" " // idmsg) + if (gotmsg /= "") call diag(" " // gotmsg) + if (expectedmsg /= "") call diag(" " // expectedmsg) + end if + contains + subroutine out(str) + character(len=*), intent(in) :: str + write (test_unit, '(A)', advance="NO") str + end + end subroutine testline + + subroutine ok(condition, msg) + logical, intent(in) :: condition + character(len=*), intent(in), optional :: msg + if (present(msg)) then + call testline(condition, msg, "", "", "") + else + call testline(condition, "", "", "", "") + end if + end + + subroutine pass(msg) + character(len=*), intent(in), optional :: msg + call ok(.true., msg) + end + + subroutine fail(msg) + character(len=*), intent(in), optional :: msg + call ok(.false., msg) + end + + subroutine todo_s_i(msg, howmany) + character(len=*), intent(in) :: msg + integer, intent(in) :: howmany + todomsg = msg + todos = howmany + end + + subroutine todo + call todo_s_i("", 1) + end + + subroutine todo_s(msg) + character(len=*), intent(in) :: msg + call todo_s_i(msg, 1) + end + + subroutine todo_i(howmany) + integer, intent(in) :: howmany + call todo_s_i("", howmany) + end + +end module test_base + +module test_planning + use test_base, only: test_unit, tests + implicit none + + integer, private :: planned = 0 + +contains + + subroutine bail_out(msg) + character(len=*), intent(in), optional :: msg + if (present(msg)) then + write (test_unit, '("Bail out! ",A)') msg + else + write (test_unit, '("Bail out!")') + end if + stop + end + + subroutine plan(tests) + integer, intent(in) :: tests + + select case (tests) + case (:-1) + call bail_out("A plan with a negative number of tests") + case (0) + write (test_unit, '("1..0")') + stop ! The same as skip_all without a given reason + case (1:) + if (planned > 0) & + & call bail_out("More than one plan in test output") + planned = tests + write (test_unit, '("1..",I0)') planned + end select + end + + subroutine done_testing(howmany) + integer, intent(in), optional :: howmany + + ! Put plan at the end of test output + if (present(howmany)) then + call plan(howmany) + else + if (planned == 0) call plan(tests) + ! else - We already have a plan + end if + end + + subroutine skip_all(msg) + character(len=*), intent(in), optional :: msg + if (present(msg)) then + write (test_unit, '("1..0 # Skipped: ",A)') msg + else + write (test_unit, '("1..0 # Skipped all")') + end if + stop + end + +end module test_planning + +! Template instances of integer kinds for "is" + +module is_i8_mod + use, intrinsic :: iso_fortran_env, only: wp => int8 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_i.inc" +end + +module is_i16_mod + use, intrinsic :: iso_fortran_env, only: wp => int16 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_i.inc" +end + +module is_i32_mod + use, intrinsic :: iso_fortran_env, only: wp => int32 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_i.inc" +end + +module is_i64_mod + use, intrinsic :: iso_fortran_env, only: wp => int64 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_i.inc" +end + +module is_i + use is_i8_mod, only: is_i8 => is + use is_i16_mod, only: is_i16 => is + use is_i32_mod, only: is_i32 => is + use is_i64_mod, only: is_i64 => is + interface is + module procedure is_i8, is_i16, is_i32, is_i64 + end interface +end + +! Template instances of real kinds for "is" + +module is_r32_mod + use, intrinsic :: iso_fortran_env, only: wp => real32 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_r.inc" +end + +module is_r64_mod + use, intrinsic :: iso_fortran_env, only: wp => real64 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_r.inc" +end + +module is_r128_mod + use, intrinsic :: iso_fortran_env, only: wp => real128 + use, non_intrinsic :: test_base, only: testline, tests +contains + include "is_r.inc" +end + +module is_r + use is_r32_mod, only: isrel_r32 => isrel, isabs_r32 => isabs, & + & isnear_r32 => isnear + use is_r64_mod, only: isrel_r64 => isrel, isabs_r64 => isabs, & + & isnear_r64 => isnear + use is_r128_mod, only: isrel_r128 => isrel, isabs_r128 => isabs, & + & isnear_r128 => isnear + interface isrel + module procedure isrel_r32, isrel_r64, isrel_r128 + end interface + + interface isabs + module procedure isabs_r32, isabs_r64, isabs_r128 + end interface + + interface isnear + module procedure isnear_r32, isnear_r64, isnear_r128 + end interface +end + +module test_more + use test_base, only: testline, tests, test_unit + use test_planning, only: bail_out ! for negative skips + use is_i, only: is, is_i8, is_i16, is_i32, is_i64 + use is_r, only: isabs, isrel, isnear, & + & isabs_r32, isrel_r32, isnear_r32, & + & isabs_r64, isrel_r64, isnear_r64, & + & isabs_r128, isrel_r128, isnear_r128 + + ! Complex numbers cannot be compared, hence no is_c module + + implicit none + + interface skip + module procedure skip_i, skip_s, skip_s_i, skip + end interface + + interface is + module procedure is_s, is_l + end interface + +contains + + subroutine skip_s_i(msg, howmany) + character(len=*), intent(in) :: msg + integer, intent(in) :: howmany + character(len=120) skipmsg + integer i + + if (howmany <= 0) then + call bail_out("Skipped non-positive number of tests") + end if + + if (msg == "") then + skipmsg = "# SKIP" + else + skipmsg = "# SKIP: " // trim(msg) + end if + + do i = 1, howmany + tests = tests + 1 + write (test_unit, '("ok ",I0," ",A)') tests, trim(skipmsg) + end do + end + + subroutine skip + call skip_s_i("", 1) + end + + subroutine skip_s(msg) + character(len=*), intent(in) :: msg + call skip_s_i(msg, 1) + end + + subroutine skip_i(howmany) + integer, intent(in) :: howmany + call skip_s_i("", howmany) + end + + ! Duplicates of is_i routines in file is_i.inc and ditto is_r + ! They are not factored any further, because it is easier + ! to see all the output together rather than in separate routines + + subroutine is_s(got, expected, msg) + character(len=*), intent(in) :: got + character(len=*), intent(in) :: expected + character(len=*), intent(in), optional :: msg + character(len=:), allocatable :: testmsg, idmsg + character(len=120) gotmsg, expectedmsg + logical good + + if (present(msg)) then + allocate(character(len=len_trim(msg)+20) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,A,A)') 'Failed test: "', trim(msg), '"' + testmsg = trim(msg) + else + allocate(character(len=30) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,I0)') 'Failed test no. ', tests + 1 + testmsg = "" + end if + write (unit=gotmsg, fmt='(A,A,A)') ' got: "', got, '"' + write (unit=expectedmsg, fmt='(A,A,A)') 'expected: "', expected, '"' + + good = got == expected + call testline(good, testmsg, idmsg, gotmsg, expectedmsg) + end + + subroutine is_l(got, expected, msg) + logical, intent(in) :: got, expected + character(len=*), intent(in), optional :: msg + character(len=:), allocatable :: testmsg, idmsg + character(len=120) gotmsg, expectedmsg + logical good + + if (present(msg)) then + allocate(character(len=len_trim(msg)+20) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,A,A)') 'Failed test: "', trim(msg), '"' + testmsg = trim(msg) + else + allocate(character(len=30) :: testmsg, idmsg) + write (unit=idmsg, fmt='(A,I0)') 'Failed test no. ', tests + 1 + testmsg = "" + end if + write (unit=gotmsg, fmt='(A,L1)') ' got: ', got + write (unit=expectedmsg, fmt='(A,L1)') 'expected: ', expected + + good = got .eqv. expected + call testline(good, testmsg, idmsg, gotmsg, expectedmsg) + end + +end module test_more + +module test + use test_base, only: test_unit, diag_unit, & + & ok, diag, note, pass, fail, todo + use test_planning, only: plan, done_testing, skip_all, bail_out + use test_more, only: is, isabs, isrel, isnear, skip +end module test + diff --git a/developer_tests/window/test_window.f90 b/developer_tests/window/test_window.f90 new file mode 100644 index 0000000000..b910d8b69d --- /dev/null +++ b/developer_tests/window/test_window.f90 @@ -0,0 +1,50 @@ +program test_window + +use mpi_utilities_mod, only : initialize_mpi_utilities, finalize_mpi_utilities, my_task_id, task_count, task_sync + +use ensemble_manager_mod, only : init_ensemble_manager, end_ensemble_manager, ensemble_type, set_num_extra_copies +use distributed_state_mod, only : create_state_window, free_state_window, get_state +use types_mod, only : i8, r8 + +use test ! fortran-testanything + +implicit none + +integer :: num_copies = 10 +integer :: real_ens_members = 3 +real(r8) :: res(3) +integer(i8) :: num_vars = 201 +type(ensemble_type) :: ens_handle + +call initialize_mpi_utilities('test_window') + +if (my_task_id() == 0 ) then + call plan(3*task_count()) +endif + +call init_ensemble_manager(ens_handle, num_copies, num_vars) +call set_num_extra_copies(ens_handle, num_copies - real_ens_members) + +ens_handle%copies(1:real_ens_members,:) = my_task_id() +ens_handle%copies(real_ens_members+1:num_copies,:) = -100 + +call create_state_window(ens_handle) + +! result should be index-1 mod task_count() for round robin distribution +res = get_state(1_i8, ens_handle) +call ok(res(1) == mod(1-1, task_count())) + +res = get_state(27_i8, ens_handle) +call ok(res(1) == mod(27-1, task_count())) + +res = get_state(198_i8, ens_handle) +call ok(res(1) == mod(198-1, task_count())) + +call free_state_window(ens_handle) + +call end_ensemble_manager(ens_handle) + +call finalize_mpi_utilities() + + +end program test_window \ No newline at end of file diff --git a/developer_tests/window/work/input.nml b/developer_tests/window/work/input.nml new file mode 100644 index 0000000000..71e9d827f1 --- /dev/null +++ b/developer_tests/window/work/input.nml @@ -0,0 +1,15 @@ +&ensemble_manager_nml +/ + +&utilities_nml + module_details = .false. + / + +&preprocess_nml + input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' + output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' + input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' + output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' + obs_type_files = '../../../observations/forward_operators/obs_def_gps_mod.f90' + quantity_files = '../../../assimilation_code/modules/observations/default_quantities_mod.f90' + / diff --git a/developer_tests/window/work/quickbuild.sh b/developer_tests/window/work/quickbuild.sh new file mode 100755 index 0000000000..5482d0f19a --- /dev/null +++ b/developer_tests/window/work/quickbuild.sh @@ -0,0 +1,40 @@ +#!/usr/bin/env bash + +# DART software - Copyright UCAR. This open source software is provided +# by UCAR, "as is", without charge, subject to all terms of use at +# http://www.image.ucar.edu/DAReS/DART/DART_download + +main() { + + +export DART=$(git rev-parse --show-toplevel) +source "$DART"/build_templates/buildfunctions.sh + +MODEL="none" +EXTRA=$DART/models/template/threed_model_mod.f90 +dev_test=1 +LOCATION="threed_sphere" +TEST="window" + +programs=( +test_window +) + +# quickbuild arguments +arguments "$@" + +# clean the directory +\rm -f -- *.o *.mod Makefile .cppdefs + +# build and run preprocess before making any other DART executables +buildpreprocess + +# build DART +buildit + +# clean up +\rm -f -- *.o *.mod + +} + +main "$@" diff --git a/models/FESOM/model_mod.f90 b/models/FESOM/model_mod.f90 index 3233cf1492..b6387431f6 100644 --- a/models/FESOM/model_mod.f90 +++ b/models/FESOM/model_mod.f90 @@ -60,7 +60,7 @@ module model_mod use obs_kind_mod, only : get_index_for_quantity -use mpi_utilities_mod, only: my_task_id, broadcast_minmax, task_count +use mpi_utilities_mod, only: my_task_id, all_reduce_min_max, task_count use fesom_modules, only: read_node, read_aux3, read_depth, read_namelist, & nCells => myDim_nod2D, & ! number of surface locations @@ -828,7 +828,7 @@ subroutine pert_model_copies(ens_handle, ens_size, pert_amp, interf_provided) enddo ! get global min/max for each variable -call broadcast_minmax(min_var, max_var, num_variables) +call all_reduce_min_max(min_var, max_var, num_variables) deallocate(within_range) call init_random_seq(random_seq, my_task_id()+1) diff --git a/models/POP/readme.rst b/models/POP/readme.rst index 6ece91e7fa..3726297166 100644 --- a/models/POP/readme.rst +++ b/models/POP/readme.rst @@ -1,3 +1,5 @@ +.. _POP: + POP === @@ -83,9 +85,10 @@ Summary To use DART and CESM POP2 on NSF NCAR's supercomputer, you will need to complete the following steps. +#. Download an intial ensemble of POP2 restart files from the `NSF NCAR Geoscience + Data Exchange `_ #. Configure the scripts for your specific experiment by editing ``DART_params.csh``. -#. Stage your initial ensemble using ``copy_POP_JRA_restarts.py``. #. Run the appropriate DART setup script to create and build the CESM case. If the DART setup script runs to completion, it will print instructions to the @@ -114,53 +117,6 @@ in subdirectories that correspond releases of CESM. For example: contains scripts that should be used with CESM releases 2.1.0-2.1.3. -copy_POP_JRA_restarts.py -~~~~~~~~~~~~~~~~~~~~~~~~ - -This script stages an intial ensemble of POP2 restart files by copying files -from a prior experiment run by *Who Kim*. Thanks Who! - -These restart files can be used as an initial ensemble of model -states. The files are kept in a directory on GLADE that is owned by the Climate -and Global Dynamics (CGD) Ocean Section: - -.. code-block:: - - /glade/campaign/cgd/oce/people/whokim/csm/g210.G_JRA.v14.gx1v7.01 - -Unless you're already a member of the CGD Ocean Section, you must be granted -access to this directory by CISL. Use the `Service Desk -`_ to request permission. If -you're unable to get permission, contact DAReS staff for assistance by emailing -dart@ucar.edu. - -Filepaths beginning with ``/glade/campaign/*`` can't be accessed from NSF NCAR's -supercomputer nodes. You must log on to NSF NCAR's data visualization computer to -copy files from ``/glade/campaign/*``. - -This python script was created by *Dan Amrhein*. Thanks Dan! - -+-------------------------------+-----------------------------------------------------------+ -| Script name | Description | -+===============================+===========================================================+ -| ``copy_POP_JRA_restarts.py`` | This script copies restart files from the | -| | g210.G_JRA.v14.gx1v7.01 experiment that are saved in | -| | campaign storage. You must be granted access to the CGD | -| | Ocean Section campaign storage directory and be logged on | -| | to NSF NCAR's data visualization computer in order to run | -| | this script. The assignment of the ``stagedir`` variable | -| | in this script should match the assignment of the | -| | ``stagedir`` variable in ``DART_params.csh``. | -+-------------------------------+-----------------------------------------------------------+ - -In order to use this script, log in to NSF NCAR's data visualization computer and -use python to run the script. For example: - -.. code-block:: - - $ cd DART/models/POP/shell_scripts/cesm2_1 - $ python copy_POP_JRA_restarts.py - DART_params.csh ~~~~~~~~~~~~~~~ @@ -175,11 +131,18 @@ It is run by the setup scripts. | | that you need to set in order to build and run cases. You | | | must read this file carefully and configure the variables | | | to match your needs. The assignment of the ``stagedir`` | -| | variable in this script should match the assignment of | -| | the ``stagedir`` variable in | -| | ``copy_POP_JRA_restarts.py``. | +| | variable in this script should match the directory path | +| | where the restarts from the GDEX were downloaded/stored. | +---------------------+-----------------------------------------------------------+ +SourceMods are required to enable POP to recompute the barotropic velocity to prevent the +barotropic solver from crashing. You can find and download the SourceMods available +for POP-DART `here `_. + +Put the SourceMods in your home directory: + +~/${cesmtag}/SourceMods + Setup scripts ~~~~~~~~~~~~~ @@ -257,11 +220,10 @@ Karspeck et al. (2013) [3]_ find that an ensemble of 1 January model states selected from a multi-decade free-running integration of POP2 can be used as an initial ensemble. -If you have access to CGD's Ocean Section directory on ``/glade/campaign`` you -can use the `copy_POP_JRA_restarts.py`_ script to stage a collection of POP -restart files from Who Kim's multi-century ``g210.G_JRA.v14.gx1v7.01`` -experiment to serve as an initial ensemble. This experiment uses the JRA-55 -dataset for atmospheric forcing (Tsujino et al. 2018 [4]_). +You can access a collection of POP restart files from Who Kim's multi-century +``g210.G_JRA.v14.gx1v7.01`` experiment to serve as an initial ensemble in the +`NSF NCAR Geoscience Data Exchange `_. This +experiment uses the JRA-55 dataset for atmospheric forcing (Tsujino et al. 2018 [4]_). Observation sequence files ~~~~~~~~~~~~~~~~~~~~~~~~~~ diff --git a/models/POP/shell_scripts/cesm2_1/copy_POP_JRA_restarts.py b/models/POP/shell_scripts/cesm2_1/copy_POP_JRA_restarts.py deleted file mode 100644 index d0667833a2..0000000000 --- a/models/POP/shell_scripts/cesm2_1/copy_POP_JRA_restarts.py +++ /dev/null @@ -1,148 +0,0 @@ -from shutil import copyfile, copyfileobj -import getpass -import os -import sys -import glob -from pathlib import Path - -# -----------------------------------------------------# -# Define needed function # -# -----------------------------------------------------# - - -def unzipNsave(old, new): - ''' - Saves *.gz file in old as unzipped file in new - ''' - import gzip - - with gzip.open(old, 'rb') as f_in: - with open(new, 'wb') as f_out: - copyfileobj(f_in, f_out) - - -# -----------------------------------------------------# -# Main body of the script # -# -----------------------------------------------------# - -# Get username # -# -----------------------------------------------------# - -USER = getpass.getuser() - -# Create stage directory if it doesn't exist # -# -----------------------------------------------------# - -case_name = 'g210.G_JRA.v14.gx1v7.01' -stagedir = '/glade/scratch/' + USER + '/' + case_name + '/rest/2010-01-01-00000/' -Path("stagedir").mkdir(parents=True, exist_ok=True) - -# First load from the directory of MONTHLY saves # -# -----------------------------------------------------# - -load_path = '/glade/campaign/cgd/oce/people/whokim/csm/' + case_name + '/rest/' - -years = list(range(41, 202, 10))+list(range(209, 251)) -# For debugging: -# years = list(range(31,32)) -inds = list(range(1, len(years)+1)) - -print('Copying from ' + load_path + ' to ' + stagedir) - -for ii, year in enumerate(years): - - # Zero pad the year and index string so that there are 4 digits - zero_filled_year = str(year).zfill(4) - zero_filled_ii = str(inds[ii]).zfill(4) - - print('Copying year ' + str(year) + ' to ensemble number ' - + zero_filled_ii) - - yearpath = load_path + '/' + zero_filled_year + '-01-01-00000/' - - # Define file paths - popold = yearpath + case_name + '.pop.r.' + zero_filled_year \ - + '-01-01-00000.nc' - ovfold = yearpath + case_name + '.pop.ro.' + zero_filled_year \ - + '-01-01-00000' - ciceold = yearpath + case_name + '.cice.r.' + zero_filled_year \ - + '-01-01-00000.nc' - - popnew = stagedir + case_name + '.pop_' + zero_filled_ii \ - + '.r.2010-01-01-00000.nc' - ovfnew = stagedir + case_name + '.pop_' + zero_filled_ii \ - + '.ro.2010-01-01-00000' - cicenew = stagedir + case_name + '.cice_' + zero_filled_ii \ - + '.r.2010-01-01-00000.nc' - - # Copy over - copyfile(popold, popnew) - copyfile(ovfold, ovfnew) - copyfile(ciceold, cicenew) - - -# Then load the rest from DAILY saves -# Some years (251-262) must be decompressed -# -----------------------------------------------------# - -case_name_dr = 'g210.G_JRA.v14.gx1v7.01.dr' -load_path = '/glade/campaign/cgd/oce/people/whokim/csm/' + case_name_dr + '/rest' - -ly = len(years) -years2 = list(range(251, 272)) - -# Start from where we left off with the daily saves with inds -inds = list(range(ly+1, ly+len(years2)+1)) - -print('Copying from ' + load_path + ' to ' + stagedir) - -for ii, year in enumerate(years2): - - # Zero pad the year and index string so that there are 4 digits - zero_filled_year = str(year).zfill(4) - zero_filled_ii = str(inds[ii]).zfill(4) - - yearpath = load_path + '/' + zero_filled_year + '-01-01-00000/' - - # Define file paths - popold = glob.glob(yearpath + case_name_dr + '.pop.r.' + zero_filled_year - + '-01-01-00000.nc*')[0] - ovfold = glob.glob(yearpath + case_name_dr + '.pop.ro.' + zero_filled_year - + '-01-01-00000*')[0] - ciceold = glob.glob(yearpath + case_name_dr + '.cice.r.' + zero_filled_year - + '-01-01-00000.nc*')[0] - - popnew = stagedir + case_name + '.pop_' + zero_filled_ii \ - + '.r.2010-01-01-00000.nc' - ovfnew = stagedir + case_name + '.pop_' + zero_filled_ii \ - + '.ro.2010-01-01-00000' - cicenew = stagedir + case_name + '.cice_' + zero_filled_ii \ - + '.r.2010-01-01-00000.nc' - - # Copy over, decompressing as necessary - - if popold[-3:] == '.gz': - print('Decompressing and copying year ' + str(year) - + ' to ensemble number ' + zero_filled_ii) - else: - print('Copying year ' + str(year) + ' to ensemble number ' - + zero_filled_ii) - - if popold[-3:] == '.gz': - unzipNsave(popold, popnew) - else: - copyfile(popold, popnew) - - # OVF - if ovfold[-3:] == '.gz': - unzipNsave(ovfold, ovfnew) - else: - copyfile(ovfold, ovfnew) - - # CICE - if ciceold[-3:] == '.gz': - unzipNsave(ciceold, cicenew) - else: - copyfile(ciceold, cicenew) - -print('Completed') diff --git a/models/mpas_atm/model_mod.f90 b/models/mpas_atm/model_mod.f90 index 29282f7a8a..02e7361c70 100644 --- a/models/mpas_atm/model_mod.f90 +++ b/models/mpas_atm/model_mod.f90 @@ -104,7 +104,7 @@ module model_mod QTY_SURFACE_TYPE, & ! for rttov QTY_CLOUD_FRACTION ! for rttov -use mpi_utilities_mod, only: my_task_id, broadcast_minmax +use mpi_utilities_mod, only: my_task_id, all_reduce_min_max use random_seq_mod, only: random_seq_type, init_random_seq, random_gaussian @@ -1812,7 +1812,7 @@ subroutine pert_model_copies(ens_handle, ens_size, pert_amp, interf_provided) enddo ! get global min/max for each variable -call broadcast_minmax(min_var, max_var, num_variables) +call all_reduce_min_max(min_var, max_var, num_variables) deallocate(within_range) call init_random_seq(random_seq, my_task_id()+1) @@ -4686,8 +4686,8 @@ subroutine convert_vert_distrib(state_handle, ens_size, location, obs_kind, ztyp fdata = 0.0_r8 do i = 1, n where (istatus == 0) - fdata(i, :) = zGridFace(k_low(i, :),c(i))*(1.0_r8 - fract(i, :)) + & - zGridFace(k_up (i, :),c(i))*fract(i, :) + fdata(i, :) = zGridCenter(k_low(i, :),c(i))*(1.0_r8 - fract(i, :)) + & + zGridCenter(k_up (i, :),c(i))*fract(i, :) end where enddo @@ -4954,7 +4954,7 @@ subroutine convert_vert_distrib_state(state_handle, ens_size, location, quantity ! we have the vert_level and cellid - no need to call find_triangle or find_vert_indices zout(:) = vert_level - + istatus(:) = 0 if (debug > 9 .and. do_output()) then write(string2,'("zout_in_level for member 1:",F10.2)') zout(1) call error_handler(E_MSG, 'convert_vert_distrib_state',string2,source, revision, revdate) @@ -5013,10 +5013,12 @@ subroutine convert_vert_distrib_state(state_handle, ens_size, location, quantity ! of the quantities should use the level centers. if ( ndim == 1 ) then zout(:) = zGridFace(1, cellid) + istatus(:) = 0 else zout(:) = zGridCenter(vert_level, cellid) if ( quantity == QTY_VERTICAL_VELOCITY ) zout(:) = zGridFace(vert_level, cellid) if ( quantity == QTY_EDGE_NORMAL_SPEED ) zout(:) = zGridEdge(vert_level, cellid) + istatus(:) = 0 endif if (debug > 9 .and. do_output()) then @@ -5045,7 +5047,7 @@ subroutine convert_vert_distrib_state(state_handle, ens_size, location, quantity ! surf F, norm F: need fullp only ! surf F, norm T: need both surfp and fullp - at_surf = (ztypein == VERTISSURFACE) + at_surf = (ztypein == VERTISSURFACE) !HK ztypin is set to VERTISLEVEL before entering this case statement do_norm = .not. no_normalization_of_scale_heights ! if normalizing pressure and we're on the surface, by definition scale height