Skip to content

Commit

Permalink
Merge branch 'main' into pop_docs
Browse files Browse the repository at this point in the history
  • Loading branch information
@hkershaw-brown
hkershaw-brown authored Nov 8, 2024
2 parents 9159690 + 89c6fc1 commit 6631c47
Show file tree
Hide file tree
Showing 7 changed files with 50 additions and 181 deletions.
45 changes: 0 additions & 45 deletions assimilation_code/modules/assimilation/normal_distribution_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
31 changes: 5 additions & 26 deletions assimilation_code/modules/utilities/mpi_utilities_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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'

Expand Down Expand Up @@ -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)
Expand All @@ -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
Expand Down
31 changes: 5 additions & 26 deletions assimilation_code/modules/utilities/mpif08_utilities_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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'

Expand Down Expand Up @@ -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)
Expand All @@ -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
Expand Down
24 changes: 4 additions & 20 deletions assimilation_code/modules/utilities/null_mpi_utilities_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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'

Expand Down Expand Up @@ -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

!-----------------------------------------------------------------------------

Expand Down
92 changes: 32 additions & 60 deletions assimilation_code/modules/utilities/time_manager_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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.

Expand All @@ -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
Expand All @@ -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

Expand All @@ -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
Expand Down Expand Up @@ -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
Expand Down
4 changes: 2 additions & 2 deletions models/FESOM/model_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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)
Expand Down
Loading

0 comments on commit 6631c47

Please sign in to comment.