diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F index 7d5b8c01ee..05f2b03e3c 100644 --- a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F +++ b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F @@ -20,7 +20,6 @@ module atm_time_integration use mpas_timer #ifdef DO_PHYSICS - use mpas_atmphys_driver_microphysics use mpas_atmphys_todynamics use mpas_atmphys_utilities #endif @@ -165,8 +164,6 @@ subroutine atm_srk3(domain, dt, itimestep) logical, pointer :: config_positive_definite logical, pointer :: config_monotonic real (kind=RKIND), pointer :: config_dt - character (len=StrKIND), pointer :: config_microp_scheme - character (len=StrKIND), pointer :: config_convection_scheme integer, pointer :: num_scalars, index_qv, nCells, nCellsSolve, nEdges, nEdgesSolve, nVertices, nVerticesSolve, nVertLevels @@ -198,8 +195,6 @@ subroutine atm_srk3(domain, dt, itimestep) real (kind=RKIND), dimension(:,:), pointer :: u, uReconstructZonal, uReconstructMeridional, uReconstructX, uReconstructY, uReconstructZ real (kind=RKIND), dimension(:,:,:), pointer :: scalars, scalars_1, scalars_2 - real (kind=RKIND), dimension(:,:), pointer :: rqvdynten - logical, parameter :: debug = .false. @@ -212,8 +207,6 @@ subroutine atm_srk3(domain, dt, itimestep) call mpas_pool_get_config(domain % blocklist % configs, 'config_positive_definite', config_positive_definite) call mpas_pool_get_config(domain % blocklist % configs, 'config_monotonic', config_monotonic) call mpas_pool_get_config(domain % blocklist % configs, 'config_dt', config_dt) - call mpas_pool_get_config(domain % blocklist % configs, 'config_microp_scheme', config_microp_scheme) - call mpas_pool_get_config(domain % blocklist % configs, 'config_convection_scheme', config_convection_scheme) call mpas_pool_get_config(domain % blocklist % configs, 'config_IAU_option', config_IAU_option) ! config variables for dynamics-transport splitting, WCS 18 November 2014 @@ -1249,64 +1242,6 @@ subroutine atm_srk3(domain, dt, itimestep) block => block % next end do - ! - ! call to parameterizations of cloud microphysics. calculation of the tendency of water vapor to horizontal and - ! vertical advection needed for the Tiedtke parameterization of convection. - ! - -#ifdef DO_PHYSICS - block => domain % blocklist - do while(associated(block)) - - call mpas_pool_get_subpool(block % structs, 'mesh', mesh) - call mpas_pool_get_subpool(block % structs, 'state', state) - call mpas_pool_get_subpool(block % structs, 'diag', diag) - call mpas_pool_get_subpool(block % structs, 'diag_physics', diag_physics) - call mpas_pool_get_subpool(block % structs, 'tend_physics', tend_physics) - call mpas_pool_get_subpool(block % structs, 'tend', tend) - call mpas_pool_get_array(state, 'scalars', scalars_1, 1) - call mpas_pool_get_array(state, 'scalars', scalars_2, 2) - call mpas_pool_get_dimension(state, 'index_qv', index_qv) - - call mpas_pool_get_dimension(block % dimensions, 'nThreads', nThreads) - - call mpas_pool_get_dimension(block % dimensions, 'cellSolveThreadStart', cellSolveThreadStart) - call mpas_pool_get_dimension(block % dimensions, 'cellSolveThreadEnd', cellSolveThreadEnd) - - if(config_convection_scheme == 'cu_grell_freitas' .or. & - config_convection_scheme == 'cu_tiedtke' .or. & - config_convection_scheme == 'cu_ntiedtke') then - - call mpas_pool_get_array(tend_physics, 'rqvdynten', rqvdynten) - - !NOTE: The calculation of the tendency due to horizontal and vertical advection for the water vapor mixing ratio - !requires that the subroutine atm_advance_scalars_mono was called on the third Runge Kutta step, so that a halo - !update for the scalars at time_levs(1) is applied. A halo update for the scalars at time_levs(2) is done above. - if (config_monotonic) then - rqvdynten(:,:) = ( scalars_2(index_qv,:,:) - scalars_1(index_qv,:,:) ) / config_dt - else - rqvdynten(:,:) = 0._RKIND - end if - end if - - !simply set to zero negative mixing ratios of different water species (for now): - where ( scalars_2(:,:,:) < 0.0) & - scalars_2(:,:,:) = 0.0 - - !call microphysics schemes: - if (trim(config_microp_scheme) /= 'off') then - call mpas_timer_start('microphysics') -!$OMP PARALLEL DO - do thread=1,nThreads - call driver_microphysics ( block % configs, mesh, state, 2, diag, diag_physics, tend, itimestep, & - cellSolveThreadStart(thread), cellSolveThreadEnd(thread)) - end do -!$OMP END PARALLEL DO - call mpas_timer_stop('microphysics') - end if - block => block % next - end do - ! ! Note: A halo exchange for 'exner' here as well as at the end of ! the first (n-1) dynamics subcycles can substitute for the exchange at @@ -1314,9 +1249,6 @@ subroutine atm_srk3(domain, dt, itimestep) ! and at the end of dynamics subcycles may in future allow for aggregation ! of the 'exner' exchange with other exchanges. ! -#endif - - call summarize_timestep(domain) end subroutine atm_srk3 diff --git a/src/core_atmosphere/mpas_atm_core.F b/src/core_atmosphere/mpas_atm_core.F index 6f2a4353fa..66fbae6133 100644 --- a/src/core_atmosphere/mpas_atm_core.F +++ b/src/core_atmosphere/mpas_atm_core.F @@ -825,7 +825,7 @@ subroutine atm_do_timestep(domain, dt, itimestep) use atm_time_integration #ifdef DO_PHYSICS use mpas_atmphys_control - use mpas_atmphys_driver + use mpas_atmphys_driver, only : physics_driver_before_dyn, physics_driver_after_dyn use mpas_atmphys_manager use mpas_atmphys_update #endif @@ -852,17 +852,39 @@ subroutine atm_do_timestep(domain, dt, itimestep) call mpas_get_time(curr_time=currTime, dateTimeString=timeStamp, ierr=ierr) - + ! + ! This is a call to all of the slow, column-oriented physics, except for resolved + ! cloud microphysical processes. The microphysics is handled after the dynamics so + ! that any supersaturation caused by various dynamical processes (advection, etc) + ! may be removed before the end of each time step. + ! #ifdef DO_PHYSICS !proceed with physics if moist_physics is set to true: if(moist_physics) then call physics_timetracker(domain,dt,clock,itimestep,xtime_s) - call physics_driver(domain,itimestep,xtime_s) + call physics_driver_before_dyn(domain,itimestep,xtime_s) endif #endif call atm_timestep(domain, dt, timeStamp, itimestep) + ! + ! This part is the microphysics, which needs to happen after the dynamics + ! time stepping. + ! +#ifdef DO_PHYSICS + !proceed with physics if moist_physics is set to true: + if(moist_physics) then + call physics_driver_after_dyn(domain,itimestep,xtime_s) + endif +#endif + + ! + ! This is the summary from the physics before the dynamics solver, + ! the accumulated processing of the dynamics, and the final physics + ! call after the dynamics completes. + call summarize_timestep(domain) + end subroutine atm_do_timestep diff --git a/src/core_atmosphere/physics/Makefile b/src/core_atmosphere/physics/Makefile index e8cb03f6f5..3b1dd144eb 100644 --- a/src/core_atmosphere/physics/Makefile +++ b/src/core_atmosphere/physics/Makefile @@ -79,6 +79,7 @@ mpas_atmphys_driver.o: \ mpas_atmphys_driver_radiation_sw.o \ mpas_atmphys_driver_sfclayer.o \ mpas_atmphys_driver_oml.o \ + mpas_atmphys_driver_microphysics.o \ mpas_atmphys_constants.o \ mpas_atmphys_interface.o \ mpas_atmphys_update.o \ diff --git a/src/core_atmosphere/physics/mpas_atmphys_driver.F b/src/core_atmosphere/physics/mpas_atmphys_driver.F index b120d0ccc4..d846b449d3 100644 --- a/src/core_atmosphere/physics/mpas_atmphys_driver.F +++ b/src/core_atmosphere/physics/mpas_atmphys_driver.F @@ -10,15 +10,20 @@ module mpas_atmphys_driver use mpas_kind_types use mpas_pool_routines + ! For physics_driver_before_dyn use mpas_atmphys_driver_cloudiness use mpas_atmphys_driver_convection use mpas_atmphys_driver_gwdo use mpas_atmphys_driver_lsm use mpas_atmphys_driver_pbl - use mpas_atmphys_driver_radiation_lw + use mpas_atmphys_driver_radiation_lw use mpas_atmphys_driver_radiation_sw use mpas_atmphys_driver_sfclayer use mpas_atmphys_driver_oml + + ! For physics_driver_after_dyn + use mpas_atmphys_driver_microphysics + use mpas_atmphys_constants use mpas_atmphys_interface use mpas_atmphys_update @@ -27,50 +32,18 @@ module mpas_atmphys_driver implicit none private - public:: physics_driver + public:: physics_driver_before_dyn,physics_driver_after_dyn !MPAS top physics driver. !Laura D. Fowler (send comments to laura@ucar.edu). !2013-05-01. ! -! subroutine physics_driver is the top physics driver from which separate drivers for all physics +! subroutine physics_driver_before_dyn is the top physics driver from which separate drivers for all physics ! parameterizations, except cloud microphysics parameterizations are called. ! -! subroutines called in mpas_atmphys_driver: -! ------------------------------------------ -! allocate_forall_physics : allocate local arrays defining atmospheric soundings (pressure,..) -! allocate_cloudiness : allocate all local arrays used in driver_cloudiness. -! allocate_convection : allocate all local arrays used in driver_convection. -! allocate_gwdo : allocate all local arrays used in driver_gwdo. -! allocate_lsm : allocate all local arrays used in driver_lsm. -! allocate_pbl : allocate all local arrays used in driver_pbl. -! allocate_radiation_lw : allocate all local arrays used in driver_radiation_lw. -! allocate_radiation_sw : allocate all local arrays used in driver_radiation_sw. -! allocate_sfclayer : allocate all local arrays used in driver_sfclayer. -! -! deallocate_forall_physics : deallocate local arrays defining atmospheric soundings. -! deallocate_cloudiness : dedeallocate all local arrays used in driver_cloudiness. -! deallocate_convection : deallocate all local arrays used in driver_convection. -! deallocate_gwdo : deallocate all local arrays used in driver_gwdo. -! deallocate_lsm : deallocate all local arrays used in driver_lsm. -! deallocate_pbl : deallocate all local arrays used in driver_pbl. -! deallocate_radiation_lw : deallocate all local arrays used in driver_radiation_lw. -! deallocate_radiation_sw : deallocate all local arrays used in driver_radiation_sw. -! deallocate_sfclayer : deallocate all local arrays used in driver_sfclayer. -! -! MPAS_to_physics : -! driver_cloudiness : driver for parameterization of fractional cloudiness. -! driver_convection : driver for parameterization of convection. -! driver_gwdo : driver for parameterization of gravity wave drag over orography. -! driver_lsm : driver for land-surface scheme. -! driver_pbl : driver for planetary boundary layer scheme. -! driver_radiation_sw : driver for short wave radiation schemes. -! driver_radiation_lw : driver for long wave radiation schemes. -! driver_sfclayer : driver for surface layer scheme. -! update_convection_step1 : updates lifetime of deep convective clouds in Kain-Fritsch scheme. -! update_convection_step2 : updates accumulated precipitation output from convection schemes. -! update_radiation_diagnostics: updates accumualted radiation diagnostics from radiation schemes. +! subroutine physics_driver_after_dyn is the top physics driver from which a separate driver for +! cloud microphysics is called. ! ! add-ons and modifications to sourcecode: ! ---------------------------------------- @@ -97,13 +70,17 @@ module mpas_atmphys_driver ! * modified call to driver_cloudiness to accomodate the calculation of the cloud fraction with the Thompson ! cloud microphysics scheme. ! Laura D. Fowler (laura@ucar.edu) / 2016-06-04. +! * added a second subroutine so now all of the physics is called from this driver module; the original +! routine has an added suffix _before_dyn, and the new routine is named similarly: _after_dyn processing; the +! cloud microphysics processing is handled in the _after_dyn routine. +! Dave Gill (gill@ucar.edu) / 2017-09-14. contains !================================================================================================================= - subroutine physics_driver(domain,itimestep,xtime_s) + subroutine physics_driver_before_dyn(domain,itimestep,xtime_s) !================================================================================================================= !input arguments: @@ -146,11 +123,46 @@ subroutine physics_driver(domain,itimestep,xtime_s) integer,pointer:: nThreads integer,dimension(:),pointer:: cellSolveThreadStart, cellSolveThreadEnd +! subroutines called +! ------------------------------------------ +! allocate_forall_physics : allocate local arrays defining atmospheric soundings (pressure,..) +! allocate_cloudiness : allocate all local arrays used in driver_cloudiness. +! allocate_convection : allocate all local arrays used in driver_convection. +! allocate_gwdo : allocate all local arrays used in driver_gwdo. +! allocate_lsm : allocate all local arrays used in driver_lsm. +! allocate_pbl : allocate all local arrays used in driver_pbl. +! allocate_radiation_lw : allocate all local arrays used in driver_radiation_lw. +! allocate_radiation_sw : allocate all local arrays used in driver_radiation_sw. +! allocate_sfclayer : allocate all local arrays used in driver_sfclayer. +! +! deallocate_forall_physics : deallocate local arrays defining atmospheric soundings. +! deallocate_cloudiness : dedeallocate all local arrays used in driver_cloudiness. +! deallocate_convection : deallocate all local arrays used in driver_convection. +! deallocate_gwdo : deallocate all local arrays used in driver_gwdo. +! deallocate_lsm : deallocate all local arrays used in driver_lsm. +! deallocate_pbl : deallocate all local arrays used in driver_pbl. +! deallocate_radiation_lw : deallocate all local arrays used in driver_radiation_lw. +! deallocate_radiation_sw : deallocate all local arrays used in driver_radiation_sw. +! deallocate_sfclayer : deallocate all local arrays used in driver_sfclayer. +! +! MPAS_to_physics : conversion of input "MPAS" variables to "physics" variables. +! driver_cloudiness : driver for parameterization of fractional cloudiness. +! driver_convection : driver for parameterization of convection. +! driver_gwdo : driver for parameterization of gravity wave drag over orography. +! driver_lsm : driver for land-surface scheme. +! driver_pbl : driver for planetary boundary layer scheme. +! driver_radiation_sw : driver for short wave radiation schemes. +! driver_radiation_lw : driver for long wave radiation schemes. +! driver_sfclayer : driver for surface layer scheme. +! update_convection_step1 : updates lifetime of deep convective clouds in Kain-Fritsch scheme. +! update_convection_step2 : updates accumulated precipitation output from convection schemes. +! update_radiation_diagnostics: updates accumualted radiation diagnostics from radiation schemes. + !================================================================================================================= !call mpas_log_write('') -!call mpas_log_write('--- enter subroutine mpas_atmphys_driver:') +!call mpas_log_write('--- enter subroutine physics_driver_before_dyn:') - call mpas_timer_start('physics driver') + call mpas_timer_start('physics_driver_before_dyn') call mpas_pool_get_config(domain%configs,'config_convection_scheme',config_convection_scheme) call mpas_pool_get_config(domain%configs,'config_gwdo_scheme' ,config_gwdo_scheme ) @@ -340,12 +352,132 @@ subroutine physics_driver(domain,itimestep,xtime_s) endif - call mpas_timer_stop('physics driver') + call mpas_timer_stop('physics_driver_before_dyn') + +!call mpas_log_write('--- exit subroutine physics_driver_before_dyn:') +!call mpas_log_write('') + + end subroutine physics_driver_before_dyn + + +!================================================================================================================= + subroutine physics_driver_after_dyn(domain,itimestep,xtime_s) +!================================================================================================================= + +!input arguments: + integer,intent(in):: itimestep + real(kind=RKIND),intent(in):: xtime_s + +!inout arguments: + type(domain_type),intent(inout):: domain + +!local pointers: + type(mpas_pool_type) ,pointer :: configs ,& + mesh ,& + state ,& + diag ,& + diag_physics ,& + tend_physics ,& + tend + + real (kind=RKIND), dimension(:,:,:), pointer :: scalars_1 ,& + scalars_2 + real (kind=RKIND), dimension(:,:) , pointer :: rqvdynten + real (kind=RKIND) , pointer :: config_dt + + logical , pointer :: config_monotonic + + character(len=StrKIND) , pointer :: config_microp_scheme ,& + config_convection_scheme + + integer , pointer :: nThreads ,& + index_qv + integer,dimension(:) , pointer :: cellSolveThreadStart ,& + cellSolveThreadEnd + + type(block_type) , pointer :: block + +!local variables: + integer:: time_lev + integer:: thread + +!================================================================================================================= +!call mpas_log_write('') +!call mpas_log_write('--- enter subroutine physics_driver_after_dyn:') + + call mpas_timer_start('physics_driver_after_dyn') + + call mpas_pool_get_config(domain%configs,'config_microp_scheme' ,config_microp_scheme ) + + if(config_microp_scheme .ne. 'off') then + + call mpas_pool_get_config(domain%configs,'config_convection_scheme',config_convection_scheme) + call mpas_pool_get_config(domain%configs,'config_dt' ,config_dt ) + call mpas_pool_get_config(domain%configs,'config_monotonic' ,config_monotonic ) + + block => domain % blocklist + do while(associated(block)) + + call mpas_pool_get_subpool(block%structs,'mesh' ,mesh ) + call mpas_pool_get_subpool(block%structs,'state' ,state ) + call mpas_pool_get_subpool(block%structs,'diag' ,diag ) + call mpas_pool_get_subpool(block%structs,'diag_physics',diag_physics) + call mpas_pool_get_subpool(block%structs,'tend_physics',tend_physics) + call mpas_pool_get_subpool(block%structs,'tend' ,tend ) + + call mpas_pool_get_array(state, 'scalars', scalars_1, 1) + call mpas_pool_get_array(state, 'scalars', scalars_2, 2) + + call mpas_pool_get_dimension(state ,'index_qv' ,index_qv ) + call mpas_pool_get_dimension(block%dimensions,'nThreads' ,nThreads ) + call mpas_pool_get_dimension(block%dimensions,'cellSolveThreadStart',cellSolveThreadStart) + call mpas_pool_get_dimension(block%dimensions,'cellSolveThreadEnd' ,cellSolveThreadEnd ) + + !physics prep step: + time_lev = 1 + + if(config_convection_scheme == 'cu_grell_freitas' .or. & + config_convection_scheme == 'cu_tiedtke' .or. & + config_convection_scheme == 'cu_ntiedtke') then + + call mpas_pool_get_array(tend_physics, 'rqvdynten', rqvdynten) + + !NOTE: The calculation of the tendency due to horizontal and vertical advection for the water vapor mixing ratio + !requires that the subroutine atm_advance_scalars_mono was called on the third Runge Kutta step, so that a halo + !update for the scalars at time_levs(1) is applied. A halo update for the scalars at time_levs(2) is done above. + if (config_monotonic) then + rqvdynten(:,:) = ( scalars_2(index_qv,:,:) - scalars_1(index_qv,:,:) ) / config_dt + else + rqvdynten(:,:) = 0._RKIND + end if + end if + + !simply set to zero negative mixing ratios of different water species (for now): + where ( scalars_2(:,:,:) < 0.0) & + scalars_2(:,:,:) = 0.0 + + !call microphysics schemes: + if (trim(config_microp_scheme) /= 'off') then + call mpas_timer_start('microphysics') +!$OMP PARALLEL DO + do thread=1,nThreads + call driver_microphysics ( block%configs, mesh, state, 2, diag, diag_physics, tend, itimestep, & + cellSolveThreadStart(thread), cellSolveThreadEnd(thread)) + end do +!$OMP END PARALLEL DO + call mpas_timer_stop('microphysics') + end if + block => block % next + end do + + endif + + call mpas_timer_stop('physics_driver_after_dyn') -!call mpas_log_write('--- enter subroutine mpas_atmphys_driver:') +!call mpas_log_write('--- exit subroutine physics_driver_after_dyn:') !call mpas_log_write('') - end subroutine physics_driver + end subroutine physics_driver_after_dyn !================================================================================================================= end module mpas_atmphys_driver