Skip to content

Commit

Permalink
Split output_for_GCM back to two arrays
Browse files Browse the repository at this point in the history
Put surface_flux_output back on the interface, and also created
interior_tendency_output. This change lets us drop field_source from the output
type, and add_output_for_GCM() does all the error checking I had introduced
into marbl_single_output_constructor (making sure the variable is registered
and available).

Note that add_output_for_GCM() now returns field_source as well, so the GCM can
keep the surface flux output and interior tendency output separate when copying
to local memory (without knowing what fields come from surface_flux_compute()
vs interior_tendency_compute())
  • Loading branch information
mnlevy1981 committed Feb 15, 2024
1 parent d572db6 commit a33bf57
Show file tree
Hide file tree
Showing 4 changed files with 173 additions and 128 deletions.
69 changes: 54 additions & 15 deletions src/marbl_interface.F90
Original file line number Diff line number Diff line change
Expand Up @@ -92,16 +92,17 @@ module marbl_interface
real (r8), allocatable , public :: interior_tendencies(:,:) ! output
type(marbl_interior_tendency_forcing_indexing_type), public :: interior_tendency_forcing_ind ! FIXME #311: should be private
type(marbl_diagnostics_type) , public :: interior_tendency_diags ! output
type(marbl_output_for_GCM_type) , public :: interior_tendency_output ! output

! public data related to computing surface fluxes
real (r8) , public, allocatable :: tracers_at_surface(:,:) ! input
type(marbl_forcing_fields_type) , public, allocatable :: surface_flux_forcings(:) ! input
type(marbl_surface_flux_forcing_indexing_type) , public :: surface_flux_forcing_ind ! FIXME #311: should be private
real (r8) , public, allocatable :: surface_fluxes(:,:) ! output
type(marbl_diagnostics_type) , public :: surface_flux_diags ! output
type(marbl_output_for_GCM_type) , public :: surface_flux_output ! output

! public data that the GCM needs to explicitly request
type(marbl_output_for_GCM_type), public :: output_for_gcm ! output

! public data - global averages
real (r8), public, allocatable :: glo_avg_fields_interior_tendency(:) ! output (nfields)
Expand Down Expand Up @@ -775,30 +776,68 @@ end function get_settings_var_cnt

!***********************************************************************

subroutine add_output_for_GCM(this, num_elements, field_name, output_id, num_levels)
! Currently, we only need this subroutine for surface fluxes.
! If we introduce this%interior_tendency_output then this function will need
! a field_source argument (either 'surface_flux' or 'interior_tendency')
subroutine add_output_for_GCM(this, num_elements, field_name, output_id, field_source, num_levels)
! Check the registry to see if field_name is provided from surface_flux_compute()
! or interior_tendency_compute(); add it to the proper output_for_GCM type, or
! return a useful error message

class (marbl_interface_class), intent(inout) :: this
character(len=*), intent(in) :: field_name
integer(int_kind), intent(in) :: num_elements
integer(int_kind), intent(out) :: output_id
character(len=*), intent(out) :: field_source
integer(int_kind), optional, intent(in) :: num_levels

character(len=*), parameter :: subname = 'marbl_interface:add_output_for_GCM'
character(len=char_len) :: log_message
integer :: m

call this%output_for_gcm%add_output(this%output_for_gcm_registry, &
num_elements, &
field_name, &
output_id, &
this%StatusLog, &
num_levels)
if (this%StatusLog%labort_marbl) then
call this%StatusLog%log_error_trace('output_for_gcm%add_output()', subname)
output_id = 0
field_source = ""

do m=1,size(this%output_for_gcm_registry%registered_outputs)
if (trim(field_name) == trim(this%output_for_gcm_registry%registered_outputs(m)%short_name)) then
! err_message will be populated if this field is unavailable in current configuration
if (len_trim(this%output_for_gcm_registry%registered_outputs(m)%err_message) > 0) then
call this%StatusLog%log_error(this%output_for_gcm_registry%registered_outputs(m)%err_message, subname)
return
end if
exit
end if
end do

! Abort if field_name was not registered
if (m > size(this%output_for_gcm_registry%registered_outputs)) then
write(log_message, "(2A)") trim(field_name), " is not a valid output field name for the GCM"
call this%StatusLog%log_error(log_message, subname)
return
end if

write(log_message, "(3A)") "Adding ", trim(field_name), " to outputs needed by the GCM"
call this%StatusLog%log_noerror(log_message, subname)

