From 19585910eaae1a596dc4fea6a20be723281e864a Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Tue, 17 Dec 2024 17:29:20 -0700 Subject: [PATCH] Implement support for log levels in MPAS dynamical core * The following log levels are supported: * log_level_quiet (a.k.a. debugout_none in CAM-SIMA) * log_level_info (a.k.a. debugout_info in CAM-SIMA) * log_level_verbose (a.k.a. debugout_verbose in CAM-SIMA) * log_level_debug (a.k.a. debugout_debug in CAM-SIMA) * Log messages that originate from CAM-SIMA side are prefixed by "MPAS Interface (N): ", where N is the MPI rank * Log messages that originate from MPAS side are prefixed by "MPAS Subdriver (N): ", where N is the MPI rank * Log messages have been improved to be more descriptive and user-friendly --- .../mpas/driver/dyn_mpas_subdriver.F90 | 304 +++++++++++------- src/dynamics/mpas/dyn_comp.F90 | 177 +++++++--- src/dynamics/mpas/dyn_coupling.F90 | 38 ++- src/dynamics/mpas/dyn_grid.F90 | 76 +++-- 4 files changed, 398 insertions(+), 197 deletions(-) diff --git a/src/dynamics/mpas/driver/dyn_mpas_subdriver.F90 b/src/dynamics/mpas/driver/dyn_mpas_subdriver.F90 index ea793dcf..63ddc6f0 100644 --- a/src/dynamics/mpas/driver/dyn_mpas_subdriver.F90 +++ b/src/dynamics/mpas/driver/dyn_mpas_subdriver.F90 @@ -34,6 +34,19 @@ subroutine model_error_if(message, file, line) end subroutine model_error_if end interface + ! The supported log levels of MPAS dynamical core. + + !> Log nothing. + integer, parameter :: log_level_quiet = 0 + !> Log plain and user-friendly information about the status of MPAS dynamical core. + !> Public procedures should start with this log level. + integer, parameter :: log_level_info = 1 + !> Same as the above, but for private procedures. + integer, parameter :: log_level_verbose = 2 + !> Log elaborate information about the inner workings of MPAS dynamical core, which may be useful for diagnosing issues. + !> However, the log volume may be very large. + integer, parameter :: log_level_debug = 3 + !> The native floating-point precision of MPAS dynamical core. integer, parameter :: mpas_dynamical_core_real_kind = rkind @@ -43,9 +56,8 @@ end subroutine model_error_if type :: mpas_dynamical_core_type private - logical, public :: debug_output = .false. - ! Initialized by `dyn_mpas_init_phase1`. + integer :: log_level = log_level_quiet integer :: log_unit = output_unit #ifdef MPAS_USE_MPI_F08 type(mpi_comm_type) :: mpi_comm = mpi_comm_null @@ -292,17 +304,17 @@ end subroutine model_error_if var_info_type('vorticity' , 'real' , 2) & ] contains - !> Print a debug message with optionally the value(s) of a variable. + !> Print a debug message at a debug level. !> If `printer` is not supplied, the MPI root rank will print. Otherwise, the designated MPI rank will print instead. !> (KCW, 2024-02-03) - subroutine dyn_mpas_debug_print(self, message, variable, printer) + subroutine dyn_mpas_debug_print(self, level, message, printer) class(mpas_dynamical_core_type), intent(in) :: self + integer, intent(in) :: level character(*), intent(in) :: message - class(*), optional, intent(in) :: variable(:) integer, optional, intent(in) :: printer - ! Bail out early if debug output is not requested. - if (.not. self % debug_output) then + ! Bail out early if the log level is less verbose than the debug level. + if (self % log_level < level) then return end if @@ -316,13 +328,7 @@ subroutine dyn_mpas_debug_print(self, message, variable, printer) end if end if - if (present(variable)) then - write(self % log_unit, '(a)') 'dyn_mpas_debug_print (' // stringify([self % mpi_rank]) // '): ' // & - message // stringify(variable) - else - write(self % log_unit, '(a)') 'dyn_mpas_debug_print (' // stringify([self % mpi_rank]) // '): ' // & - message - end if + write(self % log_unit, '(a)') 'MPAS Subdriver (' // stringify([self % mpi_rank]) // '): ' // message end subroutine dyn_mpas_debug_print !> Convert one or more values of any intrinsic data types to a character string for pretty printing. @@ -439,7 +445,7 @@ end function stringify !> Ported and refactored for CAM-SIMA. (KCW, 2024-02-02) ! !------------------------------------------------------------------------------- - subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas_log_unit) + subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_level, log_unit, mpas_log_unit) ! Module(s) from MPAS. use atm_core_interface, only: atm_setup_core, atm_setup_domain use mpas_domain_routines, only: mpas_allocate_domain @@ -452,6 +458,7 @@ subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas integer, intent(in) :: mpi_comm #endif procedure(model_error_if) :: model_error_impl + integer, intent(in) :: log_level integer, intent(in) :: log_unit integer, intent(in) :: mpas_log_unit(2) @@ -472,11 +479,12 @@ subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas end if self % mpi_rank_root = (self % mpi_rank == 0) + self % log_level = max(min(log_level, log_level_debug), log_level_quiet) self % log_unit = log_unit - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') - call self % debug_print('Allocating core') + call self % debug_print(log_level_info, 'Allocating core') allocate(self % corelist, stat=ierr) @@ -486,7 +494,7 @@ subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas nullify(self % corelist % next) - call self % debug_print('Allocating domain') + call self % debug_print(log_level_info, 'Allocating domain') allocate(self % corelist % domainlist, stat=ierr) @@ -503,20 +511,20 @@ subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas self % domain_ptr % domainid = 0 - call self % debug_print('Calling mpas_framework_init_phase1') + call self % debug_print(log_level_info, 'Initializing MPAS framework (Phase 1/2)') - ! Initialize MPAS framework with supplied MPI communicator group. + ! Initialize MPAS framework with the supplied MPI communicator group. call mpas_framework_init_phase1(self % domain_ptr % dminfo, external_comm=self % mpi_comm) - call self % debug_print('Setting up core') + call self % debug_print(log_level_info, 'Setting up core') call atm_setup_core(self % corelist) - call self % debug_print('Setting up domain') + call self % debug_print(log_level_info, 'Setting up domain') call atm_setup_domain(self % domain_ptr) - call self % debug_print('Setting up log') + call self % debug_print(log_level_info, 'Setting up log') ! Set up the log manager as early as possible so we can use it for any errors/messages during subsequent ! initialization steps. @@ -528,11 +536,12 @@ subroutine dyn_mpas_init_phase1(self, mpi_comm, model_error_impl, log_unit, mpas ierr = self % domain_ptr % core % setup_log(self % domain_ptr % loginfo, self % domain_ptr, unitnumbers=mpas_log_unit) if (ierr /= 0) then - call self % model_error('Failed to setup log for MPAS', subname, __LINE__) + call self % model_error('Log setup failed for core ' // trim(self % domain_ptr % core % corename), & + subname, __LINE__) end if ! At this point, we should be ready to read namelist in `dyn_comp::dyn_readnl`. - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_init_phase1 !------------------------------------------------------------------------------- @@ -562,12 +571,12 @@ subroutine dyn_mpas_read_namelist(self, namelist_path, & integer :: ierr logical, pointer :: config_pointer_l - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(config_pointer_c) nullify(config_pointer_l) - call self % debug_print('Reading namelist at ', [namelist_path]) + call self % debug_print(log_level_info, 'Reading namelist at "' // trim(adjustl(namelist_path)) // '"') ! Override namelist filename so that we can rely on upstream MPAS functionality for reading its own namelist. ! The case of missing namelist groups (i.e., `iostat == iostat_end` or `iostat == iostat_eor`) will be handled gracefully. @@ -585,7 +594,7 @@ subroutine dyn_mpas_read_namelist(self, namelist_path, & ! Override designated namelist variables according to information provided from CAM-SIMA. ! These include runtime settings that cannot be determined beforehand. - call self % debug_print('Overriding designated namelist variables') + call self % debug_print(log_level_info, 'Overriding designated namelist variables') ! CAM-SIMA seems to follow "NetCDF Climate and Forecast (CF) Metadata Conventions" for calendar names. See ! CF-1.11, section "4.4.1. Calendar". @@ -606,7 +615,7 @@ subroutine dyn_mpas_read_namelist(self, namelist_path, & call self % get_variable_pointer(config_pointer_c, 'cfg', 'config_calendar_type') config_pointer_c = trim(adjustl(mpas_calendar)) - call self % debug_print('config_calendar_type = ', [config_pointer_c]) + call self % debug_print(log_level_debug, 'config_calendar_type = ' // stringify([config_pointer_c])) nullify(config_pointer_c) ! MPAS represents date and time in ISO 8601 format. However, the separator between date and time is `_` @@ -615,20 +624,20 @@ subroutine dyn_mpas_read_namelist(self, namelist_path, & call self % get_variable_pointer(config_pointer_c, 'cfg', 'config_start_time') config_pointer_c = stringify(start_date_time(1:3), '-') // '_' // stringify(start_date_time(4:6), ':') - call self % debug_print('config_start_time = ', [config_pointer_c]) + call self % debug_print(log_level_debug, 'config_start_time = ' // stringify([config_pointer_c])) nullify(config_pointer_c) call self % get_variable_pointer(config_pointer_c, 'cfg', 'config_stop_time') config_pointer_c = stringify(stop_date_time(1:3), '-') // '_' // stringify(stop_date_time(4:6), ':') - call self % debug_print('config_stop_time = ', [config_pointer_c]) + call self % debug_print(log_level_debug, 'config_stop_time = ' // stringify([config_pointer_c])) nullify(config_pointer_c) ! Format in `DD_hh:mm:ss` is acceptable. call self % get_variable_pointer(config_pointer_c, 'cfg', 'config_run_duration') config_pointer_c = stringify([run_duration(1)]) // '_' // stringify(run_duration(2:4), ':') - call self % debug_print('config_run_duration = ', [config_pointer_c]) + call self % debug_print(log_level_debug, 'config_run_duration = ' // stringify([config_pointer_c])) nullify(config_pointer_c) ! Reflect current run type to MPAS. @@ -642,10 +651,10 @@ subroutine dyn_mpas_read_namelist(self, namelist_path, & config_pointer_l = .true. end if - call self % debug_print('config_do_restart = ', [config_pointer_l]) + call self % debug_print(log_level_debug, 'config_do_restart = ' // stringify([config_pointer_l])) nullify(config_pointer_l) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_read_namelist !------------------------------------------------------------------------------- @@ -676,9 +685,9 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) integer :: ierr logical :: pio_iosystem_active - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') - call self % debug_print('Checking PIO system descriptor') + call self % debug_print(log_level_info, 'Checking PIO system descriptor') if (.not. associated(pio_iosystem)) then call self % model_error('Invalid PIO system descriptor', subname, __LINE__) @@ -690,9 +699,9 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) call self % model_error('Invalid PIO system descriptor', subname, __LINE__) end if - call self % debug_print('Calling mpas_framework_init_phase2') + call self % debug_print(log_level_info, 'Initializing MPAS framework (Phase 2/2)') - ! Initialize MPAS framework with supplied PIO system descriptor. + ! Initialize MPAS framework with the supplied PIO system descriptor. call mpas_framework_init_phase2(self % domain_ptr, io_system=pio_iosystem) ! Instantiate `streaminfo` but do not actually initialize it. Any queries made to it will always return `.false.`. @@ -704,6 +713,8 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) subname, __LINE__) end if + call self % debug_print(log_level_info, 'Defining packages') + ierr = self % domain_ptr % core % define_packages(self % domain_ptr % packages) if (ierr /= 0) then @@ -711,6 +722,8 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) subname, __LINE__) end if + call self % debug_print(log_level_info, 'Setting up packages') + ierr = self % domain_ptr % core % setup_packages( & self % domain_ptr % configs, self % domain_ptr % streaminfo, & self % domain_ptr % packages, self % domain_ptr % iocontext) @@ -720,6 +733,8 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) subname, __LINE__) end if + call self % debug_print(log_level_info, 'Setting up decompositions') + ierr = self % domain_ptr % core % setup_decompositions(self % domain_ptr % decompositions) if (ierr /= 0) then @@ -727,6 +742,8 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) subname, __LINE__) end if + call self % debug_print(log_level_info, 'Setting up clock') + ierr = self % domain_ptr % core % setup_clock(self % domain_ptr % clock, self % domain_ptr % configs) if (ierr /= 0) then @@ -736,7 +753,7 @@ subroutine dyn_mpas_init_phase2(self, pio_iosystem) ! At this point, we should be ready to set up decompositions, build halos, allocate blocks, etc. ! in `dyn_grid::model_grid_init`. - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_init_phase2 !------------------------------------------------------------------------------- @@ -775,7 +792,7 @@ subroutine dyn_mpas_init_phase3(self, number_of_constituents, pio_file) integer, pointer :: num_scalars type(mpas_pool_type), pointer :: mpas_pool - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(mpas_pool) nullify(num_scalars) @@ -785,7 +802,7 @@ subroutine dyn_mpas_init_phase3(self, number_of_constituents, pio_file) ! (i.e., segmentation fault due to invalid memory access) if `qv` is not allocated. self % number_of_constituents = max(1, number_of_constituents) - call self % debug_print('Number of constituents is ', [self % number_of_constituents]) + call self % debug_print(log_level_info, 'Number of constituents is ' // stringify([self % number_of_constituents])) ! Adding a config named `cam_pcnst` with the number of constituents will indicate to MPAS that ! it is operating as a dynamical core, and therefore it needs to allocate scalars separately @@ -797,7 +814,7 @@ subroutine dyn_mpas_init_phase3(self, number_of_constituents, pio_file) mesh_filename = 'external mesh' mesh_format = mpas_io_pnetcdf - call self % debug_print('Checking PIO file descriptor') + call self % debug_print(log_level_info, 'Checking PIO file descriptor') if (.not. associated(pio_file)) then call self % model_error('Invalid PIO file descriptor', subname, __LINE__) @@ -807,12 +824,12 @@ subroutine dyn_mpas_init_phase3(self, number_of_constituents, pio_file) call self % model_error('Invalid PIO file descriptor', subname, __LINE__) end if - call self % debug_print('Calling mpas_bootstrap_framework_phase1') + call self % debug_print(log_level_info, 'Bootstrapping MPAS framework (Phase 1/2)') ! Finish setting up blocks. call mpas_bootstrap_framework_phase1(self % domain_ptr, mesh_filename, mesh_format, pio_file_desc=pio_file) - call self % debug_print('Calling mpas_bootstrap_framework_phase2') + call self % debug_print(log_level_info, 'Bootstrapping MPAS framework (Phase 2/2)') ! Finish setting up fields. call mpas_bootstrap_framework_phase2(self % domain_ptr, pio_file_desc=pio_file) @@ -837,7 +854,7 @@ subroutine dyn_mpas_init_phase3(self, number_of_constituents, pio_file) nullify(mpas_pool) nullify(num_scalars) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_init_phase3 !------------------------------------------------------------------------------- @@ -885,7 +902,7 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) type(field3dreal), pointer :: field_3d_real type(mpas_pool_type), pointer :: mpas_pool - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(field_3d_real) nullify(mpas_pool) @@ -971,6 +988,8 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) ! Create index mapping between MPAS scalars and constituent names. For example, ! MPAS scalar index `i` corresponds to constituent index `index_mpas_scalar_to_constituent(i)`. + call self % debug_print(log_level_info, 'Creating index mapping between MPAS scalars and constituents') + allocate(self % index_mpas_scalar_to_constituent(self % number_of_constituents), stat=ierr) if (ierr /= 0) then @@ -1002,6 +1021,8 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) ! Create inverse index mapping between MPAS scalars and constituent names. For example, ! Constituent index `i` corresponds to MPAS scalar index `index_constituent_to_mpas_scalar(i)`. + call self % debug_print(log_level_info, 'Creating inverse index mapping between MPAS scalars and constituents') + allocate(self % index_constituent_to_mpas_scalar(self % number_of_constituents), stat=ierr) if (ierr /= 0) then @@ -1019,17 +1040,19 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) ! Print information about constituents. do i = 1, self % number_of_constituents - call self % debug_print('Constituent index ' // stringify([i])) - call self % debug_print(' Constituent name: ' // & + call self % debug_print(log_level_info, 'Constituent index ' // stringify([i])) + call self % debug_print(log_level_info, ' Constituent name: ' // & stringify([self % constituent_name(i)])) - call self % debug_print(' Is water species: ' // & + call self % debug_print(log_level_info, ' Is water species: ' // & stringify([self % is_water_species(i)])) - call self % debug_print(' Index mapping from constituent to MPAS scalar: ' // & + call self % debug_print(log_level_info, ' Index mapping from constituent to MPAS scalar: ' // & stringify([i]) // ' -> ' // stringify([self % index_constituent_to_mpas_scalar(i)])) end do ! Define "scalars" for MPAS. + call self % debug_print(log_level_info, 'Defining MPAS scalars') + call self % get_pool_pointer(mpas_pool, 'state') call mpas_pool_add_dimension(mpas_pool, 'index_qv', index_qv) @@ -1052,12 +1075,12 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) ! Print information about MPAS scalars. Only do it once. if (i == 1) then - call self % debug_print('MPAS scalar index ' // stringify([j])) - call self % debug_print(' MPAS scalar name: ' // & + call self % debug_print(log_level_info, 'MPAS scalar index ' // stringify([j])) + call self % debug_print(log_level_info, ' MPAS scalar name: ' // & stringify([field_3d_real % constituentnames(j)])) - call self % debug_print(' Is water species: ' // & + call self % debug_print(log_level_info, ' Is water species: ' // & stringify([self % is_water_species(self % index_mpas_scalar_to_constituent(j))])) - call self % debug_print(' Index mapping from MPAS scalar to constituent: ' // & + call self % debug_print(log_level_info, ' Index mapping from MPAS scalar to constituent: ' // & stringify([j]) // ' -> ' // stringify([self % index_mpas_scalar_to_constituent(j)])) end if end do @@ -1069,6 +1092,8 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) ! Define "scalars_tend" for MPAS. + call self % debug_print(log_level_info, 'Defining MPAS scalar tendencies') + call self % get_pool_pointer(mpas_pool, 'tend') call mpas_pool_add_dimension(mpas_pool, 'index_qv', index_qv) @@ -1088,6 +1113,17 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) do j = 1, self % number_of_constituents field_3d_real % constituentnames(j) = & 'tendency_of_' // trim(adjustl(self % constituent_name(self % index_mpas_scalar_to_constituent(j)))) + + ! Print information about MPAS scalar tendencies. Only do it once. + if (i == 1) then + call self % debug_print(log_level_info, 'MPAS scalar tendency index ' // stringify([j])) + call self % debug_print(log_level_info, ' MPAS scalar tendency name: ' // & + stringify([field_3d_real % constituentnames(j)])) + call self % debug_print(log_level_info, ' Is water species: ' // & + stringify([self % is_water_species(self % index_mpas_scalar_to_constituent(j))])) + call self % debug_print(log_level_info, ' Index mapping from MPAS scalar tendency to constituent: ' // & + stringify([j]) // ' -> ' // stringify([self % index_mpas_scalar_to_constituent(j)])) + end if end do nullify(field_3d_real) @@ -1101,11 +1137,11 @@ subroutine dyn_mpas_define_scalar(self, constituent_name, is_water_species) call mpas_pool_add_dimension(self % domain_ptr % blocklist % dimensions, 'moist_start', index_water_start) call mpas_pool_add_dimension(self % domain_ptr % blocklist % dimensions, 'moist_end', index_water_end) - call self % debug_print('index_qv = ' // stringify([index_qv])) - call self % debug_print('moist_start = ' // stringify([index_water_start])) - call self % debug_print('moist_end = ' // stringify([index_water_end])) + call self % debug_print(log_level_debug, 'index_qv = ' // stringify([index_qv])) + call self % debug_print(log_level_debug, 'moist_start = ' // stringify([index_water_start])) + call self % debug_print(log_level_debug, 'moist_end = ' // stringify([index_water_end])) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_define_scalar !------------------------------------------------------------------------------- @@ -1142,12 +1178,12 @@ subroutine dyn_mpas_read_write_stream(self, pio_file, stream_mode, stream_name) type(mpas_stream_type), pointer :: mpas_stream type(var_info_type), allocatable :: var_info_list(:) - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(mpas_pool) nullify(mpas_stream) - call self % debug_print('Initializing stream "' // trim(adjustl(stream_name)) // '"') + call self % debug_print(log_level_info, 'Initializing stream "' // trim(adjustl(stream_name)) // '"') call self % init_stream_with_pool(mpas_pool, mpas_stream, pio_file, stream_mode, stream_name) @@ -1161,7 +1197,7 @@ subroutine dyn_mpas_read_write_stream(self, pio_file, stream_mode, stream_name) select case (trim(adjustl(stream_mode))) case ('r', 'read') - call self % debug_print('Reading stream "' // trim(adjustl(stream_name)) // '"') + call self % debug_print(log_level_info, 'Reading stream "' // trim(adjustl(stream_name)) // '"') call mpas_readstream(mpas_stream, 1, ierr=ierr) @@ -1180,7 +1216,7 @@ subroutine dyn_mpas_read_write_stream(self, pio_file, stream_mode, stream_name) call postread_reindex(self % domain_ptr % blocklist % allfields, self % domain_ptr % packages, & mpas_pool, mpas_pool) case ('w', 'write') - call self % debug_print('Writing stream "' // trim(adjustl(stream_name)) // '"') + call self % debug_print(log_level_info, 'Writing stream "' // trim(adjustl(stream_name)) // '"') ! WARNING: ! The `{pre,post}write_reindex` subroutines are STATEFUL because they store information inside their module @@ -1202,7 +1238,7 @@ subroutine dyn_mpas_read_write_stream(self, pio_file, stream_mode, stream_name) call self % model_error('Unsupported stream mode "' // trim(adjustl(stream_mode)) // '"', subname, __LINE__) end select - call self % debug_print('Closing stream "' // trim(adjustl(stream_name)) // '"') + call self % debug_print(log_level_info, 'Closing stream "' // trim(adjustl(stream_name)) // '"') call mpas_closestream(mpas_stream, ierr=ierr) @@ -1217,7 +1253,7 @@ subroutine dyn_mpas_read_write_stream(self, pio_file, stream_mode, stream_name) deallocate(mpas_stream) nullify(mpas_stream) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_read_write_stream !------------------------------------------------------------------------------- @@ -1281,7 +1317,7 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file type(field5dreal), pointer :: field_5d_real type(var_info_type), allocatable :: var_info_list(:) - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(field_0d_char) nullify(field_1d_char) @@ -1308,7 +1344,7 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file stream_filename = 'external stream' stream_format = mpas_io_pnetcdf - call self % debug_print('Checking PIO file descriptor') + call self % debug_print(log_level_verbose, 'Checking PIO file descriptor') if (.not. associated(pio_file)) then call self % model_error('Invalid PIO file descriptor', subname, __LINE__) @@ -1320,14 +1356,14 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file select case (trim(adjustl(stream_mode))) case ('r', 'read') - call self % debug_print('Creating "' // trim(adjustl(stream_name)) // '" stream for reading') + call self % debug_print(log_level_verbose, 'Creating stream "' // trim(adjustl(stream_name)) // '" for reading') call mpas_createstream( & mpas_stream, self % domain_ptr % iocontext, stream_filename, stream_format, mpas_io_read, & clobberrecords=.false., clobberfiles=.false., truncatefiles=.false., & precision=mpas_io_native_precision, pio_file_desc=pio_file, ierr=ierr) case ('w', 'write') - call self % debug_print('Creating "' // trim(adjustl(stream_name)) // '" stream for writing') + call self % debug_print(log_level_verbose, 'Creating stream "' // trim(adjustl(stream_name)) // '" for writing') call mpas_createstream( & mpas_stream, self % domain_ptr % iocontext, stream_filename, stream_format, mpas_io_write, & @@ -1343,15 +1379,13 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file var_info_list = parse_stream_name(stream_name) - ! Add variables to stream. - call self % debug_print('Adding variables to stream') - + ! Add variables contained in `var_info_list` to stream. do i = 1, size(var_info_list) - call self % debug_print('var_info_list(' // stringify([i]) // ') % name = ' // & + call self % debug_print(log_level_debug, 'var_info_list(' // stringify([i]) // ') % name = ' // & stringify([var_info_list(i) % name])) - call self % debug_print('var_info_list(' // stringify([i]) // ') % type = ' // & + call self % debug_print(log_level_debug, 'var_info_list(' // stringify([i]) // ') % type = ' // & stringify([var_info_list(i) % type])) - call self % debug_print('var_info_list(' // stringify([i]) // ') % rank = ' // & + call self % debug_print(log_level_debug, 'var_info_list(' // stringify([i]) // ') % rank = ' // & stringify([var_info_list(i) % rank])) if (trim(adjustl(stream_mode)) == 'r' .or. trim(adjustl(stream_mode)) == 'read') then @@ -1361,14 +1395,14 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file ! This can happen if users attempt to initialize/restart the model with data generated by ! older versions of MPAS. Print a debug message to let users decide if this is acceptable. if (.not. any(var_is_present)) then - call self % debug_print('Skipping variable "' // trim(adjustl(var_info_list(i) % name)) // & + call self % debug_print(log_level_verbose, 'Skipping variable "' // trim(adjustl(var_info_list(i) % name)) // & '" due to not present') cycle end if if (any(var_is_present .and. .not. var_is_tkr_compatible)) then - call self % debug_print('Skipping variable "' // trim(adjustl(var_info_list(i) % name)) // & + call self % debug_print(log_level_verbose, 'Skipping variable "' // trim(adjustl(var_info_list(i) % name)) // & '" due to not TKR compatible') cycle @@ -1383,6 +1417,9 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file call mpas_pool_add_config(mpas_pool, trim(adjustl(var_info_list(i) % name) // ':packages'), '') ! Add "" to stream. + call self % debug_print(log_level_verbose, 'Adding variable "' // trim(adjustl(var_info_list(i) % name)) // & + '" to stream "' // trim(adjustl(stream_name)) // '"') + select case (trim(adjustl(var_info_list(i) % type))) case ('character') select case (var_info_list(i) % rank) @@ -1559,7 +1596,6 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file if (trim(adjustl(stream_mode)) == 'w' .or. trim(adjustl(stream_mode)) == 'write') then ! Add MPAS-specific attributes to stream. - call self % debug_print('Adding attributes to stream') ! Attributes related to MPAS core (i.e., `core_type`). call add_stream_attribute('conventions', self % domain_ptr % core % conventions) @@ -1578,7 +1614,7 @@ subroutine dyn_mpas_init_stream_with_pool(self, mpas_pool, mpas_stream, pio_file call add_stream_attribute('y_period', self % domain_ptr % y_period) end if - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') contains !> Helper subroutine for adding a 0-d stream attribute by calling `mpas_writestreamatt` with error checking. !> (KCW, 2024-03-14) @@ -1589,6 +1625,9 @@ subroutine add_stream_attribute_0d(attribute_name, attribute_value) character(*), intent(in) :: attribute_name class(*), intent(in) :: attribute_value + call self % debug_print(log_level_verbose, 'Adding attribute "' // trim(adjustl(attribute_name)) // & + '" to stream "' // trim(adjustl(stream_name)) // '"') + select type (attribute_value) type is (character(*)) call mpas_writestreamatt(mpas_stream, & @@ -1629,6 +1668,9 @@ subroutine add_stream_attribute_1d(attribute_name, attribute_value) character(*), intent(in) :: attribute_name class(*), intent(in) :: attribute_value(:) + call self % debug_print(log_level_verbose, 'Adding attribute "' // trim(adjustl(attribute_name)) // & + '" to stream "' // trim(adjustl(stream_name)) // '"') + select type (attribute_value) type is (integer) call mpas_writestreamatt(mpas_stream, & @@ -1921,7 +1963,7 @@ subroutine dyn_mpas_check_variable_status(self, var_is_present, var_is_tkr_compa type(field4dreal), pointer :: field_4d_real type(field5dreal), pointer :: field_5d_real - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(field_0d_char) nullify(field_1d_char) @@ -2237,6 +2279,9 @@ subroutine dyn_mpas_check_variable_status(self, var_is_present, var_is_tkr_compa return end if + call self % debug_print(log_level_verbose, 'Checking variable "' // trim(adjustl(var_info % name)) // & + '" for presence and TKR compatibility') + do i = 1, size(var_name_list) ! Check if the variable is present on the file. ierr = pio_inq_varid(pio_file, trim(adjustl(var_name_list(i))), varid) @@ -2293,11 +2338,11 @@ subroutine dyn_mpas_check_variable_status(self, var_is_present, var_is_tkr_compa var_is_tkr_compatible(i) = .true. end do - call self % debug_print('var_name_list = ' // stringify(var_name_list)) - call self % debug_print('var_is_present = ' // stringify(var_is_present)) - call self % debug_print('var_is_tkr_compatible = ' // stringify(var_is_tkr_compatible)) + call self % debug_print(log_level_debug, 'var_name_list = ' // stringify(var_name_list)) + call self % debug_print(log_level_debug, 'var_is_present = ' // stringify(var_is_present)) + call self % debug_print(log_level_debug, 'var_is_tkr_compatible = ' // stringify(var_is_tkr_compatible)) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_check_variable_status !------------------------------------------------------------------------------- @@ -2335,7 +2380,7 @@ subroutine dyn_mpas_exchange_halo(self, field_name) type(field5dreal), pointer :: field_5d_real type(mpas_pool_field_info_type) :: mpas_pool_field_info - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(field_1d_integer) nullify(field_2d_integer) @@ -2346,7 +2391,7 @@ subroutine dyn_mpas_exchange_halo(self, field_name) nullify(field_4d_real) nullify(field_5d_real) - call self % debug_print('Inquiring field information for "' // trim(adjustl(field_name)) // '"') + call self % debug_print(log_level_info, 'Inquiring field information for "' // trim(adjustl(field_name)) // '"') call mpas_pool_get_field_info(self % domain_ptr % blocklist % allfields, & trim(adjustl(field_name)), mpas_pool_field_info) @@ -2359,12 +2404,13 @@ subroutine dyn_mpas_exchange_halo(self, field_name) ! No halo layers to exchange. This field is not decomposed. if (mpas_pool_field_info % nhalolayers == 0) then - call self % debug_print('Skipping field "' // trim(adjustl(field_name)) // '"') + call self % debug_print(log_level_info, 'Skipping field "' // trim(adjustl(field_name)) // & + '" due to not decomposed') return end if - call self % debug_print('Exchanging halo layers for "' // trim(adjustl(field_name)) // '"') + call self % debug_print(log_level_info, 'Exchanging halo layers for "' // trim(adjustl(field_name)) // '"') select case (mpas_pool_field_info % fieldtype) case (mpas_pool_integer) @@ -2479,7 +2525,7 @@ subroutine dyn_mpas_exchange_halo(self, field_name) call self % model_error('Unsupported field type (Must be one of: integer, real)', subname, __LINE__) end select - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_exchange_halo !------------------------------------------------------------------------------- @@ -2517,7 +2563,7 @@ subroutine dyn_mpas_compute_unit_vector(self) real(rkind), pointer :: east(:, :), north(:, :) type(mpas_pool_type), pointer :: mpas_pool - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(ncells) nullify(latcell, loncell) @@ -2535,6 +2581,8 @@ subroutine dyn_mpas_compute_unit_vector(self) call self % get_variable_pointer(east, 'mesh', 'east') call self % get_variable_pointer(north, 'mesh', 'north') + call self % debug_print(log_level_info, 'Computing unit vectors') + do i = 1, ncells east(1, i) = -sin(loncell(i)) east(2, i) = cos(loncell(i)) @@ -2559,7 +2607,7 @@ subroutine dyn_mpas_compute_unit_vector(self) nullify(mpas_pool) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_compute_unit_vector !------------------------------------------------------------------------------- @@ -2595,7 +2643,7 @@ subroutine dyn_mpas_compute_edge_wind(self, wind_tendency) real(rkind), pointer :: ucellzonal(:, :), ucellmeridional(:, :) real(rkind), pointer :: uedge(:, :) - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(nedges) @@ -2639,6 +2687,12 @@ subroutine dyn_mpas_compute_edge_wind(self, wind_tendency) call self % get_variable_pointer(uedge, 'state', 'u', time_level=1) end if + if (wind_tendency) then + call self % debug_print(log_level_info, 'Computing edge-normal wind tendency vectors') + else + call self % debug_print(log_level_info, 'Computing edge-normal wind vectors') + end if + do i = 1, nedges cell1 = cellsonedge(1, i) cell2 = cellsonedge(2, i) @@ -2674,7 +2728,7 @@ subroutine dyn_mpas_compute_edge_wind(self, wind_tendency) call self % exchange_halo('u') end if - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_compute_edge_wind !------------------------------------------------------------------------------- @@ -2724,7 +2778,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) type(mpas_pool_type), pointer :: mpas_pool type(mpas_time_type) :: mpas_time - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(initial_time_1, initial_time_2) nullify(xtime) @@ -2757,9 +2811,10 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) self % coupling_time_interval = coupling_time_interval self % number_of_time_steps = 0 - call self % debug_print('Coupling time interval is ' // stringify([real(self % coupling_time_interval, rkind)]) // & - ' seconds') - call self % debug_print('Time step is ' // stringify([config_dt]) // ' seconds') + call self % debug_print(log_level_info, 'Coupling time interval is ' // & + stringify([real(self % coupling_time_interval, rkind)]) // ' seconds') + call self % debug_print(log_level_info, 'Time step is ' // & + stringify([config_dt]) // ' seconds') nullify(config_dt) @@ -2767,7 +2822,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) call mpas_constants_compute_derived() ! Set up OpenMP threading. - call self % debug_print('Setting up OpenMP threading') + call self % debug_print(log_level_info, 'Setting up OpenMP threading') call mpas_atm_threading_init(self % domain_ptr % blocklist, ierr=ierr) @@ -2777,7 +2832,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) end if ! Set up inner dimensions used by arrays in optimized dynamics subroutines. - call self % debug_print('Setting up dimensions') + call self % debug_print(log_level_info, 'Setting up dimensions') call self % get_variable_pointer(nvertlevels, 'dim', 'nVertLevels') call self % get_variable_pointer(maxedges, 'dim', 'maxEdges') @@ -2808,7 +2863,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) if (.not. config_do_restart) then ! Run type is initial run. - call self % debug_print('Initializing time levels') + call self % debug_print(log_level_info, 'Initializing time levels') call self % get_pool_pointer(mpas_pool, 'state') @@ -2828,7 +2883,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) ! Initialize atmospheric variables (e.g., momentum, thermodynamic... variables in governing equations) ! as well as various aspects of time in MPAS. - call self % debug_print('Initializing atmospheric variables') + call self % debug_print(log_level_info, 'Initializing atmospheric variables') ! Controlled by `config_start_time` in namelist. mpas_time = mpas_get_clock_time(self % domain_ptr % clock, mpas_start_time, ierr=ierr) @@ -2892,7 +2947,7 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) call self % model_error('Failed to exchange halo layers for group "initialization:pv_edge,ru,rw"', subname, __LINE__) end if - call self % debug_print('Initializing dynamics') + call self % debug_print(log_level_info, 'Initializing dynamics') ! Prepare dynamics for time integration. call mpas_atm_dynamics_init(self % domain_ptr) @@ -2907,9 +2962,9 @@ subroutine dyn_mpas_init_phase4(self, coupling_time_interval) call mpas_allocate_scratch_field(field_2d_real) nullify(field_2d_real) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') - call self % debug_print('Successful initialization of MPAS dynamical core') + call self % debug_print(log_level_info, 'Successful initialization of MPAS dynamical core') contains !> Test if `a` is divisible by `b`, where `a` and `b` are both reals. !> (KCW, 2024-05-25) @@ -2996,7 +3051,7 @@ subroutine dyn_mpas_run(self) type(mpas_time_type) :: mpas_time_end, mpas_time_now ! This derived type is analogous to `ESMF_Time`. type(mpas_timeinterval_type) :: mpas_time_interval ! This derived type is analogous to `ESMF_TimeInterval`. - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(config_dt) nullify(mpas_pool_diag, mpas_pool_mesh, mpas_pool_state) @@ -3018,7 +3073,7 @@ subroutine dyn_mpas_run(self) call self % model_error('Failed to get time for "mpas_now"', subname, __LINE__) end if - call self % debug_print('Time integration of MPAS dynamical core begins at ' // trim(adjustl(date_time))) + call self % debug_print(log_level_info, 'Time integration of MPAS dynamical core begins at ' // trim(adjustl(date_time))) call mpas_set_timeinterval(mpas_time_interval, s=self % coupling_time_interval, ierr=ierr) @@ -3055,7 +3110,7 @@ subroutine dyn_mpas_run(self) call self % model_error('Failed to get time for "mpas_now"', subname, __LINE__) end if - call self % debug_print('Time step ' // stringify([self % number_of_time_steps]) // ' completed') + call self % debug_print(log_level_info, 'Time step ' // stringify([self % number_of_time_steps]) // ' completed') end do call mpas_get_time(mpas_time_now, datetimestring=date_time, ierr=ierr) @@ -3064,7 +3119,7 @@ subroutine dyn_mpas_run(self) call self % model_error('Failed to get time for "mpas_now"', subname, __LINE__) end if - call self % debug_print('Time integration of MPAS dynamical core ends at ' // trim(adjustl(date_time))) + call self % debug_print(log_level_info, 'Time integration of MPAS dynamical core ends at ' // trim(adjustl(date_time))) ! Compute diagnostic variables like `pressure`, `rho` and `theta` from time level 1 of MPAS `state` pool ! by calling upstream MPAS functionality. @@ -3073,7 +3128,7 @@ subroutine dyn_mpas_run(self) nullify(config_dt) nullify(mpas_pool_diag, mpas_pool_mesh, mpas_pool_state) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_run !------------------------------------------------------------------------------- @@ -3113,7 +3168,7 @@ subroutine dyn_mpas_final(self) integer :: ierr type(field2dreal), pointer :: field_2d_real - call self % debug_print(subname // ' entered') + call self % debug_print(log_level_debug, subname // ' entered') nullify(field_2d_real) @@ -3129,6 +3184,8 @@ subroutine dyn_mpas_final(self) call mpas_deallocate_scratch_field(field_2d_real) nullify(field_2d_real) + call self % debug_print(log_level_info, 'Finalizing dynamics') + ! Opposite to `mpas_atm_dynamics_init`. call mpas_atm_dynamics_finalize(self % domain_ptr) @@ -3141,6 +3198,8 @@ subroutine dyn_mpas_final(self) nullify(exchange_halo_group) + call self % debug_print(log_level_info, 'Cleaning up OpenMP threading') + ! Opposite to `mpas_atm_threading_init`. call mpas_atm_threading_finalize(self % domain_ptr % blocklist, ierr=ierr) @@ -3148,6 +3207,8 @@ subroutine dyn_mpas_final(self) call self % model_error('Failed to clean up OpenMP threading', subname, __LINE__) end if + call self % debug_print(log_level_info, 'Cleaning up clock') + ! Opposite to `mpas_create_clock`, which was called by `atm_simulation_clock_init`, then `atm_setup_clock`. call mpas_destroy_clock(self % domain_ptr % clock, ierr=ierr) @@ -3155,6 +3216,8 @@ subroutine dyn_mpas_final(self) call self % model_error('Failed to clean up clock', subname, __LINE__) end if + call self % debug_print(log_level_info, 'Cleaning up decompositions') + ! Opposite to `mpas_decomp_create_decomp_list`, which was called by `atm_setup_decompositions`. call mpas_decomp_destroy_decomp_list(self % domain_ptr % decompositions) @@ -3167,6 +3230,8 @@ subroutine dyn_mpas_final(self) ! Opposite to `mpas_timer_init`, which was called by `mpas_framework_init_phase2`. call mpas_timer_finalize(self % domain_ptr) + call self % debug_print(log_level_info, 'Cleaning up log') + ! Opposite to `mpas_log_init`, which was called by `atm_setup_log`. call mpas_log_finalize(ierr) @@ -3174,12 +3239,14 @@ subroutine dyn_mpas_final(self) call self % model_error('Failed to clean up log', subname, __LINE__) end if + call self % debug_print(log_level_info, 'Finalizing MPAS framework') + ! Opposite to `mpas_framework_init_phase1` and `mpas_framework_init_phase2`. call mpas_framework_finalize(self % domain_ptr % dminfo, self % domain_ptr) - call self % debug_print(subname // ' completed') + call self % debug_print(log_level_debug, subname // ' completed') - call self % debug_print('Successful finalization of MPAS dynamical core') + call self % debug_print(log_level_info, 'Successful finalization of MPAS dynamical core') ! Second, clean up this MPAS dynamical core instance. @@ -3197,6 +3264,7 @@ subroutine dyn_mpas_final(self) self % number_of_constituents = 0 ! Initialized by `dyn_mpas_init_phase1`. + self % log_level = log_level_quiet self % log_unit = output_unit self % mpi_comm = mpi_comm_null self % mpi_rank = 0 @@ -3209,8 +3277,6 @@ subroutine dyn_mpas_final(self) nullify(self % corelist) nullify(self % domain_ptr) - - self % debug_output = .false. end subroutine dyn_mpas_final !------------------------------------------------------------------------------- @@ -3377,6 +3443,8 @@ subroutine dyn_mpas_get_local_mesh_dimension(self, & integer, pointer :: nverticessolve_pointer integer, pointer :: nvertlevels_pointer + call self % debug_print(log_level_debug, subname // ' entered') + nullify(ncells_pointer) nullify(ncellssolve_pointer) nullify(nedges_pointer) @@ -3385,6 +3453,8 @@ subroutine dyn_mpas_get_local_mesh_dimension(self, & nullify(nverticessolve_pointer) nullify(nvertlevels_pointer) + call self % debug_print(log_level_info, 'Inquiring local mesh dimensions') + call self % get_variable_pointer(ncells_pointer, 'dim', 'nCells') call self % get_variable_pointer(ncellssolve_pointer, 'dim', 'nCellsSolve') call self % get_variable_pointer(nedges_pointer, 'dim', 'nEdges') @@ -3411,6 +3481,8 @@ subroutine dyn_mpas_get_local_mesh_dimension(self, & nullify(nvertices_pointer) nullify(nverticessolve_pointer) nullify(nvertlevels_pointer) + + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_get_local_mesh_dimension !------------------------------------------------------------------------------- @@ -3446,12 +3518,16 @@ subroutine dyn_mpas_get_global_mesh_dimension(self, & integer, pointer :: nverticessolve_pointer integer, pointer :: nvertlevels_pointer + call self % debug_print(log_level_debug, subname // ' entered') + nullify(maxedges_pointer) nullify(ncellssolve_pointer) nullify(nedgessolve_pointer) nullify(nverticessolve_pointer) nullify(nvertlevels_pointer) + call self % debug_print(log_level_info, 'Inquiring global mesh dimensions') + call self % get_variable_pointer(maxedges_pointer, 'dim', 'maxEdges') call self % get_variable_pointer(ncellssolve_pointer, 'dim', 'nCellsSolve') call self % get_variable_pointer(nedgessolve_pointer, 'dim', 'nEdgesSolve') @@ -3476,6 +3552,8 @@ subroutine dyn_mpas_get_global_mesh_dimension(self, & nullify(nedgessolve_pointer) nullify(nverticessolve_pointer) nullify(nvertlevels_pointer) + + call self % debug_print(log_level_debug, subname // ' completed') end subroutine dyn_mpas_get_global_mesh_dimension !> Helper subroutine for returning a pointer of `mpas_pool_type` to the named pool. diff --git a/src/dynamics/mpas/dyn_comp.F90 b/src/dynamics/mpas/dyn_comp.F90 index 453c6bf3..4c8e15a2 100644 --- a/src/dynamics/mpas/dyn_comp.F90 +++ b/src/dynamics/mpas/dyn_comp.F90 @@ -47,21 +47,21 @@ module dyn_comp integer, protected :: ncells_global, nedges_global, nvertices_global, ncells_max, nedges_max real(kind_dyn_mpas), protected :: sphere_radius contains - !> Print a debug message with optionally the value(s) of a variable. + !> Print a debug message at a debug level. !> If `printer` is not supplied, the MPI root rank will print. Otherwise, the designated MPI rank will print instead. !> (KCW, 2024-02-03) - subroutine dyn_debug_print(message, variable, printer) + subroutine dyn_debug_print(level, message, printer) ! Module(s) from CAM-SIMA. - use cam_logfile, only: debug_output, debugout_none, iulog + use cam_logfile, only: debug_output, iulog use spmd_utils, only: iam, masterproc use string_utils, only: stringify + integer, intent(in) :: level character(*), intent(in) :: message - class(*), optional, intent(in) :: variable(:) integer, optional, intent(in) :: printer - ! Bail out early if debug output is not requested. - if (.not. (debug_output > debugout_none)) then + ! Bail out early if the log level is less verbose than the debug level. + if (debug_output < level) then return end if @@ -75,16 +75,10 @@ subroutine dyn_debug_print(message, variable, printer) end if end if - if (present(variable)) then - write(iulog, '(a)') 'dyn_debug_print (' // stringify([iam]) // '): ' // & - message // stringify(variable) - else - write(iulog, '(a)') 'dyn_debug_print (' // stringify([iam]) // '): ' // & - message - end if + write(iulog, '(a)') 'MPAS Interface (' // stringify([iam]) // '): ' // message end subroutine dyn_debug_print - !> Read MPAS namelist from supplied path. + !> Read MPAS namelist from the supplied path. !> Additionally, perform early initialization of MPAS dynamical core. !> (KCW, 2024-02-09) ! @@ -94,8 +88,9 @@ subroutine dyn_readnl(namelist_path) use cam_abortutils, only: endrun use cam_control_mod, only: initial_run use cam_instance, only: atm_id - use cam_logfile, only: debug_output, debugout_none, iulog + use cam_logfile, only: debug_output, debugout_debug, debugout_info, iulog use spmd_utils, only: mpicom + use string_utils, only: stringify use time_manager, only: get_start_date, get_stop_date, get_run_duration, timemgr_get_calendar_cf ! Module(s) from CESM Share. use shr_file_mod, only: shr_file_getunit @@ -114,17 +109,19 @@ subroutine dyn_readnl(namelist_path) sec_since_midnight ! Second(s) since midnight. type(iosystem_desc_t), pointer :: pio_iosystem - nullify(pio_iosystem) + call dyn_debug_print(debugout_debug, subname // ' entered') - ! Enable/disable the debug output of MPAS dynamical core according to the debug verbosity level of CAM-SIMA. - mpas_dynamical_core % debug_output = (debug_output > debugout_none) + nullify(pio_iosystem) ! Get free units for MPAS so it can write its own log files, e.g., `log.atmosphere.0000.{out,err}`. log_unit(1) = shr_file_getunit() log_unit(2) = shr_file_getunit() - ! Initialize MPAS framework with supplied MPI communicator group and log units. - call mpas_dynamical_core % init_phase1(mpicom, endrun, iulog, log_unit) + call dyn_debug_print(debugout_info, 'Initializing MPAS dynamical core (Phase 1/4)') + + ! Initialize MPAS framework with the supplied MPI communicator group, procedure pointer to terminate the model, + ! log level, and units. + call mpas_dynamical_core % init_phase1(mpicom, endrun, debug_output, iulog, log_unit) cam_calendar = timemgr_get_calendar_cf() @@ -137,16 +134,22 @@ subroutine dyn_readnl(namelist_path) call get_run_duration(run_duration(1), sec_since_midnight) run_duration(2:4) = sec_to_hour_min_sec(sec_since_midnight) + call dyn_debug_print(debugout_info, 'Reading namelist') + ! Read MPAS-related namelist variables from `namelist_path`, e.g., `atm_in`. call mpas_dynamical_core % read_namelist(namelist_path, & cam_calendar, start_date_time, stop_date_time, run_duration, initial_run) pio_iosystem => shr_pio_getiosys(atm_id) - ! Initialize MPAS framework with supplied PIO system descriptor. + call dyn_debug_print(debugout_info, 'Initializing MPAS dynamical core (Phase 2/4)') + + ! Initialize MPAS framework with the supplied PIO system descriptor. call mpas_dynamical_core % init_phase2(pio_iosystem) nullify(pio_iosystem) + + call dyn_debug_print(debugout_debug, subname // ' completed') contains !> Convert second(s) to hour(s), minute(s), and second(s). !> (KCW, 2024-02-07) @@ -180,11 +183,13 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) use cam_constituents, only: const_name, const_is_water_species, num_advected, readtrace use cam_control_mod, only: initial_run use cam_initfiles, only: initial_file_get_id, topo_file_get_id + use cam_logfile, only: debugout_debug, debugout_info use cam_pio_utils, only: clean_iodesc_list use cam_thermo_formula, only: energy_formula_dycore, energy_formula_dycore_mpas use inic_analytic, only: analytic_ic_active use physics_types, only: dycore_energy_consistency_adjust use runtime_obj, only: runtime_options + use string_utils, only: stringify use time_manager, only: get_step_size ! Module(s) from CCPP. use phys_vars_init_check, only: std_name_len @@ -204,6 +209,8 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) type(file_desc_t), pointer :: pio_init_file type(file_desc_t), pointer :: pio_topo_file + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(pio_init_file) nullify(pio_topo_file) @@ -226,6 +233,8 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) is_water_species(i) = const_is_water_species(i) end do + call dyn_debug_print(debugout_info, 'Defining MPAS scalars and scalar tendencies') + ! Inform MPAS about constituent names and their corresponding waterness. call mpas_dynamical_core % define_scalar(constituent_name, is_water_species) @@ -248,26 +257,40 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) mpas_dynamical_core % map_mpas_scalar_index(thermodynamic_active_species_ice_idx(i)) end do + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_num = ' // & + stringify([thermodynamic_active_species_num])) + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_liq_num = ' // & + stringify([thermodynamic_active_species_liq_num])) + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_ice_num = ' // & + stringify([thermodynamic_active_species_ice_num])) + + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_idx_dycore = ' // & + stringify(thermodynamic_active_species_idx_dycore)) + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_liq_idx_dycore = ' // & + stringify(thermodynamic_active_species_liq_idx_dycore)) + call dyn_debug_print(debugout_debug, 'thermodynamic_active_species_ice_idx_dycore = ' // & + stringify(thermodynamic_active_species_ice_idx_dycore)) + pio_init_file => initial_file_get_id() pio_topo_file => topo_file_get_id() if (initial_run) then ! Run type is initial run. - call dyn_debug_print('Calling check_topography_data') + call dyn_debug_print(debugout_info, 'Checking for consistency in topography data') call check_topography_data(pio_topo_file) if (analytic_ic_active()) then - call dyn_debug_print('Calling set_analytic_initial_condition') + call dyn_debug_print(debugout_info, 'Initializing MPAS state variables by setting analytic initial condition') call set_analytic_initial_condition() else + call dyn_debug_print(debugout_info, 'Initializing MPAS state variables by reading initial condition') + ! Perform default initialization for all constituents. ! Subsequently, they can be overridden depending on the namelist option (below) and ! the actual availability (checked and handled by MPAS). - call dyn_debug_print('Calling dyn_exchange_constituent_states') - call dyn_exchange_constituent_states(direction='e', exchange=.true., conversion=.false.) ! Namelist option that controls if constituents are to be read from the file. @@ -282,6 +305,8 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) else ! Run type is branch or restart run. + call dyn_debug_print(debugout_info, 'Initializing MPAS state variables by restarting') + ! Read variables that belong to the "input" and "restart" streams in MPAS. call mpas_dynamical_core % read_write_stream(pio_init_file, 'r', 'input+restart') end if @@ -297,8 +322,12 @@ subroutine dyn_init(cam_runtime_opts, dyn_in, dyn_out) ! then yield control back to the caller. coupling_time_interval = get_step_size() + call dyn_debug_print(debugout_info, 'Initializing MPAS dynamical core (Phase 4/4)') + ! Finish MPAS dynamical core initialization. After this point, MPAS dynamical core is ready for time integration. call mpas_dynamical_core % init_phase4(coupling_time_interval) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine dyn_init !> Check for consistency in topography data. The presence of topography file is inferred from the `pio_file` pointer. @@ -310,7 +339,7 @@ subroutine check_topography_data(pio_file) ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate, endrun use cam_field_read, only: cam_read_field - use cam_logfile, only: debug_output, debugout_none + use cam_logfile, only: debug_output, debugout_debug, debugout_none, debugout_verbose use dynconst, only: constant_g => gravit ! Module(s) from external libraries. use pio, only: file_desc_t, pio_file_is_open @@ -325,12 +354,14 @@ subroutine check_topography_data(pio_file) real(kind_r8), allocatable :: surface_geopotential(:) ! Read from topography file. real(kind_dyn_mpas), pointer :: zgrid(:, :) ! From MPAS. Geometric height (meters) at layer interfaces. + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(zgrid) call mpas_dynamical_core % get_variable_pointer(zgrid, 'mesh', 'zgrid') if (associated(pio_file)) then - call dyn_debug_print('Topography file is used') + call dyn_debug_print(debugout_verbose, 'Topography file is used for consistency check') if (.not. pio_file_is_open(pio_file)) then call endrun('Invalid PIO file descriptor', subname, __LINE__) @@ -362,7 +393,7 @@ subroutine check_topography_data(pio_file) deallocate(surface_geopotential) deallocate(surface_geometric_height) else - call dyn_debug_print('Topography file is not used') + call dyn_debug_print(debugout_verbose, 'Topography file is not used for consistency check') ! Surface geometric height in MPAS should be zero. if (any(abs(real(zgrid(1, 1:ncells_solve), kind_r8)) > error_tolerance)) then @@ -371,11 +402,16 @@ subroutine check_topography_data(pio_file) end if nullify(zgrid) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine check_topography_data !> Set analytic initial condition for MPAS. !> (KCW, 2024-05-22) subroutine set_analytic_initial_condition() + ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_debug + character(*), parameter :: subname = 'dyn_comp::set_analytic_initial_condition' integer, allocatable :: global_grid_index(:) real(kind_r8), allocatable :: buffer_2d_real(:, :), buffer_3d_real(:, :, :) @@ -385,6 +421,8 @@ subroutine set_analytic_initial_condition() real(kind_dyn_mpas), pointer :: zgrid(:, :) ! Geometric height (meters) at layer interfaces. ! Dimension and vertical index orders follow MPAS convention. + call dyn_debug_print(debugout_debug, subname // ' entered') + call init_shared_variables() call set_mpas_state_u() @@ -394,6 +432,8 @@ subroutine set_analytic_initial_condition() call set_mpas_state_rho_base_theta_base() call final_shared_variables() + + call dyn_debug_print(debugout_debug, subname // ' completed') contains !> Initialize variables that are shared and repeatedly used by the `set_mpas_state_*` internal subroutines. !> (KCW, 2024-05-13) @@ -401,6 +441,7 @@ subroutine init_shared_variables() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate, endrun use cam_grid_support, only: cam_grid_get_latvals, cam_grid_get_lonvals, cam_grid_id + use cam_logfile, only: debugout_verbose use dynconst, only: deg_to_rad use vert_coord, only: pverp @@ -410,7 +451,7 @@ subroutine init_shared_variables() integer, pointer :: indextocellid(:) real(kind_r8), pointer :: lat_deg(:), lon_deg(:) - call dyn_debug_print('Preparing to set analytic initial condition') + call dyn_debug_print(debugout_verbose, 'Preparing to set analytic initial condition') nullify(zgrid) nullify(indextocellid) @@ -476,6 +517,7 @@ end subroutine final_shared_variables subroutine set_mpas_state_u() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate + use cam_logfile, only: debugout_verbose use dyn_tests_utils, only: vc_height use inic_analytic, only: dyn_set_inic_col use vert_coord, only: pver @@ -485,7 +527,7 @@ subroutine set_mpas_state_u() integer :: ierr real(kind_dyn_mpas), pointer :: ucellzonal(:, :), ucellmeridional(:, :) - call dyn_debug_print('Setting MPAS state "u"') + call dyn_debug_print(debugout_verbose, 'Setting MPAS state "u"') nullify(ucellzonal, ucellmeridional) @@ -523,10 +565,13 @@ end subroutine set_mpas_state_u !> Set MPAS state `w` (i.e., vertical velocity at cell interfaces). !> (KCW, 2024-05-13) subroutine set_mpas_state_w() + ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_verbose + character(*), parameter :: subname = 'dyn_comp::set_analytic_initial_condition::set_mpas_state_w' real(kind_dyn_mpas), pointer :: w(:, :) - call dyn_debug_print('Setting MPAS state "w"') + call dyn_debug_print(debugout_verbose, 'Setting MPAS state "w"') nullify(w) @@ -546,6 +591,7 @@ subroutine set_mpas_state_scalars() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate use cam_constituents, only: num_advected + use cam_logfile, only: debugout_verbose use dyn_tests_utils, only: vc_height use inic_analytic, only: dyn_set_inic_col use vert_coord, only: pver @@ -561,7 +607,7 @@ subroutine set_mpas_state_scalars() integer, pointer :: index_qv real(kind_dyn_mpas), pointer :: scalars(:, :, :) - call dyn_debug_print('Setting MPAS state "scalars"') + call dyn_debug_print(debugout_verbose, 'Setting MPAS state "scalars"') nullify(index_qv) nullify(scalars) @@ -593,10 +639,10 @@ subroutine set_mpas_state_scalars() if (mpas_dynamical_core % get_constituent_name(mpas_dynamical_core % map_constituent_index(index_qv)) == & constituent_qv_standard_name) then ! The definition of `qv` matches exactly what MPAS wants. No conversion is needed. - call dyn_debug_print('No conversion is needed for water vapor mixing ratio') + call dyn_debug_print(debugout_verbose, 'No conversion is needed for water vapor mixing ratio') else ! The definition of `qv` actually represents specific humidity. Conversion is needed. - call dyn_debug_print('Conversion is needed and applied for water vapor mixing ratio') + call dyn_debug_print(debugout_verbose, 'Conversion is needed and applied for water vapor mixing ratio') ! Convert specific humidity to water vapor mixing ratio. scalars(index_qv, :, 1:ncells_solve) = & @@ -618,6 +664,7 @@ end subroutine set_mpas_state_scalars subroutine set_mpas_state_rho_theta() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate + use cam_logfile, only: debugout_verbose use dyn_tests_utils, only: vc_height use dynconst, only: constant_p0 => pref, constant_rd => rair, constant_rv => rh2o use inic_analytic, only: dyn_set_inic_col @@ -640,7 +687,7 @@ subroutine set_mpas_state_rho_theta() real(kind_dyn_mpas), pointer :: theta(:, :) real(kind_dyn_mpas), pointer :: scalars(:, :, :) - call dyn_debug_print('Setting MPAS state "rho" and "theta"') + call dyn_debug_print(debugout_verbose, 'Setting MPAS state "rho" and "theta"') nullify(index_qv) nullify(rho) @@ -737,6 +784,7 @@ end subroutine set_mpas_state_rho_theta subroutine set_mpas_state_rho_base_theta_base() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate + use cam_logfile, only: debugout_verbose use dynconst, only: constant_p0 => pref, constant_rd => rair use vert_coord, only: pver @@ -750,7 +798,7 @@ subroutine set_mpas_state_rho_base_theta_base() real(kind_dyn_mpas), pointer :: theta_base(:, :) real(kind_dyn_mpas), pointer :: zz(:, :) - call dyn_debug_print('Setting MPAS state "rho_base" and "theta_base"') + call dyn_debug_print(debugout_verbose, 'Setting MPAS state "rho_base" and "theta_base"') nullify(rho_base) nullify(theta_base) @@ -840,6 +888,7 @@ subroutine dyn_exchange_constituent_states(direction, exchange, conversion) ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate, endrun use cam_constituents, only: const_is_dry, const_is_water_species, num_advected + use cam_logfile, only: debugout_debug, debugout_info use physics_types, only: phys_state use vert_coord, only: pver ! Module(s) from CCPP. @@ -860,22 +909,24 @@ subroutine dyn_exchange_constituent_states(direction, exchange, conversion) real(kind_r8), allocatable :: sigma_all_q(:) ! Summation of all water species mixing ratios. real(kind_dyn_mpas), pointer :: scalars(:, :, :) ! This points to MPAS memory. + call dyn_debug_print(debugout_debug, subname // ' entered') + select case (trim(adjustl(direction))) case ('e', 'export') if (exchange) then - call dyn_debug_print('Setting MPAS state "scalars" from physics state "constituents"') + call dyn_debug_print(debugout_info, 'Setting MPAS state "scalars" from physics state "constituents"') end if if (conversion) then - call dyn_debug_print('Converting MPAS state "scalars"') + call dyn_debug_print(debugout_info, 'Converting MPAS state "scalars"') end if case ('i', 'import') if (exchange) then - call dyn_debug_print('Setting physics state "constituents" from MPAS state "scalars"') + call dyn_debug_print(debugout_info, 'Setting physics state "constituents" from MPAS state "scalars"') end if if (conversion) then - call dyn_debug_print('Converting physics state "constituents"') + call dyn_debug_print(debugout_info, 'Converting physics state "constituents"') end if case default call endrun('Unsupported exchange direction "' // trim(adjustl(direction)) // '"', subname, __LINE__) @@ -982,17 +1033,22 @@ subroutine dyn_exchange_constituent_states(direction, exchange, conversion) ! Because we are injecting data directly into MPAS memory, halo layers need to be updated manually. call mpas_dynamical_core % exchange_halo('scalars') end if + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine dyn_exchange_constituent_states !> Inquire local and global mesh dimensions. Save them as protected module variables. !> (KCW, 2024-11-21) subroutine dyn_inquire_mesh_dimensions() ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_debug, debugout_info use string_utils, only: stringify character(*), parameter :: subname = 'dyn_comp::dyn_inquire_mesh_dimensions' - call dyn_debug_print('Inquiring local and global mesh dimensions') + call dyn_debug_print(debugout_debug, subname // ' entered') + + call dyn_debug_print(debugout_info, 'Inquiring local and global mesh dimensions') call mpas_dynamical_core % get_local_mesh_dimension( & ncells, ncells_solve, nedges, nedges_solve, nvertices, nvertices_solve, nvertlevels) @@ -1001,13 +1057,15 @@ subroutine dyn_inquire_mesh_dimensions() ncells_global, nedges_global, nvertices_global, nvertlevels, ncells_max, nedges_max, & sphere_radius) - call dyn_debug_print('ncells_global = ' // stringify([ncells_global])) - call dyn_debug_print('nedges_global = ' // stringify([nedges_global])) - call dyn_debug_print('nvertices_global = ' // stringify([nvertices_global])) - call dyn_debug_print('nvertlevels = ' // stringify([nvertlevels])) - call dyn_debug_print('ncells_max = ' // stringify([ncells_max])) - call dyn_debug_print('nedges_max = ' // stringify([nedges_max])) - call dyn_debug_print('sphere_radius = ' // stringify([sphere_radius])) + call dyn_debug_print(debugout_debug, 'ncells_global = ' // stringify([ncells_global])) + call dyn_debug_print(debugout_debug, 'nedges_global = ' // stringify([nedges_global])) + call dyn_debug_print(debugout_debug, 'nvertices_global = ' // stringify([nvertices_global])) + call dyn_debug_print(debugout_debug, 'nvertlevels = ' // stringify([nvertlevels])) + call dyn_debug_print(debugout_debug, 'ncells_max = ' // stringify([ncells_max])) + call dyn_debug_print(debugout_debug, 'nedges_max = ' // stringify([nedges_max])) + call dyn_debug_print(debugout_debug, 'sphere_radius = ' // stringify([sphere_radius])) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine dyn_inquire_mesh_dimensions !> Mark everything in the `physics_types` module along with constituents as initialized @@ -1016,12 +1074,15 @@ end subroutine dyn_inquire_mesh_dimensions subroutine mark_variables_as_initialized() ! Module(s) from CAM-SIMA. use cam_constituents, only: const_name, num_advected + use cam_logfile, only: debugout_debug ! Module(s) from CCPP. use phys_vars_init_check, only: mark_as_initialized character(*), parameter :: subname = 'dyn_comp::mark_variables_as_initialized' integer :: i + call dyn_debug_print(debugout_debug, subname // ' entered') + ! The variables below are managed by dynamics interface. ! We are responsible for initializing and updating them. @@ -1074,28 +1135,48 @@ subroutine mark_variables_as_initialized() call mark_as_initialized('vertically_integrated_total_energy_using_physics_energy_formula_at_start_of_physics_timestep') call mark_as_initialized('vertically_integrated_total_water') call mark_as_initialized('vertically_integrated_total_water_at_start_of_physics_timestep') + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine mark_variables_as_initialized !> Run MPAS dynamical core to integrate the dynamical states with time. !> (KCW, 2024-07-11) subroutine dyn_run() + ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_debug, debugout_info + character(*), parameter :: subname = 'dyn_comp::dyn_run' + call dyn_debug_print(debugout_debug, subname // ' entered') + + call dyn_debug_print(debugout_info, 'Running MPAS dynamical core') + ! MPAS dynamical core will run until the coupling time interval is reached. call mpas_dynamical_core % run() + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine dyn_run !> Finalize MPAS dynamical core as well as its framework. !> (KCW, 2024-10-04) subroutine dyn_final() + ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_debug, debugout_info + character(*), parameter :: subname = 'dyn_comp::dyn_final' + call dyn_debug_print(debugout_debug, subname // ' entered') + + call dyn_debug_print(debugout_info, 'Finalizing MPAS dynamical core') + ! Quick hack for dumping variables from MPAS dynamical core. ! Remove it once history and restart are wired up in CAM-SIMA. call dyn_variable_dump() ! After this point, do not access anything under MPAS dynamical core or runtime errors will ensue. call mpas_dynamical_core % final() + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine dyn_final subroutine dyn_variable_dump() diff --git a/src/dynamics/mpas/dyn_coupling.F90 b/src/dynamics/mpas/dyn_coupling.F90 index 831bcd02..eeca26e2 100644 --- a/src/dynamics/mpas/dyn_coupling.F90 +++ b/src/dynamics/mpas/dyn_coupling.F90 @@ -17,6 +17,7 @@ module dyn_coupling !> (KCW, 2024-07-31) subroutine dynamics_to_physics_coupling() ! Module(s) from CAM-SIMA. + use cam_logfile, only: debugout_debug, debugout_info use dyn_comp, only: dyn_debug_print, dyn_exchange_constituent_states, ncells_solve character(*), parameter :: subname = 'dyn_coupling::dynamics_to_physics_coupling' @@ -53,11 +54,13 @@ subroutine dynamics_to_physics_coupling() real(kind_dyn_mpas), pointer :: zgrid(:, :) real(kind_dyn_mpas), pointer :: zz(:, :) + call dyn_debug_print(debugout_debug, subname // ' entered') + call init_shared_variables() call dyn_exchange_constituent_states(direction='i', exchange=.true., conversion=.false.) - call dyn_debug_print('Setting physics state variables column by column') + call dyn_debug_print(debugout_info, 'Setting physics state variables column by column') ! Set variables in the `physics_state` derived type column by column. ! This way, peak memory usage can be reduced. @@ -69,6 +72,8 @@ subroutine dynamics_to_physics_coupling() call set_physics_state_external() call final_shared_variables() + + call dyn_debug_print(debugout_debug, subname // ' completed') contains !> Initialize variables that are shared and repeatedly used by the `update_shared_variables` and !> `set_physics_state_column` internal subroutines. @@ -85,7 +90,7 @@ subroutine init_shared_variables() integer :: ierr logical, allocatable :: is_water_species(:) - call dyn_debug_print('Preparing for dynamics-physics coupling') + call dyn_debug_print(debugout_info, 'Preparing for dynamics-physics coupling') nullify(index_qv) nullify(exner) @@ -317,7 +322,7 @@ subroutine set_physics_state_external() real(kind_phys), pointer :: constituents(:, :, :) type(ccpp_constituent_prop_ptr_t), pointer :: constituent_properties(:) - call dyn_debug_print('Setting physics state variables externally') + call dyn_debug_print(debugout_info, 'Setting physics state variables externally') nullify(constituents) nullify(constituent_properties) @@ -410,7 +415,8 @@ end subroutine dynamics_to_physics_coupling !> (KCW, 2024-09-20) subroutine physics_to_dynamics_coupling() ! Module(s) from CAM-SIMA. - use dyn_comp, only: dyn_exchange_constituent_states + use cam_logfile, only: debugout_debug + use dyn_comp, only: dyn_debug_print, dyn_exchange_constituent_states character(*), parameter :: subname = 'dyn_coupling::physics_to_dynamics_coupling' integer, pointer :: index_qv @@ -420,6 +426,8 @@ subroutine physics_to_dynamics_coupling() real(kind_dyn_mpas), pointer :: scalars(:, :, :) real(kind_dyn_mpas), pointer :: zz(:, :) + call dyn_debug_print(debugout_debug, subname // ' entered') + call init_shared_variables() call dyn_exchange_constituent_states(direction='e', exchange=.true., conversion=.true.) @@ -429,19 +437,22 @@ subroutine physics_to_dynamics_coupling() call set_mpas_physics_tendency_rtheta() call final_shared_variables() + + call dyn_debug_print(debugout_debug, subname // ' completed') contains !> Initialize variables that are shared and repeatedly used by the `set_mpas_physics_tendency_*` internal subroutines. !> (KCW, 2024-09-13) subroutine init_shared_variables() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate - use dyn_comp, only: dyn_debug_print, mpas_dynamical_core, ncells_solve + use cam_logfile, only: debugout_info + use dyn_comp, only: mpas_dynamical_core, ncells_solve use vert_coord, only: pver character(*), parameter :: subname = 'dyn_coupling::physics_to_dynamics_coupling::init_shared_variables' integer :: ierr - call dyn_debug_print('Preparing for physics-dynamics coupling') + call dyn_debug_print(debugout_info, 'Preparing for physics-dynamics coupling') nullify(index_qv) nullify(rho_zz) @@ -481,14 +492,15 @@ end subroutine final_shared_variables !> (KCW, 2024-09-11) subroutine set_mpas_physics_tendency_ru() ! Module(s) from CAM-SIMA. - use dyn_comp, only: dyn_debug_print, reverse, mpas_dynamical_core, ncells_solve + use cam_logfile, only: debugout_info + use dyn_comp, only: reverse, mpas_dynamical_core, ncells_solve use physics_types, only: phys_tend character(*), parameter :: subname = 'dyn_coupling::physics_to_dynamics_coupling::set_mpas_physics_tendency_ru' integer :: i real(kind_dyn_mpas), pointer :: u_tendency(:, :), v_tendency(:, :) - call dyn_debug_print('Setting MPAS physics tendency "tend_ru_physics"') + call dyn_debug_print(debugout_info, 'Setting MPAS physics tendency "tend_ru_physics"') nullify(u_tendency, v_tendency) @@ -512,12 +524,13 @@ end subroutine set_mpas_physics_tendency_ru !> (KCW, 2024-09-11) subroutine set_mpas_physics_tendency_rho() ! Module(s) from CAM-SIMA. - use dyn_comp, only: dyn_debug_print, mpas_dynamical_core, ncells_solve + use cam_logfile, only: debugout_info + use dyn_comp, only: mpas_dynamical_core, ncells_solve character(*), parameter :: subname = 'dyn_coupling::physics_to_dynamics_coupling::set_mpas_physics_tendency_rho' real(kind_dyn_mpas), pointer :: rho_tendency(:, :) - call dyn_debug_print('Setting MPAS physics tendency "tend_rho_physics"') + call dyn_debug_print(debugout_info, 'Setting MPAS physics tendency "tend_rho_physics"') nullify(rho_tendency) @@ -538,7 +551,8 @@ end subroutine set_mpas_physics_tendency_rho subroutine set_mpas_physics_tendency_rtheta() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate - use dyn_comp, only: dyn_debug_print, reverse, mpas_dynamical_core, ncells_solve + use cam_logfile, only: debugout_info + use dyn_comp, only: reverse, mpas_dynamical_core, ncells_solve use dynconst, only: constant_rd => rair, constant_rv => rh2o use physics_types, only: dtime_phys, phys_tend use vert_coord, only: pver @@ -558,7 +572,7 @@ subroutine set_mpas_physics_tendency_rtheta() real(kind_dyn_mpas), pointer :: theta_m(:, :) real(kind_dyn_mpas), pointer :: theta_m_tendency(:, :) - call dyn_debug_print('Setting MPAS physics tendency "tend_rtheta_physics"') + call dyn_debug_print(debugout_info, 'Setting MPAS physics tendency "tend_rtheta_physics"') nullify(theta_m) nullify(theta_m_tendency) diff --git a/src/dynamics/mpas/dyn_grid.F90 b/src/dynamics/mpas/dyn_grid.F90 index 2b4e3468..8a21b2a7 100644 --- a/src/dynamics/mpas/dyn_grid.F90 +++ b/src/dynamics/mpas/dyn_grid.F90 @@ -37,6 +37,7 @@ subroutine model_grid_init() use cam_abortutils, only: endrun use cam_constituents, only: num_advected use cam_initfiles, only: initial_file_get_id + use cam_logfile, only: debugout_debug, debugout_info use dyn_comp, only: dyn_debug_print, dyn_inquire_mesh_dimensions, mpas_dynamical_core, nvertlevels use dynconst, only: dynconst_init use string_utils, only: stringify @@ -47,27 +48,31 @@ subroutine model_grid_init() character(*), parameter :: subname = 'dyn_grid::model_grid_init' type(file_desc_t), pointer :: pio_file + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(pio_file) ! Initialize mathematical and physical constants for dynamics. - call dyn_debug_print('Calling dynconst_init') - call dynconst_init() - ! Initialize vertical coordinates. + call dyn_debug_print(debugout_info, 'Initializing vertical coordinate indexes') + + ! Initialize vertical coordinate indexes. ! `pver` comes from CAM-SIMA namelist. ! `pverp` is set only after this call. Call it as soon as possible. - call dyn_debug_print('Calling vert_coord_init') - call vert_coord_init(1, pver) ! Get PIO file descriptor for initial file. pio_file => initial_file_get_id() + call dyn_debug_print(debugout_info, 'Initializing MPAS dynamical core (Phase 3/4)') + ! Finish MPAS framework initialization, including ! the allocation of all blocks and fields managed by MPAS. call mpas_dynamical_core % init_phase3(num_advected, pio_file) + call dyn_debug_print(debugout_info, 'Initializing MPAS mesh variables') + ! Read time-invariant (e.g., grid/mesh) variables. call mpas_dynamical_core % read_write_stream(pio_file, 'r', 'invariant') @@ -81,26 +86,30 @@ subroutine model_grid_init() ! Check for consistency in numbers of vertical layers. if (nvertlevels /= pver) then - call dyn_debug_print('Number of vertical layers in CAM-SIMA namelist, pver = ' // stringify([pver])) - call dyn_debug_print('Number of vertical layers in initial file, nvertlevels = ' // stringify([nvertlevels])) + call dyn_debug_print(debugout_debug, 'Number of vertical layers in CAM-SIMA namelist, pver = ' // & + stringify([pver])) + call dyn_debug_print(debugout_debug, 'Number of vertical layers in initial file, nvertlevels = ' // & + stringify([nvertlevels])) call endrun('Numbers of vertical layers mismatch', subname, __LINE__) end if - ! Initialize reference pressure for use by physics. - call dyn_debug_print('Calling init_reference_pressure') + call dyn_debug_print(debugout_info, 'Initializing reference pressure') + ! Initialize reference pressure for use by physics. call init_reference_pressure() - ! Initialize physics grid. - call dyn_debug_print('Calling init_physics_grid') + call dyn_debug_print(debugout_info, 'Initializing physics grid') + ! Initialize physics grid. call init_physics_grid() - ! Initialize CAM-SIMA grid. - call dyn_debug_print('Calling define_cam_grid') + call dyn_debug_print(debugout_info, 'Registering dynamics grid with CAM-SIMA') + ! Registering dynamics grid with CAM-SIMA. call define_cam_grid() + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine model_grid_init !> Initialize reference pressure for use by physics. @@ -109,6 +118,7 @@ subroutine init_reference_pressure() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate use cam_history_support, only: add_vert_coord + use cam_logfile, only: debugout_debug, debugout_verbose use dyn_comp, only: dyn_debug_print, mpas_dynamical_core use dynconst, only: constant_p0 => pref use ref_pres, only: ref_pres_init @@ -136,6 +146,8 @@ subroutine init_reference_pressure() real(kind_r8), pointer :: zw(:) ! CANNOT be safely deallocated because `add_vert_coord` ! just uses pointers to point at it internally. + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(rdzw) nullify(zu) nullify(zw) @@ -185,23 +197,23 @@ subroutine init_reference_pressure() p_ref_mid(:) = 0.5_kind_r8 * (p_ref_int(1:pver) + p_ref_int(2:pverp)) ! Print a nice table of reference height and pressure. - call dyn_debug_print('Reference layer information:') - call dyn_debug_print('----- | -------------- | --------------') - call dyn_debug_print('Index | Height (m) | Pressure (hPa)') + call dyn_debug_print(debugout_verbose, 'Reference layer information:') + call dyn_debug_print(debugout_verbose, '----- | -------------- | --------------') + call dyn_debug_print(debugout_verbose, 'Index | Height (m) | Pressure (hPa)') do k = 1, pver - call dyn_debug_print(repeat('-', 5) // ' | ' // stringify([zw(k)]) // & + call dyn_debug_print(debugout_verbose, repeat('-', 5) // ' | ' // stringify([zw(k)]) // & ' | ' // stringify([p_ref_int(k) / 100.0_kind_r8])) l = len(stringify([k])) - call dyn_debug_print(repeat(' ', 5 - l) // stringify([k]) // ' | ' // stringify([zu(k)]) // & + call dyn_debug_print(debugout_verbose, repeat(' ', 5 - l) // stringify([k]) // ' | ' // stringify([zu(k)]) // & ' | ' // stringify([p_ref_mid(k) / 100.0_kind_r8])) end do k = pverp - call dyn_debug_print(repeat('-', 5) // ' | ' // stringify([zw(k)]) // & + call dyn_debug_print(debugout_verbose, repeat('-', 5) // ' | ' // stringify([zw(k)]) // & ' | ' // stringify([p_ref_int(k) / 100.0_kind_r8])) call ref_pres_init(p_ref_int, p_ref_mid, num_pure_p_lev) @@ -211,6 +223,8 @@ subroutine init_reference_pressure() nullify(zu) nullify(zw) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine init_reference_pressure !> Initialize physics grid in terms of dynamics decomposition. @@ -219,7 +233,8 @@ end subroutine init_reference_pressure subroutine init_physics_grid() ! Module(s) from CAM-SIMA. use cam_abortutils, only: check_allocate - use dyn_comp, only: mpas_dynamical_core, ncells_global, ncells_solve, sphere_radius + use cam_logfile, only: debugout_debug + use dyn_comp, only: dyn_debug_print, mpas_dynamical_core, ncells_global, ncells_solve, sphere_radius use dynconst, only: constant_pi => pi, rad_to_deg use physics_column_type, only: physics_column_t use physics_grid, only: phys_grid_init @@ -237,6 +252,8 @@ subroutine init_physics_grid() real(kind_dyn_mpas), pointer :: loncell(:) ! Cell center longitudes (radians). type(physics_column_t), allocatable :: dyn_column(:) ! Grid and mapping information between global and local indexes. + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(areacell) nullify(indextocellid) nullify(latcell) @@ -297,6 +314,8 @@ subroutine init_physics_grid() deallocate(dyn_column) deallocate(dyn_attribute_name) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine init_physics_grid !> This subroutine defines and registers four variants of dynamics grids in terms of dynamics decomposition. @@ -311,6 +330,7 @@ subroutine define_cam_grid() use cam_abortutils, only: check_allocate use cam_grid_support, only: cam_grid_attribute_register, cam_grid_register, & horiz_coord_create, horiz_coord_t + use cam_logfile, only: debugout_debug, debugout_verbose use dyn_comp, only: dyn_debug_print, mpas_dynamical_core, & ncells_global, nedges_global, nvertices_global, & ncells_solve, nedges_solve, nvertices_solve, & @@ -353,6 +373,8 @@ subroutine define_cam_grid() ! just uses pointers to point at it internally. type(horiz_coord_t), pointer :: lon_coord + call dyn_debug_print(debugout_debug, subname // ' entered') + nullify(indextocellid, indextoedgeid, indextovertexid) nullify(areacell) nullify(latcell, loncell) @@ -397,7 +419,8 @@ subroutine define_cam_grid() global_grid_map(3, i) = global_grid_index(i) end do - call dyn_debug_print('Registering grid "mpas_cell" with id ' // stringify([dyn_grid_id('mpas_cell')])) + call dyn_debug_print(debugout_verbose, 'Registering dynamics grid "mpas_cell" with id ' // & + stringify([dyn_grid_id('mpas_cell')])) call cam_grid_register('mpas_cell', dyn_grid_id('mpas_cell'), lat_coord, lon_coord, global_grid_map, & unstruct=.true., block_indexed=.false.) @@ -418,7 +441,8 @@ subroutine define_cam_grid() lon_coord => horiz_coord_create('lon', 'ncol', ncells_global, 'longitude', 'degrees_east', & 1, ncells_solve, real(loncell, kind_r8) * rad_to_deg, map=global_grid_index) - call dyn_debug_print('Registering grid "cam_cell" with id ' // stringify([dyn_grid_id('cam_cell')])) + call dyn_debug_print(debugout_verbose, 'Registering dynamics grid "cam_cell" with id ' // & + stringify([dyn_grid_id('cam_cell')])) call cam_grid_register('cam_cell', dyn_grid_id('cam_cell'), lat_coord, lon_coord, global_grid_map, & unstruct=.true., block_indexed=.false.) @@ -456,7 +480,8 @@ subroutine define_cam_grid() global_grid_map(3, i) = global_grid_index(i) end do - call dyn_debug_print('Registering grid "mpas_edge" with id ' // stringify([dyn_grid_id('mpas_edge')])) + call dyn_debug_print(debugout_verbose, 'Registering dynamics grid "mpas_edge" with id ' // & + stringify([dyn_grid_id('mpas_edge')])) call cam_grid_register('mpas_edge', dyn_grid_id('mpas_edge'), lat_coord, lon_coord, global_grid_map, & unstruct=.true., block_indexed=.false.) @@ -494,7 +519,8 @@ subroutine define_cam_grid() global_grid_map(3, i) = global_grid_index(i) end do - call dyn_debug_print('Registering grid "mpas_vertex" with id ' // stringify([dyn_grid_id('mpas_vertex')])) + call dyn_debug_print(debugout_verbose, 'Registering dynamics grid "mpas_vertex" with id ' // & + stringify([dyn_grid_id('mpas_vertex')])) call cam_grid_register('mpas_vertex', dyn_grid_id('mpas_vertex'), lat_coord, lon_coord, global_grid_map, & unstruct=.true., block_indexed=.false.) @@ -505,6 +531,8 @@ subroutine define_cam_grid() deallocate(global_grid_index) nullify(global_grid_index, global_grid_map) nullify(lat_coord, lon_coord) + + call dyn_debug_print(debugout_debug, subname // ' completed') end subroutine define_cam_grid !> Helper function for returning grid id given its name.