! Set field source, and then add output to appropriate output_for_GCM_type
field_source = trim(this%output_for_gcm_registry%registered_outputs(m)%field_source)
if ( trim(field_source) == "surface_flux") then
call this%surface_flux_output%add_output(this%output_for_gcm_registry%registered_outputs(m)%short_name, &
this%output_for_gcm_registry%registered_outputs(m)%long_name, &
this%output_for_gcm_registry%registered_outputs(m)%units, &
num_elements, &
output_id, &
num_levels)
end if
if (trim(field_source) == "interior_tendency") then
call this%interior_tendency_output%add_output(this%output_for_gcm_registry%registered_outputs(m)%short_name, &
this%output_for_gcm_registry%registered_outputs(m)%long_name, &
this%output_for_gcm_registry%registered_outputs(m)%units, &
num_elements, &
output_id, &
num_levels)
end if

! %id is a pointer to a member of either sfo_ind or ito_ind
this%output_for_gcm_registry%registered_outputs(m)%id = output_id

end subroutine add_output_for_GCM

!***********************************************************************
Expand Down Expand Up @@ -958,7 +997,7 @@ subroutine interior_tendency_compute(this)
zooplankton_local = this%zooplankton_local, &
zooplankton_share = this%zooplankton_share, &
saved_state = this%interior_tendency_saved_state, &
output_for_gcm = this%output_for_gcm, &
output_for_gcm = this%interior_tendency_output, &
marbl_timers = this%timers, &
interior_tendency_share = this%interior_tendency_share, &
marbl_particulate_share = this%particulate_share, &
Expand Down Expand Up @@ -1008,7 +1047,7 @@ subroutine surface_flux_compute(this)
marbl_tracer_indices = this%tracer_indices, &
saved_state = this%surface_flux_saved_state, &
saved_state_ind = this%surf_state_ind, &
output_for_gcm = this%output_for_gcm, &
output_for_gcm = this%surface_flux_output, &
surface_flux_internal = this%surface_flux_internal, &
surface_flux_share = this%surface_flux_share, &
surface_flux_diags = this%surface_flux_diags, &
Expand Down
78 changes: 20 additions & 58 deletions src/marbl_interface_public_types.F90
Original file line number Diff line number Diff line change
Expand Up @@ -150,7 +150,6 @@ module marbl_interface_public_types
character(len=char_len) :: long_name
character(len=char_len) :: short_name
character(len=char_len) :: units
character(len=char_len) :: field_source
real(r8), allocatable, dimension(:) :: forcing_field_0d
real(r8), allocatable, dimension(:,:) :: forcing_field_1d
contains
Expand Down Expand Up @@ -462,62 +461,32 @@ end subroutine marbl_single_diag_init

!*****************************************************************************

subroutine marbl_single_output_constructor(this, output_registry, num_elements, num_levels, field_name, id, marbl_status_log)
subroutine marbl_single_output_constructor(this, short_name, long_name, units, num_elements, num_levels)

class(marbl_single_output_type), intent(out) :: this
type(marbl_output_for_GCM_registry_type), intent(in) :: output_registry
character(len=*), intent(in) :: short_name
character(len=*), intent(in) :: long_name
character(len=*), intent(in) :: units
integer(int_kind), intent(in) :: num_elements
integer(int_kind), intent(in) :: num_levels
character(len=*), intent(in) :: field_name
integer(int_kind), intent(in) :: id
type(marbl_log_type), intent(inout) :: marbl_status_log

character(len=*), parameter :: subname = 'marbl_interface_public_types:marbl_single_output_constructor'
character(len=char_len) :: log_message
integer :: m

do m=1,size(output_registry%registered_outputs)
if (trim(field_name) == trim(output_registry%registered_outputs(m)%short_name)) then
! err_message will be populated if this field is unavailable in current configuration
if (len_trim(output_registry%registered_outputs(m)%err_message) > 0) then
call marbl_status_log%log_error(output_registry%registered_outputs(m)%err_message, subname)
return
end if

write(log_message, "(3A)") "Adding ", trim(field_name), " to outputs needed by the GCM"
call marbl_status_log%log_noerror(log_message, subname)

this%short_name = output_registry%registered_outputs(m)%short_name
this%long_name = output_registry%registered_outputs(m)%long_name
this%units = output_registry%registered_outputs(m)%units
this%field_source = output_registry%registered_outputs(m)%field_source

if (num_levels .eq. 0) then
allocate(this%forcing_field_0d(num_elements))
this%forcing_field_0d = c0
else
allocate(this%forcing_field_1d(num_elements, num_levels))
this%forcing_field_1d = c0
end if

! Set ofg_ind index for field_name (via pointer)
output_registry%registered_outputs(m)%id = id
exit
end if
end do
this%short_name = short_name
this%long_name = long_name
this%units = units

! Abort if field_name was not registered
if (m > size(output_registry%registered_outputs)) then
write(log_message, "(2A)") trim(field_name), " is not a valid output field name for the GCM"
call marbl_status_log%log_error(log_message, subname)
return
if (num_levels .eq. 0) then
allocate(this%forcing_field_0d(num_elements))
this%forcing_field_0d = c0
else
allocate(this%forcing_field_1d(num_elements, num_levels))
this%forcing_field_1d = c0
end if

end subroutine marbl_single_output_constructor

!*****************************************************************************

subroutine marbl_output_add(this, output_registry, num_elements, field_name, output_id, marbl_status_log, num_levels)
subroutine marbl_output_add(this, short_name, long_name, units, num_elements, output_id, num_levels)

! MARBL uses pointers to create an extensible allocatable array. The output
! fields (part of the intent(out) of this routine) are stored in
Expand All @@ -526,23 +495,21 @@ subroutine marbl_output_add(this, output_registry, num_elements, field_name, out
!
! 1) allocate new_output to be size N (one element larger than this%outputs_for_GCM)
! 2) copy this%outputs_for_GCM into first N-1 elements of new_output
! 3) newest surface flux output (field_name) is Nth element of new_output
! 3) newest surface flux output or interior tendency output (short_name) is Nth element of new_output
! 4) deallocate / nullify this%outputs_for_GCM
! 5) point this%outputs_for_GCM => new_output
!
! If the number of possible surface flux output fields grows, this workflow
! may need to be replaced with something that is not O(N^2).

class(marbl_output_for_GCM_type), intent(inout) :: this
type(marbl_output_for_GCM_registry_type), intent(in) :: output_registry
character(len=*), intent(in) :: short_name
character(len=*), intent(in) :: long_name
character(len=*), intent(in) :: units
integer(int_kind), intent(in) :: num_elements
character(len=*), intent(in) :: field_name
integer(int_kind), intent(out) :: output_id
type(marbl_log_type), intent(inout) :: marbl_status_log
integer(int_kind), optional, intent(in) :: num_levels

character(len=*), parameter :: subname = 'marbl_interface_public_types:marbl_output_add'

type(marbl_single_output_type), dimension(:), pointer :: new_output
integer :: n, old_size, dim1_loc, dim2_loc, num_levels_loc

Expand All @@ -568,7 +535,6 @@ subroutine marbl_output_add(this, output_registry, num_elements, field_name, out
new_output(n)%long_name = this%outputs_for_GCM(n)%long_name
new_output(n)%short_name = this%outputs_for_GCM(n)%short_name
new_output(n)%units = this%outputs_for_GCM(n)%units
new_output(n)%field_source = this%outputs_for_GCM(n)%field_source
if (allocated(this%outputs_for_GCM(n)%forcing_field_0d)) then
dim1_loc = size(this%outputs_for_GCM(n)%forcing_field_0d)
allocate(new_output(n)%forcing_field_0d(dim1_loc))
Expand All @@ -585,11 +551,7 @@ subroutine marbl_output_add(this, output_registry, num_elements, field_name, out
end do

! 3) newest surface flux output (field_name) is Nth element of new_output
call new_output(output_id)%construct(output_registry, num_elements, num_levels_loc, field_name, output_id, marbl_status_log)
if (marbl_status_log%labort_marbl) then
call marbl_status_log%log_error_trace('new_output%construct()', subname)
return
end if
call new_output(output_id)%construct(short_name, long_name, units, num_elements, num_levels_loc)

! 4) deallocate / nullify this%outputs_for_GCM
if (old_size .gt. 0) then
Expand Down Expand Up @@ -857,7 +819,7 @@ subroutine create_registry(this, base_bio_on, conc_flux_units)
allocate(ofg_ind%flux_o2_id, source=0)
this%registered_outputs(ofg_ind_loc)%id => ofg_ind%flux_o2_id
if (.not. (base_bio_on .and. lflux_gas_o2)) &
write(this%registered_outputs(ofg_ind_loc)%err_message, "(A,1X,A)") "Can not add flux_co2 to outputs without", &
write(this%registered_outputs(ofg_ind_loc)%err_message, "(A,1X,A)") "Can not add flux_o2 to outputs without", &
"base biotic tracers and lflux_gas_o2"

ofg_ind_loc = ofg_ind_loc + 1
Expand Down
Loading

0 comments on commit a33bf57

Please sign in to comment.