diff --git a/CONTRIBUTORS.md b/CONTRIBUTORS.md index c64b5adb..69a39ad0 100644 --- a/CONTRIBUTORS.md +++ b/CONTRIBUTORS.md @@ -9,3 +9,4 @@ | ScottWales | Scott Wales | Bureau of Meteorology | 2026-02-16 | | t00sa | Sam Clarke-Green | Met Office | 2026-02-27 | | Pierre-siddall | Pierre Siddall| Met Office | 2026-02-06 | +| RichardEssery | Richard Essery | University of Edinburgh | 2026-04-18 | diff --git a/src/control/shared/fluxes.F90 b/src/control/shared/fluxes.F90 index 478b4909..45d674e5 100644 --- a/src/control/shared/fluxes.F90 +++ b/src/control/shared/fluxes.F90 @@ -95,6 +95,8 @@ MODULE fluxes_mod ! Surface FQW for sea-ice REAL(KIND=real_jlslsm), ALLOCATABLE :: fsmc_pft(:,:) ! Moisture availability factor. + REAL(KIND=real_jlslsm), ALLOCATABLE :: fsnow(:,:) + ! Snow cover fraction on tiles REAL(KIND=real_jlslsm), ALLOCATABLE :: ftl_sicat(:,:,:) ! Surface FTL for sea-ice REAL(KIND=real_jlslsm), ALLOCATABLE :: ftl_surft(:,:) @@ -215,6 +217,8 @@ MODULE fluxes_mod ! Surface FQW for sea-ice REAL(KIND=real_jlslsm), POINTER :: fsmc_pft(:,:) ! Moisture availability factor. + REAL(KIND=real_jlslsm), POINTER :: fsnow(:,:) + ! Snow cover fraction on tiles REAL(KIND=real_jlslsm), POINTER :: ftl_sicat(:,:,:) ! Surface FTL for sea-ice REAL(KIND=real_jlslsm), POINTER :: ftl_surft(:,:) @@ -340,6 +344,7 @@ SUBROUTINE fluxes_alloc(land_pts, t_i_length, t_j_length, & ALLOCATE(fluxes_data%sub_surf_roff_gb(land_pts)) ALLOCATE(fluxes_data%surf_roff_gb(land_pts)) ALLOCATE(fluxes_data%alb_surft(land_pts,nsurft,4) ) +ALLOCATE(fluxes_data%fsnow(land_pts,nsurft)) ALLOCATE(fluxes_data%tstar_ij(t_i_length,t_j_length)) ALLOCATE(fluxes_data%e_sea_ij(t_i_length,t_j_length)) ALLOCATE(fluxes_data%fsmc_pft(land_pts,npft)) @@ -383,6 +388,7 @@ SUBROUTINE fluxes_alloc(land_pts, t_i_length, t_j_length, & fluxes_data%sub_surf_roff_gb(:) = 0.0 fluxes_data%surf_roff_gb(:) = 0.0 fluxes_data%alb_surft(:,:,:) = 0.0 +fluxes_data%fsnow(:,:) = 0.0 fluxes_data%tstar_ij(:,:) = 0.0 fluxes_data%e_sea_ij(:,:) = 0.0 fluxes_data%fsmc_pft(:,:) = 0.0 @@ -459,6 +465,7 @@ SUBROUTINE fluxes_dealloc(fluxes_data) DEALLOCATE(fluxes_data%sub_surf_roff_gb) DEALLOCATE(fluxes_data%surf_roff_gb) DEALLOCATE(fluxes_data%alb_surft) +DEALLOCATE(fluxes_data%fsnow) DEALLOCATE(fluxes_data%tstar_ij) DEALLOCATE(fluxes_data%e_sea_ij) DEALLOCATE(fluxes_data%fsmc_pft) @@ -541,6 +548,7 @@ SUBROUTINE fluxes_assoc(fluxes,fluxes_data) fluxes%sub_surf_roff_gb => fluxes_data%sub_surf_roff_gb fluxes%surf_roff_gb => fluxes_data%surf_roff_gb fluxes%alb_surft => fluxes_data%alb_surft +fluxes%fsnow => fluxes_data%fsnow fluxes%tstar_ij => fluxes_data%tstar_ij fluxes%e_sea_ij => fluxes_data%e_sea_ij fluxes%fsmc_pft => fluxes_data%fsmc_pft @@ -619,6 +627,7 @@ SUBROUTINE fluxes_nullify(fluxes) NULLIFY(fluxes%sub_surf_roff_gb) NULLIFY(fluxes%surf_roff_gb) NULLIFY(fluxes%alb_surft) +NULLIFY(fluxes%fsnow) NULLIFY(fluxes%tstar_ij) NULLIFY(fluxes%e_sea_ij) NULLIFY(fluxes%fsmc_pft) diff --git a/src/control/shared/jules_radiation_mod.F90 b/src/control/shared/jules_radiation_mod.F90 index 160f931c..0bf6232f 100644 --- a/src/control/shared/jules_radiation_mod.F90 +++ b/src/control/shared/jules_radiation_mod.F90 @@ -116,6 +116,7 @@ SUBROUTINE check_jules_radiation() USE ereport_mod, ONLY: ereport USE jules_print_mgr, ONLY: jules_message +USE jules_snow_mod, ONLY: i_snow_tile !----------------------------------------------------------------------------- ! Description: @@ -157,6 +158,17 @@ SUBROUTINE check_jules_radiation() CALL ereport("check_jules_radiation", errorstatus, jules_message) END IF +! Required options if a separate snow tile is being used +IF ( ANY(i_snow_tile == 1) ) THEN + IF ( (.NOT. l_spec_albedo) .OR. l_embedded_snow .OR. & + (.NOT. l_snow_albedo) ) THEN + errorstatus = 1003 + CALL ereport("check_jules_radiation", errorstatus, & + "l_spec_albedo = T, l_embedded_snow = F and l_snow_albedo = T" // & + "required if ANY i_snow_tile = 1") + END IF +END IF + !Can set the size of rad_nband IF (l_spec_albedo) THEN rad_nband = 2 diff --git a/src/control/shared/jules_snow_mod.F90 b/src/control/shared/jules_snow_mod.F90 index 8b35d2ab..6c5306ca 100644 --- a/src/control/shared/jules_snow_mod.F90 +++ b/src/control/shared/jules_snow_mod.F90 @@ -71,6 +71,12 @@ MODULE jules_snow_mod INTEGER :: & i_snow_cond_parm = 0 ! Parametrization scheme for snow conductivity + +INTEGER :: & + i_snow_tile(nsurft_max) + ! Tiles without (0) or with (1) separate energy balance + ! for snow. Must be 0 on ice tiles +DATA i_snow_tile / nsurft_max * 0 / !----------------------------------------------------------------------------- ! Parametrization of the rate of growth of snow grains @@ -307,7 +313,7 @@ MODULE jules_snow_mod NAMELIST / jules_snow/ & ! Switches nsmax, l_snowdep_surf, l_rho_snow_corr, frac_snow_subl_melt, & - graupel_options, & + graupel_options,i_snow_tile, & ! Equitemperature metamorphism l_et_metamorph, a_snow_et, b_snow_et, c_snow_et, rho_snow_et_crit, & ! Thermal conductivity of snow @@ -333,7 +339,7 @@ SUBROUTINE check_jules_snow() USE water_constants_mod, ONLY: tm -USE jules_surface_types_mod, ONLY: npft +USE jules_surface_types_mod, ONLY: ice, npft, ntype USE jules_surface_mod, ONLY: l_aggregate USE jules_vegetation_mod, ONLY: can_model @@ -388,6 +394,25 @@ SUBROUTINE check_jules_snow() END IF END IF +! Check for inconsistencies if a separate snow tile is used +IF ( ANY(i_snow_tile == 1) ) THEN +! Tiles cannot be aggregated + IF ( l_aggregate ) THEN + CALL ereport(RoutineName, errorstatus, & + 'l_aggregate = F required if ANY i_snow_tile = 1') + END IF +! The last surface type must be an ice tile + IF ( ice .NE. ntype ) THEN + CALL ereport(RoutineName, errorstatus, & + 'ice = ntype required if ANY i_snow_tile = 1') + END IF +! The ice tile cannot be selected + IF ( i_snow_tile(ntype) == 1 ) THEN + CALL ereport(RoutineName, errorstatus, & + 'i_snow_tile(ntype) = 1 is not permitted') + END IF +END IF + END SUBROUTINE check_jules_snow #endif diff --git a/src/control/shared/surf_couple_explicit_mod.F90 b/src/control/shared/surf_couple_explicit_mod.F90 index 5052725b..6f922834 100644 --- a/src/control/shared/surf_couple_explicit_mod.F90 +++ b/src/control/shared/surf_couple_explicit_mod.F90 @@ -554,8 +554,8 @@ SUBROUTINE surf_couple_explicit( & !IN soil/vegetation/land surface data : ainfo%land_index,nsurft,sm_levels,progs%canopy_surft,psparms%catch_surft, & psparms%catch_snow_surft, psparms%hcon_soilt,jules_vars%ho2r2_orog_gb, & - flandg_tmp, & - progs%snow_surft,jules_vars%sil_orog_land_gb,psparms%smvccl_soilt, & + flandg_tmp, fluxes%fsnow, & + jules_vars%sil_orog_land_gb,psparms%smvccl_soilt, & psparms%smvcst_soilt,psparms%smvcwt_soilt, & psparms%sthf_soilt, psparms%sthu_soilt,psparms%z0_surft, & psparms%z0h_bare_surft, psparms%z0m_soil_gb, & @@ -636,9 +636,11 @@ SUBROUTINE surf_couple_explicit( & urban_param%wrr_gb, & !Fluxes (IN OUT) fluxes%anthrop_heat_surft, & - !prognostics (IN) + !snow prognostics (IN) progs%nsnow_surft, progs%sice_surft, progs%sliq_surft, & - progs%snowdepth_surft, progs%tsnow_surft, progs%ds_surft, & + progs%tsnow_surft, progs%ds_surft, & + !snow prognostics (IN OUT) + progs%snowdepth_surft, progs%snow_surft, & !c_elevate (OUT) jules_vars%surf_hgt_surft, jules_vars%lw_down_elevcorr_surft, & !jules_vars_mod (OUT) diff --git a/src/control/shared/surf_couple_extra_mod.F90 b/src/control/shared/surf_couple_extra_mod.F90 index ca35dc00..2fb70b3a 100644 --- a/src/control/shared/surf_couple_extra_mod.F90 +++ b/src/control/shared/surf_couple_extra_mod.F90 @@ -703,9 +703,9 @@ SUBROUTINE surf_couple_extra( & fluxes%melt_surft,fluxes%snowinc_surft, & progs%smcl_soilt(:,:,1),psparms%sthf_soilt(:,:,1), & fluxes%surf_htf_surft, progs%t_soil_soilt(:,:,1), & - progs%tsurf_elev_surft, progs%tstar_surft, & + progs%tsurf_elev_surft, progs%tstar_surft, psparms%z0_surft, & psparms%smvcst_soilt(:,:,1), con_snow_wtrac, & - wtrac_jls%ei_surft, progs%rgrain_surft, & + wtrac_jls%ei_surft, fluxes%fsnow, progs%rgrain_surft, & progs%rgrainl_surft, progs%rho_snow_grnd_surft, & progs%sice_surft,progs%sliq_surft,progs%snow_grnd_surft, & progs%snow_surft, progs%snowdepth_surft, progs%tsnow_surft, & diff --git a/src/control/shared/surf_couple_implicit_mod.F90 b/src/control/shared/surf_couple_implicit_mod.F90 index 461271a4..b917440d 100644 --- a/src/control/shared/surf_couple_implicit_mod.F90 +++ b/src/control/shared/surf_couple_implicit_mod.F90 @@ -446,7 +446,7 @@ SUBROUTINE surf_couple_implicit( & alpha1,ashtf_surft, & jules_vars%dtrdz_charney_grid_1_ij,fracaero_t,fracaero_s, & resfs,resft, rhokh_surft, & - fluxes%emis_surft,progs%snow_surft,dtstar_surft, & + fluxes%emis_surft,fluxes%fsnow,progs%snow_surft,dtstar_surft, & ! INOUT data : progs%tstar_surft,fluxes%fqw_surft,fqw_1,ftl_1,fluxes%ftl_surft, & sf_diag, & diff --git a/src/control/shared/surf_couple_radiation_mod.F90 b/src/control/shared/surf_couple_radiation_mod.F90 index 13422ff8..b6260df8 100644 --- a/src/control/shared/surf_couple_radiation_mod.F90 +++ b/src/control/shared/surf_couple_radiation_mod.F90 @@ -248,7 +248,7 @@ SUBROUTINE surf_couple_radiation( & progs%rgrain_surft, snow_surft, progs%tstar_surft, psparms%z0_surft, & ainfo%frac_surft, & !INTENT(OUT) - fluxes%alb_surft,albobs_sc_ij,fluxes%land_albedo_ij, & + fluxes%alb_surft,albobs_sc_ij,fluxes%fsnow,fluxes%land_albedo_ij, & !New arguments replacing USE statements !jules_vars_mod (IN OUT) jules_vars%albobs_scaling_surft, & diff --git a/src/control/shared/tilepts_jls.F90 b/src/control/shared/tilepts_jls.F90 index ebd056c8..9f4c1e53 100644 --- a/src/control/shared/tilepts_jls.F90 +++ b/src/control/shared/tilepts_jls.F90 @@ -16,6 +16,8 @@ MODULE tilepts_mod SUBROUTINE tilepts(land_pts,frac,surft_pts,surft_index,l_lice_point, & l_lice_surft) +USE jules_snow_mod, ONLY: i_snow_tile + USE jules_surface_mod, ONLY: all_tiles, l_elev_land_ice USE jules_surface_types_mod, ONLY: ntype, ice, elev_ice, elev_rock @@ -41,6 +43,9 @@ SUBROUTINE tilepts(land_pts,frac,surft_pts,surft_index,l_lice_point, & LOGICAL, INTENT(IN) :: l_lice_point(land_pts) LOGICAL, INTENT(OUT) :: l_lice_surft(ntype) +LOGICAL :: snow_tile(land_pts) + ! Indicates if a separate snow tile is used on a land point + LOGICAL :: use_tile ! Indicates if we will model the tile for the current ! land point @@ -62,6 +67,7 @@ SUBROUTINE tilepts(land_pts,frac,surft_pts,surft_index,l_lice_point, & !$OMP elev_ice, elev_rock, l_lice_point, l_lice_surft, & !$OMP l_elev_land_ice, surft_pts ) & !$OMP PRIVATE( l, n, c, use_tile ) +snow_tile(:) = .FALSE. DO n = 1,ntype c = 0 l_lice_surft(n) = ( n == ice .OR. ANY(elev_ice == n) ) @@ -75,6 +81,14 @@ SUBROUTINE tilepts(land_pts,frac,surft_pts,surft_index,l_lice_point, & ! as this may be needed for on- or off-line icesheet coupling IF (l_elev_land_ice .AND. l_lice_point(l) .AND. l_lice_surft(n)) & use_tile = .TRUE. + ! Ice tile can be used as a separate snow tile + IF (ANY(i_snow_tile == 1)) THEN + IF (i_snow_tile(n) == 1 .AND. frac(l,n) > 0.0) THEN + snow_tile(l) = .TRUE. + ELSE IF (n == ice .AND. snow_tile(l)) THEN + use_tile = .TRUE. + END IF + END IF ELSE ! with all_tiles we do ! * All tiles except the ice/elevated tiles on non-land-ice points diff --git a/src/io/model_interface/extract_var.inc b/src/io/model_interface/extract_var.inc index 49d27599..e82946c4 100644 --- a/src/io/model_interface/extract_var.inc +++ b/src/io/model_interface/extract_var.inc @@ -76,7 +76,7 @@ USE sf_diags_mod, ONLY: & sf_diag USE jules_snow_mod, ONLY: & - canSnowTile + canSnowTile, i_snow_tile USE overbank_inundation_mod, ONLY: & frac_fplain_lp, frac_fplain_rp @@ -1100,15 +1100,20 @@ CASE ( 'snow_depth_gb' ) cube = map_from_land(cube_land) CASE ( 'snow_frac' ) - ! Sum frac over tiles with snow. workspace_land(:) = 0.0 - IF ( l_aggregate ) THEN - WHERE ( progs%snow_surft(:,1) + progs%snow_grnd_surft(:,1) > EPSILON(1.0) ) - workspace_land(:) = 1.0 - END WHERE - ELSE - workspace_land(:) = SUM(ainfo%frac_surft, 2, progs%snow_surft + & - progs%snow_grnd_surft > EPSILON(1.0)) + IF (ANY(i_snow_tile == 1)) THEN + ! Sum snow cover fractions on tiles + workspace_land(:) = SUM(fluxes%fsnow * ainfo%frac_surft(:,:),2) + ELSE + ! Sum frac over tiles with snow. + IF ( l_aggregate ) THEN + WHERE ( progs%snow_surft(:,1) + progs%snow_grnd_surft(:,1) > EPSILON(1.0) ) + workspace_land(:) = 1.0 + END WHERE + ELSE + workspace_land(:) = SUM(ainfo%frac_surft, 2, progs%snow_surft + & + progs%snow_grnd_surft > EPSILON(1.0)) + END IF END IF cube_land = cube_from_array(workspace_land) cube = map_from_land(cube_land) diff --git a/src/io/model_interface/fluxes.F90 b/src/io/model_interface/fluxes.F90 new file mode 100644 index 00000000..45d674e5 --- /dev/null +++ b/src/io/model_interface/fluxes.F90 @@ -0,0 +1,675 @@ +! *****************************COPYRIGHT******************************* +! (C) Crown copyright Met Office. All rights reserved. +! For further details please refer to the file COPYRIGHT.txt +! which you should have received as part of this distribution. +! *****************************COPYRIGHT******************************* +! +! Module containing surface fluxes. +! +! Code Description: +! Language: FORTRAN 90 +! This code is written to UMDP3 v8.2 programming standards. +! +! Code Owner: Please refer to ModuleLeaders.txt and UM file CodeOwners.txt +! This file belongs in section: Land +! +! Changes to variable names to enable soil tiling +! Anything named _tile is now ambiguous, so the following convention is +! adopted: +! _gb for variables on land points +! _ij for variables with i and j indices +! Surface heterogeneity... +! _surft for surface tiled variables (generally size n_surft) +! _pft for plant functional type surface tiled variables (generally sized n_pft) +! _sicat for sea ice catergories (generally size nice_use) +! Sub-surface heterogeneity... +! _soilt for soil tiled variables + +! Implementation: +! Each variable is declared in both the 'data' TYPE and the 'pointer' type. +! Instances of these types are declared at at high level as required +! This is to facilitate advanced memory management features, which are generally +! not visible in the science code. +! Checklist for adding a new variable: +! -add to data_type +! -add to pointer_type +! -add to the allocate routine, passing in any new dimension sizes required +! by argument (not via USE statement) +! -add to the deallocate routine +! -add to the assoc and nullify routines + +MODULE fluxes_mod + +USE um_types, ONLY: real_jlslsm +USE missing_data_mod, ONLY: rmdi + +IMPLICIT NONE + +TYPE :: fluxes_data_type + + ! T if varying grey emissivity is used for surface tile type + LOGICAL, ALLOCATABLE :: l_emis_surft_set(:) + ! anthrop_heat is required by both the UM and standalone configurations + + REAL(KIND=real_jlslsm), ALLOCATABLE :: anthrop_heat_surft(:,:) + ! Additional heat source on surface tiles used for anthropgenic urban heat + ! source (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: surf_ht_store_surft(:,:) + ! Diagnostic to store values of C*(dT/dt) during calculation of energy + ! balance + REAL(KIND=real_jlslsm), ALLOCATABLE :: sw_sicat(:,:) + ! Net SW on sea ice categories + REAL(KIND=real_jlslsm), ALLOCATABLE :: alb_sicat(:,:,:) + ! Albedo of sea ice categories (ij point, sicat, band- see below) + REAL(KIND=real_jlslsm), ALLOCATABLE :: penabs_rad_frac(:,:,:) + ! Fraction of downward solar that penetrates the sea ice and is absorbed + REAL(KIND=real_jlslsm), ALLOCATABLE :: sw_sea(:) + ! Net SW on open sea + REAL(KIND=real_jlslsm), ALLOCATABLE :: alb_surft(:,:,:) + ! Albedo for surface tiles + ! (:,:,1) direct beam visible + ! (:,:,2) diffuse visible + ! (:,:,3) direct beam near-IR + ! (:,:,4) diffuse near-IR + REAL(KIND=real_jlslsm), ALLOCATABLE :: e_sea_ij(:,:) + ! Evaporation from sea times leads fraction. Zero over land + ! (kg per square metre per sec) + REAL(KIND=real_jlslsm), ALLOCATABLE :: ecan_ij(:,:) + ! Gridbox mean evaporation from canopy/surface store (kg/m2/s) + ! Zero over sea + REAL(KIND=real_jlslsm), ALLOCATABLE :: ecan_surft(:,:) + ! Canopy evaporation from for snow-free land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: ei_ij(:,:) + ! Sublimation from lying snow or sea-ice (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: ei_surft(:,:) + ! EI for land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: esoil_ij_soilt(:,:,:) + ! Surface evapotranspiration from soil moisture store (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: esoil_surft(:,:) + ! ESOIL for snow-free land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: ext_soilt(:,:,:) + ! Extraction of water from each soil layer (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: fqw_surft(:,:) + ! Surface FQW for land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: fqw_sicat(:,:,:) + ! Surface FQW for sea-ice + REAL(KIND=real_jlslsm), ALLOCATABLE :: fsmc_pft(:,:) + ! Moisture availability factor. + REAL(KIND=real_jlslsm), ALLOCATABLE :: fsnow(:,:) + ! Snow cover fraction on tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: ftl_sicat(:,:,:) + ! Surface FTL for sea-ice + REAL(KIND=real_jlslsm), ALLOCATABLE :: ftl_surft(:,:) + ! Surface FTL for land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: h_sea_ij(:,:) + ! Surface sensible heat flux over sea times leads fraction (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: hf_snow_melt_gb(:) + ! Gridbox snowmelt heat flux (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: land_albedo_ij(:,:,:) + ! GBM albedo + ! (:,:,1) direct beam visible + ! (:,:,2) diffuse visible + ! (:,:,3) direct beam near-IR + ! (:,:,4) diffuse near-IR + REAL(KIND=real_jlslsm), ALLOCATABLE :: le_surft(:,:) + ! Surface latent heat flux for land tiles + REAL(KIND=real_jlslsm), ALLOCATABLE :: melt_surft(:,:) + ! Snowmelt on land tiles (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: snowinc_surft(:,:) + ! Increment to snow from sublimation and melting on + ! surface tiles (kg m-2 TS-1) + REAL(KIND=real_jlslsm), ALLOCATABLE :: tot_tfall_surft(:,:) + ! Canopy througfall rate on land tiles (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: sice_melt(:,:,:) + ! Sea ice top melt on sea ice categories (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: ei_sice(:,:,:) + ! Sea ice sublimation on sea ice categories (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: surf_ht_flux_sice(:,:,:) + ! Heat flux through sea-ice (W/m2, positive downwards) + ! used in models coupled to ocean/sea-ice + REAL(KIND=real_jlslsm), ALLOCATABLE :: sea_ice_htf_sicat(:,:,:) + ! Heat flux through sea-ice (W/m2, positive downwards) + ! used in models NOT coupled to ocean/sea-ice + REAL(KIND=real_jlslsm), ALLOCATABLE :: snomlt_sub_htf_gb(:) + ! Sub-canopy snowmelt heat flux (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: snow_melt_gb(:) + ! snowmelt on land points (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: sub_surf_roff_gb(:) + ! Sub-surface runoff (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: surf_ht_flux_ij(:,:) + ! Net downward heat flux at surface over land and sea-ice fraction of + !gridbox (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: snow_soil_htf(:,:) + ! Heat flux under snow to subsurface on tiles (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: surf_htf_surft(:,:) + ! Surface heat flux on land tiles (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: surf_roff_gb(:) + ! Surface runoff (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: radnet_surft(:,:) + ! Surface net radiation on tiles ( W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: tot_tfall_gb(:) + ! Total throughfall (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: tstar_ij(:,:) + ! GBM surface temperature (K) + REAL(KIND=real_jlslsm), ALLOCATABLE :: emis_surft(:,:) + ! Tile emissivity + REAL(KIND=real_jlslsm), ALLOCATABLE :: sw_surft(:,:) + ! Surface net shortwave on tiles (W/m2) + REAL(KIND=real_jlslsm), ALLOCATABLE :: rflow_gb(:) + ! River outflow on model grid (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: rrun_gb(:) + ! Runoff after river routing on model grid (kg/m2/s) + REAL(KIND=real_jlslsm), ALLOCATABLE :: z0m_surft(:,:) + ! Tile roughness lengths for momentum. + REAL(KIND=real_jlslsm), ALLOCATABLE :: z0h_surft(:,:) + ! Tile roughness lengths for heat and moisture (m). + REAL(KIND=real_jlslsm), ALLOCATABLE :: lake_evap(:) + ! Evaporation from lakes (kg m-2 s-1). +END TYPE + +TYPE :: fluxes_type + + ! T if varying grey emissivity is used for surface tile type + LOGICAL, POINTER :: l_emis_surft_set(:) + + ! anthrop_heat is required by both the UM and standalone configurations + REAL(KIND=real_jlslsm), POINTER :: anthrop_heat_surft(:,:) + ! Additional heat source on surface tiles used for anthropgenic urban heat + ! source (W/m2) + REAL(KIND=real_jlslsm), POINTER :: surf_ht_store_surft(:,:) + ! Diagnostic to store values of C*(dT/dt) during calculation of energy + ! balance + REAL(KIND=real_jlslsm), POINTER :: sw_sicat(:,:) + ! Net SW on sea ice categories + REAL(KIND=real_jlslsm), POINTER :: alb_sicat(:,:,:) + ! Albedo of sea ice categories (ij point, sicat, band- see below) + REAL(KIND=real_jlslsm), POINTER :: penabs_rad_frac(:,:,:) + ! Fraction of downward solar that penetrates the sea ice and is absorbed + REAL(KIND=real_jlslsm), POINTER :: sw_sea(:) + ! Net SW on open sea + REAL(KIND=real_jlslsm), POINTER :: alb_surft(:,:,:) + ! Albedo for surface tiles + ! (:,:,1) direct beam visible + ! (:,:,2) diffuse visible + ! (:,:,3) direct beam near-IR + ! (:,:,4) diffuse near-IR + REAL(KIND=real_jlslsm), POINTER :: e_sea_ij(:,:) + ! Evaporation from sea times leads fraction. Zero over land + ! (kg per square metre per sec) + REAL(KIND=real_jlslsm), POINTER :: ecan_ij(:,:) + ! Gridbox mean evaporation from canopy/surface store (kg/m2/s) + ! Zero over sea + REAL(KIND=real_jlslsm), POINTER :: ecan_surft(:,:) + ! Canopy evaporation from for snow-free land tiles + REAL(KIND=real_jlslsm), POINTER :: ei_ij(:,:) + ! Sublimation from lying snow or sea-ice (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: ei_surft(:,:) + ! EI for land tiles + REAL(KIND=real_jlslsm), POINTER :: esoil_ij_soilt(:,:,:) + ! Surface evapotranspiration from soil moisture store (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: esoil_surft(:,:) + ! ESOIL for snow-free land tiles + REAL(KIND=real_jlslsm), POINTER :: ext_soilt(:,:,:) + ! Extraction of water from each soil layer (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: fqw_surft(:,:) + ! Surface FQW for land tiles + REAL(KIND=real_jlslsm), POINTER :: fqw_sicat(:,:,:) + ! Surface FQW for sea-ice + REAL(KIND=real_jlslsm), POINTER :: fsmc_pft(:,:) + ! Moisture availability factor. + REAL(KIND=real_jlslsm), POINTER :: fsnow(:,:) + ! Snow cover fraction on tiles + REAL(KIND=real_jlslsm), POINTER :: ftl_sicat(:,:,:) + ! Surface FTL for sea-ice + REAL(KIND=real_jlslsm), POINTER :: ftl_surft(:,:) + ! Surface FTL for land tiles + REAL(KIND=real_jlslsm), POINTER :: h_sea_ij(:,:) + ! Surface sensible heat flux over sea times leads fraction (W/m2) + REAL(KIND=real_jlslsm), POINTER :: hf_snow_melt_gb(:) + ! Gridbox snowmelt heat flux (W/m2) + REAL(KIND=real_jlslsm), POINTER :: land_albedo_ij(:,:,:) + ! GBM albedo + ! (:,:,1) direct beam visible + ! (:,:,2) diffuse visible + ! (:,:,3) direct beam near-IR + ! (:,:,4) diffuse near-IR + REAL(KIND=real_jlslsm), POINTER :: le_surft(:,:) + ! Surface latent heat flux for land tiles + REAL(KIND=real_jlslsm), POINTER :: melt_surft(:,:) + ! Snowmelt on land tiles (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: snowinc_surft(:,:) + ! Increment to snow from sublimation and melting on + ! surface tiles (kg m-2 TS-1) + REAL(KIND=real_jlslsm), POINTER :: tot_tfall_surft(:,:) + ! Canopy througfall rate on land tiles (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: sice_melt(:,:,:) + ! Sea ice top melt on sea ice categories (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: ei_sice(:,:,:) + ! Sea ice sublimation on sea ice categories (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: surf_ht_flux_sice(:,:,:) + ! Heat flux through sea-ice (W/m2, positive downwards) + ! used in models coupled to ocean/sea-ice + REAL(KIND=real_jlslsm), POINTER :: sea_ice_htf_sicat(:,:,:) + ! Heat flux through sea-ice (W/m2, positive downwards) + ! used in models NOT coupled to ocean/sea-ice + REAL(KIND=real_jlslsm), POINTER :: snomlt_sub_htf_gb(:) + ! Sub-canopy snowmelt heat flux (W/m2) + REAL(KIND=real_jlslsm), POINTER :: snow_melt_gb(:) + ! snowmelt on land points (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: sub_surf_roff_gb(:) + ! Sub-surface runoff (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: surf_ht_flux_ij(:,:) + ! Net downward heat flux at surface over land and sea-ice fraction of + !gridbox (W/m2) + REAL(KIND=real_jlslsm), POINTER :: snow_soil_htf(:,:) + ! Heat flux under snow to subsurface on tiles (W/m2) + REAL(KIND=real_jlslsm), POINTER :: surf_htf_surft(:,:) + ! Surface heat flux on land tiles (W/m2) + REAL(KIND=real_jlslsm), POINTER :: surf_roff_gb(:) + ! Surface runoff (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: radnet_surft(:,:) + ! Surface net radiation on tiles ( W/m2) + REAL(KIND=real_jlslsm), POINTER :: tot_tfall_gb(:) + ! Total throughfall (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: tstar_ij(:,:) + ! GBM surface temperature (K) + REAL(KIND=real_jlslsm), POINTER :: emis_surft(:,:) + ! Tile emissivity + REAL(KIND=real_jlslsm), POINTER :: sw_surft(:,:) + ! Surface net shortwave on tiles (W/m2) + REAL(KIND=real_jlslsm), POINTER :: rflow_gb(:) + ! River outflow on model grid (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: rrun_gb(:) + ! Runoff after river routing on model grid (kg/m2/s) + REAL(KIND=real_jlslsm), POINTER :: z0m_surft(:,:) + ! Tile roughness lengths for momentum. + REAL(KIND=real_jlslsm), POINTER :: z0h_surft(:,:) + ! Tile roughness lengths for heat and moisture (m) + REAL(KIND=real_jlslsm), POINTER :: lake_evap(:) + ! Evaporation from lakes (kg m-2 s-1). +END TYPE + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='FLUXES_MOD' + +CONTAINS + +SUBROUTINE fluxes_alloc(land_pts, t_i_length, t_j_length, & + nsurft, npft, nsoilt, sm_levels, & + nice, nice_use, & + fluxes_data) + +!No USE statements other than error reporting and Dr Hook +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +IMPLICIT NONE + +!Arguments +INTEGER, INTENT(IN) :: land_pts, t_i_length, t_j_length, & + nsurft, npft, nsoilt, sm_levels, & + nice, nice_use + +TYPE(fluxes_data_type), INTENT(IN OUT) :: fluxes_data + +!Local variables + +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='FLUXES_ALLOC' + +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +ALLOCATE(fluxes_data%l_emis_surft_set(nsurft)) + +fluxes_data%l_emis_surft_set(:) = .FALSE. + +ALLOCATE(fluxes_data%surf_ht_store_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%anthrop_heat_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%sw_sicat(t_i_length * t_j_length, nice_use)) +ALLOCATE(fluxes_data%alb_sicat(t_i_length * t_j_length, nice_use, 4)) +ALLOCATE(fluxes_data%penabs_rad_frac(t_i_length * t_j_length, nice_use, 4)) +ALLOCATE(fluxes_data%sw_sea(t_i_length * t_j_length)) + +fluxes_data%surf_ht_store_surft(:,:) = 0.0 +fluxes_data%anthrop_heat_surft(:,:) = 0.0 +fluxes_data%sw_sicat(:,:) = 0.0 +fluxes_data%alb_sicat(:,:,:) = 0.0 +fluxes_data%penabs_rad_frac(:,:,:) = 0.0 +fluxes_data%sw_sea(:) = 0.0 + +ALLOCATE(fluxes_data%sub_surf_roff_gb(land_pts)) +ALLOCATE(fluxes_data%surf_roff_gb(land_pts)) +ALLOCATE(fluxes_data%alb_surft(land_pts,nsurft,4) ) +ALLOCATE(fluxes_data%fsnow(land_pts,nsurft)) +ALLOCATE(fluxes_data%tstar_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%e_sea_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%fsmc_pft(land_pts,npft)) +ALLOCATE(fluxes_data%ftl_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%le_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%h_sea_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%fqw_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%fqw_sicat(t_i_length,t_j_length,nice_use)) +ALLOCATE(fluxes_data%ftl_sicat(t_i_length,t_j_length,nice_use)) +ALLOCATE(fluxes_data%ecan_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%esoil_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%surf_ht_flux_sice(t_i_length,t_j_length,nice_use)) +ALLOCATE(fluxes_data%sea_ice_htf_sicat(t_i_length,t_j_length,nice)) +ALLOCATE(fluxes_data%surf_ht_flux_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%surf_htf_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%snow_soil_htf(land_pts,nsurft)) +ALLOCATE(fluxes_data%land_albedo_ij(t_i_length,t_j_length,4) ) +ALLOCATE(fluxes_data%ei_ij(t_i_length,t_j_length)) +ALLOCATE(fluxes_data%ei_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%ecan_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%esoil_ij_soilt(t_i_length,t_j_length,nsoilt)) +ALLOCATE(fluxes_data%ext_soilt(land_pts,nsoilt,sm_levels)) +ALLOCATE(fluxes_data%hf_snow_melt_gb(land_pts)) +ALLOCATE(fluxes_data%radnet_surft(land_pts,nsurft) ) +ALLOCATE(fluxes_data%sw_surft(land_pts,nsurft) ) +ALLOCATE(fluxes_data%emis_surft(land_pts,nsurft) ) +ALLOCATE(fluxes_data%snow_melt_gb(land_pts)) +ALLOCATE(fluxes_data%snomlt_sub_htf_gb(land_pts)) +ALLOCATE(fluxes_data%tot_tfall_gb(land_pts)) +ALLOCATE(fluxes_data%melt_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%snowinc_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%tot_tfall_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%sice_melt(t_i_length,t_j_length,nice_use)) +ALLOCATE(fluxes_data%ei_sice(t_i_length,t_j_length,nice_use)) +ALLOCATE(fluxes_data%z0m_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%z0h_surft(land_pts,nsurft)) +ALLOCATE(fluxes_data%rflow_gb(land_pts)) +ALLOCATE(fluxes_data%rrun_gb(land_pts)) +ALLOCATE(fluxes_data%lake_evap(land_pts)) + +fluxes_data%sub_surf_roff_gb(:) = 0.0 +fluxes_data%surf_roff_gb(:) = 0.0 +fluxes_data%alb_surft(:,:,:) = 0.0 +fluxes_data%fsnow(:,:) = 0.0 +fluxes_data%tstar_ij(:,:) = 0.0 +fluxes_data%e_sea_ij(:,:) = 0.0 +fluxes_data%fsmc_pft(:,:) = 0.0 +fluxes_data%ftl_surft(:,:) = 0.0 +fluxes_data%le_surft(:,:) = 0.0 +fluxes_data%h_sea_ij(:,:) = 0.0 +fluxes_data%fqw_surft(:,:) = 0.0 +fluxes_data%fqw_sicat(:,:,:) = 0.0 +fluxes_data%ftl_sicat(:,:,:) = 0.0 +fluxes_data%ecan_ij(:,:) = 0.0 +fluxes_data%esoil_surft(:,:) = 273.15 +fluxes_data%surf_ht_flux_sice(:,:,:) = 0.0 +fluxes_data%sea_ice_htf_sicat(:,:,:) = 0.0 +fluxes_data%surf_ht_flux_ij(:,:) = 0.0 +fluxes_data%surf_htf_surft(:,:) = 0.0 +fluxes_data%snow_soil_htf(:,:) = 0.0 +fluxes_data%land_albedo_ij(:,:,:) = 0.0 +fluxes_data%ei_ij(:,:) = 0.0 +fluxes_data%ei_surft(:,:) = 0.0 +fluxes_data%ecan_surft(:,:) = 0.0 +fluxes_data%esoil_ij_soilt(:,:,:) = 0.0 +fluxes_data%ext_soilt(:,:,:) = 0.0 +fluxes_data%hf_snow_melt_gb(:) = 0.0 +fluxes_data%radnet_surft(:,:) = 0.0 +fluxes_data%sw_surft(:,:) = 0.0 +fluxes_data%emis_surft(:,:) = 0.0 +fluxes_data%snow_melt_gb(:) = 0.0 +fluxes_data%snomlt_sub_htf_gb(:) = 0.0 +fluxes_data%tot_tfall_gb(:) = 0.0 +fluxes_data%melt_surft(:,:) = 0.0 +fluxes_data%snowinc_surft(:,:) = 0.0 +fluxes_data%tot_tfall_surft(:,:) = 0.0 +fluxes_data%sice_melt(:,:,:) = 0.0 +fluxes_data%ei_sice(:,:,:) = 0.0 +fluxes_data%z0m_surft(:,:) = 0.0 +fluxes_data%z0h_surft(:,:) = 0.0 +fluxes_data%rflow_gb(:) = 0.0 +fluxes_data%rrun_gb(:) = 0.0 +fluxes_data%lake_evap(:) = rmdi + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE fluxes_alloc + +!=============================================================================== +SUBROUTINE fluxes_dealloc(fluxes_data) + +!No USE statements other than Dr Hook +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +IMPLICIT NONE + +!Arguments +TYPE(fluxes_data_type), INTENT(IN OUT) :: fluxes_data + +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='FLUXES_DEALLOC' + +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +DEALLOCATE(fluxes_data%l_emis_surft_set) +DEALLOCATE(fluxes_data%surf_ht_store_surft) +DEALLOCATE(fluxes_data%anthrop_heat_surft) +DEALLOCATE(fluxes_data%sw_sicat) +DEALLOCATE(fluxes_data%alb_sicat) +DEALLOCATE(fluxes_data%penabs_rad_frac) +DEALLOCATE(fluxes_data%sw_sea) +DEALLOCATE(fluxes_data%sub_surf_roff_gb) +DEALLOCATE(fluxes_data%surf_roff_gb) +DEALLOCATE(fluxes_data%alb_surft) +DEALLOCATE(fluxes_data%fsnow) +DEALLOCATE(fluxes_data%tstar_ij) +DEALLOCATE(fluxes_data%e_sea_ij) +DEALLOCATE(fluxes_data%fsmc_pft) +DEALLOCATE(fluxes_data%ftl_surft) +DEALLOCATE(fluxes_data%le_surft) +DEALLOCATE(fluxes_data%h_sea_ij) +DEALLOCATE(fluxes_data%fqw_surft) +DEALLOCATE(fluxes_data%fqw_sicat) +DEALLOCATE(fluxes_data%ftl_sicat) +DEALLOCATE(fluxes_data%ecan_ij) +DEALLOCATE(fluxes_data%esoil_surft) +DEALLOCATE(fluxes_data%surf_ht_flux_sice) +DEALLOCATE(fluxes_data%sea_ice_htf_sicat) +DEALLOCATE(fluxes_data%surf_ht_flux_ij) +DEALLOCATE(fluxes_data%surf_htf_surft) +DEALLOCATE(fluxes_data%snow_soil_htf) +DEALLOCATE(fluxes_data%land_albedo_ij) +DEALLOCATE(fluxes_data%ei_ij) +DEALLOCATE(fluxes_data%ei_surft) +DEALLOCATE(fluxes_data%ecan_surft) +DEALLOCATE(fluxes_data%esoil_ij_soilt) +DEALLOCATE(fluxes_data%ext_soilt) +DEALLOCATE(fluxes_data%hf_snow_melt_gb) +DEALLOCATE(fluxes_data%radnet_surft) +DEALLOCATE(fluxes_data%sw_surft) +DEALLOCATE(fluxes_data%emis_surft) +DEALLOCATE(fluxes_data%snow_melt_gb) +DEALLOCATE(fluxes_data%snomlt_sub_htf_gb) +DEALLOCATE(fluxes_data%tot_tfall_gb) +DEALLOCATE(fluxes_data%melt_surft) +DEALLOCATE(fluxes_data%snowinc_surft) +DEALLOCATE(fluxes_data%tot_tfall_surft) +DEALLOCATE(fluxes_data%sice_melt) +DEALLOCATE(fluxes_data%ei_sice) +DEALLOCATE(fluxes_data%z0m_surft) +DEALLOCATE(fluxes_data%z0h_surft) +DEALLOCATE(fluxes_data%rflow_gb) +DEALLOCATE(fluxes_data%rrun_gb) +DEALLOCATE(fluxes_data%lake_evap) + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE fluxes_dealloc + +!============================================================================== +SUBROUTINE fluxes_assoc(fluxes,fluxes_data) + +!No USE statements other than Dr Hook +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +IMPLICIT NONE + +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes + !Instance of the pointer type we are associating + +TYPE(fluxes_data_type), INTENT(IN OUT), TARGET :: fluxes_data + !Instance of the data type we are associating to + +!Local variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='FLUXES_ASSOC' + +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +CALL fluxes_nullify(fluxes) + +fluxes%l_emis_surft_set => fluxes_data%l_emis_surft_set +fluxes%surf_ht_store_surft => fluxes_data%surf_ht_store_surft +fluxes%anthrop_heat_surft => fluxes_data%anthrop_heat_surft +fluxes%sw_sicat => fluxes_data%sw_sicat +fluxes%alb_sicat => fluxes_data%alb_sicat +fluxes%penabs_rad_frac => fluxes_data%penabs_rad_frac +fluxes%sw_sea => fluxes_data%sw_sea +fluxes%sub_surf_roff_gb => fluxes_data%sub_surf_roff_gb +fluxes%surf_roff_gb => fluxes_data%surf_roff_gb +fluxes%alb_surft => fluxes_data%alb_surft +fluxes%fsnow => fluxes_data%fsnow +fluxes%tstar_ij => fluxes_data%tstar_ij +fluxes%e_sea_ij => fluxes_data%e_sea_ij +fluxes%fsmc_pft => fluxes_data%fsmc_pft +fluxes%ftl_surft => fluxes_data%ftl_surft +fluxes%le_surft => fluxes_data%le_surft +fluxes%h_sea_ij => fluxes_data%h_sea_ij +fluxes%fqw_surft => fluxes_data%fqw_surft +fluxes%fqw_sicat => fluxes_data%fqw_sicat +fluxes%ftl_sicat => fluxes_data%ftl_sicat +fluxes%ecan_ij => fluxes_data%ecan_ij +fluxes%esoil_surft => fluxes_data%esoil_surft +fluxes%surf_ht_flux_sice => fluxes_data%surf_ht_flux_sice +fluxes%sea_ice_htf_sicat => fluxes_data%sea_ice_htf_sicat +fluxes%surf_ht_flux_ij => fluxes_data%surf_ht_flux_ij +fluxes%surf_htf_surft => fluxes_data%surf_htf_surft +fluxes%snow_soil_htf => fluxes_data%snow_soil_htf +fluxes%land_albedo_ij => fluxes_data%land_albedo_ij +fluxes%ei_ij => fluxes_data%ei_ij +fluxes%ei_surft => fluxes_data%ei_surft +fluxes%ecan_surft => fluxes_data%ecan_surft +fluxes%esoil_ij_soilt => fluxes_data%esoil_ij_soilt +fluxes%ext_soilt => fluxes_data%ext_soilt +fluxes%hf_snow_melt_gb => fluxes_data%hf_snow_melt_gb +fluxes%radnet_surft => fluxes_data%radnet_surft +fluxes%sw_surft => fluxes_data%sw_surft +fluxes%emis_surft => fluxes_data%emis_surft +fluxes%snow_melt_gb => fluxes_data%snow_melt_gb +fluxes%snomlt_sub_htf_gb => fluxes_data%snomlt_sub_htf_gb +fluxes%tot_tfall_gb => fluxes_data%tot_tfall_gb +fluxes%melt_surft => fluxes_data%melt_surft +fluxes%snowinc_surft => fluxes_data%snowinc_surft +fluxes%tot_tfall_surft => fluxes_data%tot_tfall_surft +fluxes%sice_melt => fluxes_data%sice_melt +fluxes%ei_sice => fluxes_data%ei_sice +fluxes%z0m_surft => fluxes_data%z0m_surft +fluxes%z0h_surft => fluxes_data%z0h_surft +fluxes%rflow_gb => fluxes_data%rflow_gb +fluxes%rrun_gb => fluxes_data%rrun_gb +fluxes%lake_evap => fluxes_data%lake_evap + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE fluxes_assoc + +!============================================================================== + +SUBROUTINE fluxes_nullify(fluxes) + +!No USE statements other than Dr Hook +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +IMPLICIT NONE + +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes + !Instance of the pointer type we are nullifying + +!Local variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='FLUXES_NULLIFY' + +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +NULLIFY(fluxes%l_emis_surft_set) +NULLIFY(fluxes%surf_ht_store_surft) +NULLIFY(fluxes%anthrop_heat_surft) +NULLIFY(fluxes%sw_sicat) +NULLIFY(fluxes%alb_sicat) +NULLIFY(fluxes%penabs_rad_frac) +NULLIFY(fluxes%sw_sea) +NULLIFY(fluxes%sub_surf_roff_gb) +NULLIFY(fluxes%surf_roff_gb) +NULLIFY(fluxes%alb_surft) +NULLIFY(fluxes%fsnow) +NULLIFY(fluxes%tstar_ij) +NULLIFY(fluxes%e_sea_ij) +NULLIFY(fluxes%fsmc_pft) +NULLIFY(fluxes%ftl_surft) +NULLIFY(fluxes%le_surft) +NULLIFY(fluxes%h_sea_ij) +NULLIFY(fluxes%fqw_surft) +NULLIFY(fluxes%fqw_sicat) +NULLIFY(fluxes%ftl_sicat) +NULLIFY(fluxes%ecan_ij) +NULLIFY(fluxes%esoil_surft) +NULLIFY(fluxes%surf_ht_flux_sice) +NULLIFY(fluxes%sea_ice_htf_sicat) +NULLIFY(fluxes%surf_ht_flux_ij) +NULLIFY(fluxes%surf_htf_surft) +NULLIFY(fluxes%snow_soil_htf) +NULLIFY(fluxes%land_albedo_ij) +NULLIFY(fluxes%ei_ij) +NULLIFY(fluxes%ei_surft) +NULLIFY(fluxes%ecan_surft) +NULLIFY(fluxes%esoil_ij_soilt) +NULLIFY(fluxes%ext_soilt) +NULLIFY(fluxes%hf_snow_melt_gb) +NULLIFY(fluxes%radnet_surft) +NULLIFY(fluxes%sw_surft) +NULLIFY(fluxes%emis_surft) +NULLIFY(fluxes%snow_melt_gb) +NULLIFY(fluxes%snomlt_sub_htf_gb) +NULLIFY(fluxes%tot_tfall_gb) +NULLIFY(fluxes%melt_surft) +NULLIFY(fluxes%snowinc_surft) +NULLIFY(fluxes%tot_tfall_surft) +NULLIFY(fluxes%sice_melt) +NULLIFY(fluxes%ei_sice) +NULLIFY(fluxes%z0m_surft) +NULLIFY(fluxes%z0h_surft) +NULLIFY(fluxes%rflow_gb) +NULLIFY(fluxes%rrun_gb) +NULLIFY(fluxes%lake_evap) + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE fluxes_nullify + +END MODULE fluxes_mod diff --git a/src/io/model_interface/jules_radiation_mod.F90 b/src/io/model_interface/jules_radiation_mod.F90 new file mode 100644 index 00000000..2f7b9fb6 --- /dev/null +++ b/src/io/model_interface/jules_radiation_mod.F90 @@ -0,0 +1,393 @@ +! *****************************COPYRIGHT******************************* +! (C) Crown copyright Met Office. All rights reserved. +! For further details please refer to the file COPYRIGHT.txt +! which you should have received as part of this distribution. +! *****************************COPYRIGHT******************************* + +MODULE jules_radiation_mod + +!----------------------------------------------------------------------------- +! Description: +! Contains radiation options and a namelist for setting them +! +! Code Owner: Please refer to ModuleLeaders.txt +! This file belongs in TECHNICAL +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +USE missing_data_mod, ONLY: imdi, rmdi + +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +!----------------------------------------------------------------------------- +! Switches +!----------------------------------------------------------------------------- +LOGICAL :: & + l_spec_albedo = .FALSE., & + ! Switch spectrally varying land albedo + l_spec_alb_bs = .FALSE., & + ! Switch to have only a bluesky albedo when using + ! spectrally varying albedo + l_niso_direct = .FALSE., & + ! Switch to use the true non-isotropic form of direct + ! back-scattering for the canopy albedo + l_snow_albedo = .FALSE., & + ! Switch for prognostic snow albedo (on land) + l_embedded_snow = .FALSE., & + ! Switch for calculation of albedo with snow in + ! the canopy (exclusive of l_snow_albedo) + l_mask_snow_orog = .FALSE., & + ! Switch for orographic masking of snow + l_albedo_obs = .FALSE., & + ! scale the albedo on tiles to agree with obs + l_dolr_land_black = .FALSE., & + ! Do not use the surface emissivity in adjusting the OLR + ! at land points. + ! This flag is introduced for historical compatibility + ! only. There is no equivalent choice at sea points. + l_spec_sea_alb = .FALSE., & + ! Switch spectrally varying open sea albedo + l_sea_alb_var_chl = .FALSE., & + ! Switch varying chlorophyll in open sea albedo + l_partition_albsoil = .FALSE., & + ! Impose a spectral partition on the bare soil albedo + l_hapke_soil = .FALSE. + ! Calculate the direct albedo using Hapke's two-stream + ! approximation of isotropic scattering by soil + ! particles, but omitting the opposition effect. + ! DOI=10.1029/JB086iB04p03039 + +LOGICAL :: & + l_cosz = .TRUE. + ! Switch for turning on calculations of cosz + ! Used in standalone JULES only. This is effectively .true. when & + ! coupled to the UM as it is calculated therein. + +INTEGER :: & + i_sea_alb_method = imdi + ! Method of diagnosing the Ocean Surface Albedo + ! 1 - Briegleb and Ramanathan, 1982, J. Appl. Met. + ! (doi:10.1175/1520-0450(1982)021<1160:SADVIC>2.0.CO;2) + ! 2 - Modified Barker and Li, 1995, J. Climate, + ! (doi:10.1175/1520-0442(1995)008<2213:ISOCSS>2.0.CO;2) + ! 3 - Jin et al. 2011, Optics Express + ! (doi:10.1364/OE.19.026429) + +REAL(KIND=real_jlslsm) :: & + wght_alb(4) = [ 0.0, 0.5, 0.0, 0.5 ] + ! Weights to form broad-band albedo from components +REAL :: & + ratio_albsoil = rmdi, & + ! Ratio of NIR soil albedo to VIS soil albedo + swdn_frac_albsoil = rmdi + ! Fraction of downward SW in NIR used to partition soil albedo + +REAL(KIND=real_jlslsm) :: & + fixed_sea_albedo = rmdi + ! Value for open sea albedo if using simple fixed value method + + +!----------------------------------------------------------------------------- +! Single namelist definition for UM and standalone +!----------------------------------------------------------------------------- +NAMELIST / jules_radiation/ & + l_spec_albedo, l_spec_alb_bs, l_albedo_obs, l_niso_direct, & + l_snow_albedo, l_embedded_snow, l_mask_snow_orog, & + l_dolr_land_black, l_spec_sea_alb, l_sea_alb_var_chl, i_sea_alb_method, & + wght_alb, fixed_sea_albedo, & + l_partition_albsoil, ratio_albsoil, swdn_frac_albsoil, l_hapke_soil, & +! Standalone only switches + l_cosz + + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='JULES_RADIATION_MOD' + +CONTAINS + +#if !defined(RIVERS_ONLY) +SUBROUTINE check_jules_radiation() + +USE ancil_info, ONLY: rad_nband + +USE ereport_mod, ONLY: ereport +USE jules_print_mgr, ONLY: jules_message +USE jules_snow_mod, ONLY: i_snow_tile + +!----------------------------------------------------------------------------- +! Description: +! Checks JULES_RADIATION namelist for consistency +! +! Code Owner: Please refer to ModuleLeaders.txt +! This file belongs in TECHNICAL +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +IMPLICIT NONE + +INTEGER :: errorstatus + +! The snow scheme only calculates spectrally varying albedos +IF ( l_snow_albedo .AND. .NOT. l_spec_albedo ) THEN + errorstatus = 101 + CALL ereport("check_jules_radiation", errorstatus, & + "If l_snow_albedo = T then l_spec_albedo must also be T") +END IF + +! The embedded snow scheme only calculates spectrally varying albedos +IF ( l_embedded_snow .AND. .NOT. l_spec_albedo ) THEN + errorstatus = 1001 + CALL ereport("check_jules_radiation", errorstatus, & + "If l_embedded_snow = T then l_spec_albedo must also be T") +END IF + +! The embedded snow and the original snow albedo scheme are mutually +! exclusive +IF ( l_embedded_snow .AND. l_snow_albedo ) THEN + errorstatus = 1002 + WRITE(jules_message,'(A,2(1x,L1))') & + "l_embedded_snow and l_snow_albedo are mutually exclusive: " // & + "l_embedded_snow, l_snow_albedo = ", l_embedded_snow, l_snow_albedo + CALL ereport("check_jules_radiation", errorstatus, jules_message) +END IF + +! Required options if a separate snow tile is being used +IF ( ANY(i_snow_tile == 1) ) THEN + IF ( (.NOT. l_spec_albedo) .OR. l_embedded_snow .OR. & + (.NOT. l_snow_albedo) ) THEN + errorstatus = 1003 + CALL ereport("check_jules_radiation", errorstatus, & + "l_spec_albedo = T, l_embedded_snow = F and l_snow_albedo = T" // & + "required if ANY i_snow_tile = 1") + END IF +END IF + +!Can set the size of rad_nband +IF (l_spec_albedo) THEN + rad_nband = 2 +ELSE + rad_nband = 1 +END IF + +! Check that the spectral partitioning of the soil albedo is sensible. +IF (l_partition_albsoil) THEN + IF ( (swdn_frac_albsoil < 0.0) .OR. & + (swdn_frac_albsoil > 1.0) ) THEN + errorstatus = 2001 + CALL ereport("check_jules_radiation", errorstatus, & + "swdn_frac_albsoil cannot lie outside the range [0, 1].") + END IF + IF ( (ratio_albsoil < 1.0) .OR. & + (ratio_albsoil > 10.0) ) THEN + errorstatus = 2002 + CALL ereport("check_jules_radiation", errorstatus, & + "ratio_albsoil is not permoitted to lie outside " // & + "the range [1, 10].") + END IF +END IF + +END SUBROUTINE check_jules_radiation +#endif + +SUBROUTINE print_nlist_jules_radiation() + +USE jules_print_mgr, ONLY: jules_print + +IMPLICIT NONE + +CHARACTER(LEN=50000) :: lineBuffer + +CALL jules_print('jules_radiation', & + 'Contents of namelist jules_radiation') + +WRITE(lineBuffer, *) ' l_spec_albedo = ', l_spec_albedo +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_spec_alb_bs = ', l_spec_alb_bs +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_niso_direct = ', l_niso_direct +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_snow_albedo = ', l_snow_albedo +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_embedded_snow = ', l_embedded_snow +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_mask_snow_orog = ', l_mask_snow_orog +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_albedo_obs = ', l_albedo_obs +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_dolr_land_black = ', l_dolr_land_black +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_spec_sea_alb = ', l_spec_sea_alb +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_sea_alb_var_chl = ', l_sea_alb_var_chl +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' i_sea_alb_method = ', i_sea_alb_method +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' fixed_sea_albedo = ', fixed_sea_albedo +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' wght_alb = ', wght_alb +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_partition_albsoil = ', l_partition_albsoil +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' ratio_albsoil = ', ratio_albsoil +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' swdn_frac_albsoil = ', swdn_frac_albsoil +CALL jules_print('jules_radiation', lineBuffer) + +WRITE(lineBuffer, *) ' l_hapke_soil = ', l_hapke_soil +CALL jules_print('jules_radiation', lineBuffer) + +CALL jules_print('jules_radiation', & + '- - - - - - end of namelist - - - - - -') + +END SUBROUTINE print_nlist_jules_radiation + +#if defined(UM_JULES) && !defined(LFRIC) +SUBROUTINE read_nml_jules_radiation (unitnumber) + +! Description: +! Read the JULES_RADIATION namelist + +USE setup_namelist, ONLY: setup_nml_type +USE check_iostat_mod, ONLY: check_iostat +USE UM_parcore, ONLY: mype + +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +USE errormessagelength_mod, ONLY: errormessagelength + +IMPLICIT NONE + +! Subroutine arguments +INTEGER, INTENT(IN) :: unitnumber + +INTEGER :: my_comm +INTEGER :: mpl_nml_type +INTEGER :: ErrorStatus +INTEGER :: icode +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='READ_NML_JULES_RADIATION' + +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 + +CHARACTER(LEN=errormessagelength) :: iomessage + +! set number of each type of variable in my_namelist type +INTEGER, PARAMETER :: no_of_types = 3 +INTEGER, PARAMETER :: n_int = 1 +INTEGER, PARAMETER :: n_log = 13 +INTEGER, PARAMETER :: n_real = 7 + +TYPE :: my_namelist + SEQUENCE + INTEGER :: i_sea_alb_method + LOGICAL :: l_spec_albedo + LOGICAL :: l_spec_alb_bs + LOGICAL :: l_niso_direct + LOGICAL :: l_snow_albedo + LOGICAL :: l_embedded_snow + LOGICAL :: l_mask_snow_orog + LOGICAL :: l_albedo_obs + LOGICAL :: l_dolr_land_black + LOGICAL :: l_spec_sea_alb + LOGICAL :: l_sea_alb_var_chl + LOGICAL :: l_partition_albsoil + LOGICAL :: l_hapke_soil + LOGICAL :: l_cosz + REAL(KIND=real_jlslsm) :: wght_alb(4) + REAL(KIND=real_jlslsm) :: fixed_sea_albedo + REAL :: ratio_albsoil + REAL :: swdn_frac_albsoil +END TYPE my_namelist + +TYPE (my_namelist) :: my_nml + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +CALL gc_get_communicator(my_comm, icode) + +CALL setup_nml_type(no_of_types, mpl_nml_type, n_int_in = n_int, & + n_log_in = n_log, n_real_in = n_real) + +IF (mype == 0) THEN + + READ (UNIT = unitnumber, NML = jules_radiation, IOSTAT = errorstatus, & + IOMSG = iomessage) + CALL check_iostat(errorstatus, "namelist jules_radiation", iomessage) + + my_nml % i_sea_alb_method = i_sea_alb_method + my_nml % l_spec_albedo = l_spec_albedo + my_nml % l_spec_alb_bs = l_spec_alb_bs + my_nml % l_niso_direct = l_niso_direct + my_nml % l_snow_albedo = l_snow_albedo + my_nml % l_embedded_snow = l_embedded_snow + my_nml % l_mask_snow_orog = l_mask_snow_orog + my_nml % l_albedo_obs = l_albedo_obs + my_nml % l_dolr_land_black = l_dolr_land_black + my_nml % l_spec_sea_alb = l_spec_sea_alb + my_nml % l_sea_alb_var_chl = l_sea_alb_var_chl + my_nml % l_cosz = l_cosz + my_nml % wght_alb = wght_alb + my_nml % fixed_sea_albedo = fixed_sea_albedo + my_nml % l_partition_albsoil = l_partition_albsoil + my_nml % l_hapke_soil = l_hapke_soil + my_nml % ratio_albsoil = ratio_albsoil + my_nml % swdn_frac_albsoil = swdn_frac_albsoil +END IF + +CALL mpl_bcast(my_nml,1,mpl_nml_type,0,my_comm,icode) + +IF (mype /= 0) THEN + + i_sea_alb_method = my_nml % i_sea_alb_method + l_spec_albedo = my_nml % l_spec_albedo + l_spec_alb_bs = my_nml % l_spec_alb_bs + l_snow_albedo = my_nml % l_snow_albedo + l_niso_direct = my_nml % l_niso_direct + l_embedded_snow = my_nml % l_embedded_snow + l_masK_snow_orog = my_nml % l_masK_snow_orog + l_albedo_obs = my_nml % l_albedo_obs + l_dolr_land_black = my_nml % l_dolr_land_black + l_spec_sea_alb = my_nml % l_spec_sea_alb + l_sea_alb_var_chl = my_nml % l_sea_alb_var_chl + l_cosz = my_nml % l_cosz + wght_alb = my_nml % wght_alb + fixed_sea_albedo = my_nml % fixed_sea_albedo + l_partition_albsoil = my_nml % l_partition_albsoil + l_hapke_soil = my_nml % l_hapke_soil + ratio_albsoil = my_nml % ratio_albsoil + swdn_frac_albsoil = my_nml % swdn_frac_albsoil +END IF + +CALL mpl_type_free(mpl_nml_type,icode) + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE read_nml_jules_radiation +#endif + +END MODULE jules_radiation_mod diff --git a/src/io/model_interface/jules_snow_mod.F90 b/src/io/model_interface/jules_snow_mod.F90 new file mode 100644 index 00000000..583a6bfa --- /dev/null +++ b/src/io/model_interface/jules_snow_mod.F90 @@ -0,0 +1,749 @@ +! *****************************COPYRIGHT******************************* +! (C) Crown copyright Met Office. All rights reserved. +! For further details please refer to the file COPYRIGHT.txt +! which you should have received as part of this distribution. +! *****************************COPYRIGHT******************************* + +MODULE jules_snow_mod + +USE max_dimensions, ONLY: snow_layers_max, npft_max, nsurft_max + +!----------------------------------------------------------------------------- +! Description: +! Contains snow options and a namelist for setting them +! +! Code Owner: Please refer to ModuleLeaders.txt +! This file belongs in TECHNICAL +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +USE missing_data_mod, ONLY: rmdi +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +INTEGER :: & + nsmax = 0 ! Maximum number of snow layers + +!----------------------------------------------------------------------------- +! Switches +!----------------------------------------------------------------------------- +LOGICAL :: & + l_snowdep_surf = .FALSE., & + ! use equivalent canopy snow depth for surface + ! calculations on tiles with a snow canopy + l_rho_snow_corr = .TRUE., & + ! Switch for using a correction to the density of the + ! snow pack when nsnow=0 when relayering in the new snow + ! scheme + ! Has no effect for nsmax < 1 + l_et_metamorph = .FALSE., & + ! Switch to include equitemperature metamorphism + l_snow_infilt = .FALSE., & + ! Include infiltration of rain water into snow + ! (Only applies to the multilayer scheme) + l_snow_nocan_hc = .FALSE. + ! Flag to negelect canopy heat capacity in the presence + ! of snow if not using an explicit canopy model + +INTEGER :: & + frac_snow_subl_melt = 0 + ! Switch for use of snow-cover fraction in the calculation + ! of sublimation and melting + ! 0 = off + ! 1 = on +! Parametrizations of thermal conductivity + +INTEGER, PARAMETER :: & + ip_snow_cond_yen81 = 0, & + ! Parametrization of snow conductivity following + ! Yen (1981), "Review of the thermal properties of snow, + ! ice and sea ice", Technical Report 81-10 + ! Cold Regions Research and Engineering Laboratory, + ! Hanover, NH + ip_snow_cond_calonne11 = 1 + ! Parametrization of snow conductivity following + ! Calonne et al. (2011), GRL, 38, L23501 + +INTEGER :: & + i_snow_cond_parm = 0 + ! Parametrization scheme for snow conductivity + +INTEGER :: & + i_snow_tile(nsurft_max) + ! Tiles without (0) or with (1) separate energy balance + ! for snow. Must be 0 on ice tiles +DATA i_snow_tile / nsurft_max * 0 / + +!----------------------------------------------------------------------------- +! Parametrization of the rate of growth of snow grains +INTEGER, PARAMETER :: & + ip_grain_growth_marshall89 = 0, & + ! Original scheme following Marshall (1989) + ip_grain_growth_taillandier_et = 1 + ! Rate of growth of grains at temperatures below + ! freezing follows the equitemperature scheme of + ! Taillandier et al. (2007), J. Geophys. Res., 112, + ! F03003 + +INTEGER :: & + i_grain_growth_opt = 0 + ! Option for rate of growth of snow grains + +!----------------------------------------------------------------------------- +! Options for relayering the snow pack. +INTEGER, PARAMETER :: & + ip_relayer_linear = 0, & + ! Relayer based on thickness + ip_relayer_rgrain_inv = 1 + ! As the linear scheme, but relayering the inverse of + ! the grain size + +INTEGER :: & + i_relayer_opt = 0 + ! Option for relayering the snow pack. + +!----------------------------------------------------------------------------- +! Options for melting at the base of thin snow packs. +INTEGER, PARAMETER :: & + ip_basal_melting_off = 0, & + ! Basal melting not considered in the zero-layer scheme + ip_basal_melting_diag_0 = 1 + ! Basal melting is diagnosed for thin snow packs + ! at points where the zero-layer scheme is active. + ! The multilayer scheme is unaffected. + +INTEGER :: & + i_basal_melting_opt = 0 + ! Option for basal melting + +! ---------------------------------------------------------------------------- +! Parametrisations for treatment of graupel +INTEGER :: & + graupel_options = 0 + ! Switch for treatment of graupel in the snow scheme + +INTEGER, PARAMETER :: & + graupel_as_snow = 0, & + ! Include graupel in the surface snowfall + ignore_graupel = 1 + ! Ignore graupel in the surface snowfall +!----------------------------------------------------------------------------- +! Parameters for equitemperature metamorphism +! +REAL(KIND=real_jlslsm) :: & + a_snow_et = 0.0, & + b_snow_et = 0.0, & + c_snow_et = 0.0 + ! Parameters in rate equation for ET metamorphism +REAL(KIND=real_jlslsm) :: & + rho_snow_et_crit = 0.0 + ! Critical density + +!----------------------------------------------------------------------------- +! Radiation parameters for namelist +! +! Defaults taken from JULES examples +!----------------------------------------------------------------------------- +! Parameters for snow grain size. +REAL(KIND=real_jlslsm) :: & + r0 = 50.0, & + ! Grain size for fresh snow (microns) + rmax = 2000.0 + ! Maximum snow grain size (microns) + + +! Parameters for snow grain size growth (HCTN30.16) +! Values are for melting snow, cold fresh snow and cold aged snow respectively +REAL(KIND=real_jlslsm) :: & + snow_ggr(3) ! snow grain area growth rates (microns**2 s-1) +DATA snow_ggr / 0.6, 0.06, 0.23e6 / + + +! Parameters for prognostic, spectral snow albedo. +REAL(KIND=real_jlslsm) :: & + amax(2) ! Maximum albedo for fresh snow (values for VIS and NIR) +DATA amax / 0.98, 0.7 / + +REAL(KIND=real_jlslsm) :: & + aicemax(2) ! Maximum albedo for bare land ice (values for VIS and NIR) + ! Values hand tuned by Sarah Shannon to match a range of + ! glaciers in standalone JULES (pers comm.) +DATA aicemax / 0.78, 0.36 / + +!----------------------------------------------------------------------- +! Parameters for prognostic, spectral snow albedo using the embedded +! canopy scheme. +!----------------------------------------------------------------------- +REAL(KIND=real_jlslsm) :: dce = 2.4848e-7 + ! Effective dimension of black carbon +REAL(KIND=real_jlslsm) :: cnr_g(2, 2) = RESHAPE( & + [ 0.900396, 0.00143416, 0.947363, 0.00570307 ], & + [ 2, 2 ] ) + ! Parameters fitting the asymmetry +REAL(KIND=real_jlslsm) :: cnr_om(3, 2, 2) = RESHAPE( & + [ 0.000869870, 0.000194185, 1.09046e-05, & + 0.0, 0.2709248, 0.0, & + 0.103851, 0.0201497, 0.00101678, & + 0.0, 0.4306217, 0.0 ], & + [ 3, 2, 2 ] ) + ! Parameters fitting the single scattering albedo + +! Parameters for (diagnostic) all-band snow albedo. +REAL(KIND=real_jlslsm) :: & + dtland = 2.0, & + ! Degrees C below freezing point at which snow albedo equals + ! cold deep snow albedo. This is 2 in HCTN30.4. + ! Must not be zero! + kland_numerator = 0.3, & + ! KLAND is calculated as KLAND_NUMERATOR / DTLAND once the + ! namelist has been read + maskd = 50.0 + ! Used in the exponent of equation weighting snow and snow-free + ! albedo to get tile albedo + ! This is 0.2 in HCTN30.5, where it is used as maskd*snowMass, + ! assuming snow density=250kg/m3 + ! It is now used as maskd*snowDepth, so maskd=50 gives the same + ! relationship as HCTN30.5 + +!----------------------------------------------------------------------------- +! Calculated radiation parameters +!----------------------------------------------------------------------------- +REAL(KIND=real_jlslsm) :: kland + ! Used in snow-ageing effect on snow albedo + ! This is 0.3 in HCTN30.4, although note that the last term + ! in that eqn should be divided by dtland + ! KLAND is calculated as KLAND_NUMERATOR / DTLAND once the + ! namelist has been read +REAL(KIND=real_jlslsm) :: tcland + ! Temperature below which snow albedo equals cold deep snow + ! albedo (HCTN30.4) + ! This is initialised as TM - DTLAND once the namelist has + ! been read + +!----------------------------------------------------------------------------- +! Other snow parameters for namelist +!----------------------------------------------------------------------------- +REAL(KIND=real_jlslsm) :: & + rho_snow_const = 250.0, & + ! constant density of lying snow (kg per m**3) + ! This is used: + ! (a) as the snow density when nsmax=0 + ! (b) with nsmax>0 and l_snowdep_surf=.TRUE. this + ! is the density of snow on the canopy + ! (not on the ground) + rho_snow_fresh = 100.0, & + ! Density of fresh snow (kg per m**3) + ! Only used with nsmax>0. + snow_hcon = 0.265, & + ! Thermal conductivity of lying snow (Watts per m per K) + snow_hcap = 0.63e6, & + ! Thermal capacity of lying snow (J/K/m3) + snowliqcap = 0.05, & + ! Liquid water holding capacity of lying snow as a fraction + ! of snow mass + snowinterceptfact = 0.7, & + ! Constant in relationship between mass of intercepted snow + ! and snowfall rate + snowloadlai = 4.4, & + ! Ratio of maximum canopy snow load to leaf area index(kg m-2) + snowunloadfact = 0.4, & + ! Constant in relationship between canopy snow unloading and + ! canopy snow melt rate + rho_firn_albedo = 550.0, & + ! Threshold surface density where albedo parameterisation + ! switches to scaling with surface density not grain size + ! (kg per m**3) + rho_firn_pore_restrict = 450.0, & + ! Threshold snowpack density where holding capacity/ability + ! to percolate meltwater starts to be reduced + ! (kg per m**3) + rho_firn_pore_closure = 850.0 + ! Threshold snowpack density where holding capacity/ability + ! to percolate meltwater becomes 0 + ! (kg per m**3) + +!----------------------------------------------------------------------------- +! Allocatable versions of namelist arrays +!----------------------------------------------------------------------------- +REAL(KIND=real_jlslsm) :: dzsnow(snow_layers_max) = rmdi + ! Prescribed thickness of snow layers (m) + ! This is the thickness of each snow layer when it is not + ! the bottom layer (note that dzSnow(nsMax) is not used + ! because that is always the bottom layer) + +LOGICAL :: cansnowtile(nsurft_max) + ! Switch for canopy snow model on each tile + ! Must be false for non-PFT tiles +DATA cansnowtile / nsurft_max * .FALSE. / + +LOGICAL :: cansnowpft(npft_max) + ! This is what appears in the namelist, to ensure that true + ! values can only given for pfts +DATA cansnowpft / npft_max * .FALSE. / + +REAL(KIND=real_jlslsm) :: unload_rate_cnst(npft_max) +! ! Constant canopy unloading rate for the tile +REAL(KIND=real_jlslsm) :: unload_rate_u(npft_max) +! ! Wind-speed dependent canopy unloading rate +! ! for the tile + +DATA unload_rate_cnst / npft_max * 0.0 / +DATA unload_rate_u / npft_max * 0.0 / + +REAL(KIND=real_jlslsm) :: can_clump(npft_max) + ! Clumping factor for snow in the canopy +REAL(KIND=real_jlslsm) :: lai_alb_lim_sn(npft_max) + ! Lower limit on permitted LAI in albedo with snow +REAL(KIND=real_jlslsm) :: n_lai_exposed(npft_max) + ! Shape parameter for exposed canopy with + ! embedded snow + +DATA can_clump / npft_max * 0.0 / +DATA lai_alb_lim_sn / npft_max * 0.0 / +DATA n_lai_exposed / npft_max * 0.0 / + + +!----------------------------------------------------------------------------- +! Single namelist definition for UM and standalone +!----------------------------------------------------------------------------- +NAMELIST / jules_snow/ & +! Switches + nsmax, l_snowdep_surf, l_rho_snow_corr, frac_snow_subl_melt, & + graupel_options,i_snow_tile, & +! Equitemperature metamorphism + l_et_metamorph, a_snow_et, b_snow_et, c_snow_et, rho_snow_et_crit, & +! Thermal conductivity of snow + i_snow_cond_parm, & +! Radiation parameters + r0, rmax, snow_ggr, amax, aicemax, maskd, dtland, kland_numerator, & +! Other snow parameters + rho_snow_const, rho_snow_fresh, rho_firn_albedo, & + snow_hcon, snow_hcap, snowliqcap, & + snowinterceptfact, snowloadlai, snowunloadfact, dzsnow, cansnowpft, & + unload_rate_cnst, unload_rate_u, can_clump, lai_alb_lim_sn, & + n_lai_exposed, l_snow_infilt, l_snow_nocan_hc, & + i_grain_growth_opt, i_relayer_opt, i_basal_melting_opt + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='JULES_SNOW_MOD' + +CONTAINS + +#if !defined(RIVERS_ONLY) +SUBROUTINE check_jules_snow() + +USE ereport_mod, ONLY: ereport + +USE water_constants_mod, ONLY: tm + +USE jules_surface_types_mod, ONLY: ice, npft, ntype +USE jules_surface_mod, ONLY: l_aggregate +USE jules_vegetation_mod, ONLY: can_model + +!----------------------------------------------------------------------------- +! Description: +! Checks JULES_SNOW namelist for consistency +! +! Code Owner: Please refer to ModuleLeaders.txt +! This file belongs in TECHNICAL +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +IMPLICIT NONE + +INTEGER :: errorstatus +CHARACTER(LEN=*), PARAMETER :: RoutineName='CHECK_JULES_SNOW' + +!----------------------------------------------------------------------------- + + +! Check that we have a suitable number of snow layers +IF ( nsmax > snow_layers_max ) THEN + errorstatus = 101 + CALL ereport(RoutineName, errorstatus, & + "Too many snow layers specified - increase snow_layers_max and recompile") +END IF + +! Copy the values for the given number of pfts across +IF ( .NOT. l_aggregate .AND. can_model == 4 ) THEN + canSnowTile(1:npft) = cansnowpft(1:npft) +END IF + +! Derive kland and dtland +kland = kland_numerator / dtland +tcland = tm - dtland + +! If variables are set check that they are in the correct range +IF ( ANY( ABS( aicemax(:) - rmdi ) > EPSILON(1.0) ) ) THEN + IF ( ANY( aicemax(:) < 0.01 ) .OR. ANY( aicemax(:) > 0.99 ) ) THEN + CALL ereport(RoutineName, errorstatus, & + 'aicemax must lie in the range 0.01 to 0.99') + END IF +END IF + +IF ( ABS( rho_firn_albedo - rmdi ) > EPSILON(1.0) ) THEN + IF ( rho_firn_albedo < 0.01 .OR. rho_firn_albedo > 1000.0 ) THEN + CALL ereport(RoutineName, errorstatus, & + 'rho_firn_albedo must lie in the range 0.01 to 1000.0') + END IF +END IF + +! Check for inconsistencies if a separate snow tile is used +IF ( ANY(i_snow_tile == 1) ) THEN + ! Tiles cannot be aggregated + IF ( l_aggregate ) THEN + CALL ereport(RoutineName, errorstatus, & + 'l_aggregate = F required if ANY i_snow_tile = 1') + END IF + ! The last surface type must be an ice tile + IF ( ice /= ntype ) THEN + CALL ereport(RoutineName, errorstatus, & + 'ice = ntype required if ANY i_snow_tile = 1') + END IF + ! The ice tile cannot be selected + IF ( i_snow_tile(ntype) == 1 ) THEN + CALL ereport(RoutineName, errorstatus, & + 'i_snow_tile(ntype) = 1 is not permitted') + END IF +END IF + +END SUBROUTINE check_jules_snow +#endif + +SUBROUTINE print_nlist_jules_snow() + +USE jules_print_mgr, ONLY: jules_print + +IMPLICIT NONE + +CHARACTER(LEN=50000) :: lineBuffer + + +!----------------------------------------------------------------------------- + + +CALL jules_print('jules_snow', & + 'Contents of namelist jules_snow') + +WRITE(lineBuffer, *) ' nsmax = ', nsmax +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' l_snowdep_surf = ', l_snowdep_surf +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' l_rho_snow_corr = ', l_rho_snow_corr +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' graupel_options = ', graupel_options +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' frac_snow_subl_melt = ', frac_snow_subl_melt +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' l_snow_infilt = ', l_snow_infilt +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' l_snow_nocan_hc = ', l_snow_nocan_hc +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' l_et_metamorph = ', l_et_metamorph +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' a_snow_et = ', a_snow_et +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' b_snow_et = ', b_snow_et +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' c_snow_et = ', c_snow_et +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' rho_snow_et_crit = ', rho_snow_et_crit +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' i_snow_cond_parm = ', i_snow_cond_parm +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' r0 = ', r0 +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' rmax = ', rmax +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snow_ggr = ', snow_ggr +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' amax = ', amax +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' aicemax = ', aicemax +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' maskd = ', maskd +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' dtland = ', dtland +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' kland_numerator = ', kland_numerator +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' rho_snow_const = ', rho_snow_const +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' rho_firn_albedo = ', rho_firn_albedo +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' rho_snow_fresh = ', rho_snow_fresh +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snow_hcon = ', snow_hcon +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snow_hcap = ', snow_hcap +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snowliqcap = ', snowliqcap +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snowinterceptfact = ', snowinterceptfact +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snowloadlai = ', snowloadlai +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' snowunloadfact = ', snowunloadfact +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' dzsnow = ', dzsnow +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer, *) ' cansnowpft = ', cansnowpft +CALL jules_print('jules_snow', lineBuffer) +WRITE(lineBuffer,*)' unload_rate_cnst = ', unload_rate_cnst +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,*)' unload_rate_u = ', unload_rate_u +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,*)' can_clump = ', can_clump +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,*)' lai_alb_lim_sn = ', lai_alb_lim_sn +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,*)' n_lai_exposed = ', n_lai_exposed +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,'(a,i3)')' i_grain_growth_opt = ', i_grain_growth_opt +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,'(a,i3)')' i_relayer_opt = ', i_relayer_opt +CALL jules_print('jules_snow', lineBuffer) + +WRITE(lineBuffer,'(a,i3)')' i_basal_melting_opt = ', i_basal_melting_opt +CALL jules_print('jules_snow', lineBuffer) + +CALL jules_print('jules_snow', & + '- - - - - - end of namelist - - - - - -') + +END SUBROUTINE print_nlist_jules_snow + +#if defined(UM_JULES) && !defined(LFRIC) +SUBROUTINE read_nml_jules_snow (unitnumber) + +! Description: +! Read the JULES_SNOW namelist + +USE setup_namelist, ONLY: setup_nml_type +USE check_iostat_mod, ONLY: check_iostat +USE UM_parcore, ONLY: mype +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook +USE errormessagelength_mod, ONLY: errormessagelength + +IMPLICIT NONE + +! Subroutine arguments +INTEGER, INTENT(IN) :: unitnumber + +INTEGER :: my_comm +INTEGER :: mpl_nml_type +INTEGER :: ErrorStatus +INTEGER :: icode +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='READ_NML_JULES_SNOW' +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 + +CHARACTER(LEN=errormessagelength) :: iomessage + +! set number of each type of variable in my_namelist type +INTEGER, PARAMETER :: no_of_types = 3 +INTEGER, PARAMETER :: n_int = 7 +INTEGER, PARAMETER :: n_real = 25 + snow_layers_max + 5 * npft_max +INTEGER, PARAMETER :: n_log = 5 + npft_max + +TYPE :: my_namelist + SEQUENCE + INTEGER :: nsmax + INTEGER :: frac_snow_subl_melt + INTEGER :: i_snow_cond_parm + INTEGER :: i_grain_growth_opt + INTEGER :: i_relayer_opt + INTEGER :: i_basal_melting_opt + INTEGER :: graupel_options + REAL(KIND=real_jlslsm) :: r0 + REAL(KIND=real_jlslsm) :: rmax + REAL(KIND=real_jlslsm) :: snow_ggr(3) + REAL(KIND=real_jlslsm) :: amax(2) + REAL(KIND=real_jlslsm) :: aicemax(2) + REAL(KIND=real_jlslsm) :: maskd + REAL(KIND=real_jlslsm) :: dtland + REAL(KIND=real_jlslsm) :: kland_numerator + REAL(KIND=real_jlslsm) :: rho_snow_const + REAL(KIND=real_jlslsm) :: rho_snow_fresh + REAL(KIND=real_jlslsm) :: rho_firn_albedo + REAL(KIND=real_jlslsm) :: snow_hcon + REAL(KIND=real_jlslsm) :: snow_hcap + REAL(KIND=real_jlslsm) :: snowliqcap + REAL(KIND=real_jlslsm) :: snowinterceptfact + REAL(KIND=real_jlslsm) :: snowloadlai + REAL(KIND=real_jlslsm) :: snowunloadfact + REAL(KIND=real_jlslsm) :: a_snow_et + REAL(KIND=real_jlslsm) :: b_snow_et + REAL(KIND=real_jlslsm) :: c_snow_et + REAL(KIND=real_jlslsm) :: rho_snow_et_crit + REAL(KIND=real_jlslsm) :: dzsnow(snow_layers_max) + REAL(KIND=real_jlslsm) :: unload_rate_cnst(npft_max) + REAL(KIND=real_jlslsm) :: unload_rate_u(npft_max) + REAL(KIND=real_jlslsm) :: can_clump(npft_max) + REAL(KIND=real_jlslsm) :: lai_alb_lim_sn(npft_max) + REAL(KIND=real_jlslsm) :: n_lai_exposed(npft_max) + LOGICAL :: l_snowdep_surf + LOGICAL :: l_rho_snow_corr + LOGICAL :: cansnowpft(npft_max) + LOGICAL :: l_snow_infilt + LOGICAL :: l_snow_nocan_hc + LOGICAL :: l_et_metamorph +END TYPE my_namelist + +TYPE (my_namelist) :: my_nml + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +CALL gc_get_communicator(my_comm, icode) + +CALL setup_nml_type(no_of_types, mpl_nml_type, n_int_in = n_int, & + n_real_in = n_real, n_log_in = n_log) + +IF (mype == 0) THEN + + READ (UNIT = unitnumber, NML = jules_snow, IOSTAT = errorstatus, & + IOMSG = iomessage) + CALL check_iostat(errorstatus, "namelist jules_snow", iomessage) + + my_nml % nsmax = nsmax + my_nml % frac_snow_subl_melt = frac_snow_subl_melt + my_nml % i_snow_cond_parm = i_snow_cond_parm + my_nml % graupel_options = graupel_options + my_nml % r0 = r0 + my_nml % rmax = rmax + my_nml % snow_ggr = snow_ggr + my_nml % amax = amax + my_nml % aicemax = aicemax + my_nml % maskd = maskd + my_nml % dtland = dtland + my_nml % kland_numerator = kland_numerator + my_nml % rho_snow_const = rho_snow_const + my_nml % rho_snow_fresh = rho_snow_fresh + my_nml % rho_firn_albedo = rho_firn_albedo + my_nml % snow_hcon = snow_hcon + my_nml % snow_hcap = snow_hcap + my_nml % snowliqcap = snowliqcap + my_nml % snowinterceptfact = snowinterceptfact + my_nml % snowloadlai = snowloadlai + my_nml % snowunloadfact = snowunloadfact + my_nml % a_snow_et = a_snow_et + my_nml % b_snow_et = b_snow_et + my_nml % c_snow_et = c_snow_et + my_nml % rho_snow_et_crit = rho_snow_et_crit + my_nml % dzsnow = dzsnow + my_nml % unload_rate_cnst = unload_rate_cnst + my_nml % unload_rate_u = unload_rate_u + my_nml % can_clump = can_clump + my_nml % lai_alb_lim_sn = lai_alb_lim_sn + my_nml % n_lai_exposed = n_lai_exposed + my_nml % l_snowdep_surf = l_snowdep_surf + my_nml % l_rho_snow_corr = l_rho_snow_corr + my_nml % cansnowpft = cansnowpft + my_nml % l_snow_infilt = l_snow_infilt + my_nml % l_snow_nocan_hc = l_snow_nocan_hc + my_nml % l_et_metamorph = l_et_metamorph + my_nml % i_grain_growth_opt = i_grain_growth_opt + my_nml % i_relayer_opt = i_relayer_opt + my_nml % i_basal_melting_opt = i_basal_melting_opt + +END IF + +CALL mpl_bcast(my_nml,1,mpl_nml_type,0,my_comm,icode) + +IF (mype /= 0) THEN + + nsmax = my_nml % nsmax + frac_snow_subl_melt = my_nml % frac_snow_subl_melt + i_snow_cond_parm = my_nml % i_snow_cond_parm + graupel_options = my_nml % graupel_options + r0 = my_nml % r0 + rmax = my_nml % rmax + snow_ggr = my_nml % snow_ggr + amax = my_nml % amax + aicemax = my_nml % aicemax + maskd = my_nml % maskd + dtland = my_nml % dtland + kland_numerator = my_nml % kland_numerator + rho_snow_const = my_nml % rho_snow_const + rho_snow_fresh = my_nml % rho_snow_fresh + rho_firn_albedo = my_nml % rho_firn_albedo + snow_hcon = my_nml % snow_hcon + snow_hcap = my_nml % snow_hcap + snowliqcap = my_nml % snowliqcap + snowinterceptfact = my_nml % snowinterceptfact + snowloadlai = my_nml % snowloadlai + snowunloadfact = my_nml % snowunloadfact + a_snow_et = my_nml % a_snow_et + b_snow_et = my_nml % b_snow_et + c_snow_et = my_nml % c_snow_et + rho_snow_et_crit = my_nml % rho_snow_et_crit + dzsnow = my_nml % dzsnow + unload_rate_cnst = my_nml % unload_rate_cnst + unload_rate_u = my_nml % unload_rate_u + can_clump = my_nml % can_clump + lai_alb_lim_sn = my_nml % lai_alb_lim_sn + n_lai_exposed = my_nml % n_lai_exposed + l_snowdep_surf = my_nml % l_snowdep_surf + l_rho_snow_corr = my_nml % l_rho_snow_corr + cansnowpft = my_nml % cansnowpft + l_snow_infilt = my_nml % l_snow_infilt + l_snow_nocan_hc = my_nml % l_snow_nocan_hc + l_et_metamorph = my_nml % l_et_metamorph + i_grain_growth_opt = my_nml % i_grain_growth_opt + i_relayer_opt = my_nml % i_relayer_opt + i_basal_melting_opt = my_nml % i_basal_melting_opt +END IF + +CALL mpl_type_free(mpl_nml_type,icode) + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE read_nml_jules_snow +#endif + +END MODULE jules_snow_mod diff --git a/src/io/model_interface/surf_couple_explicit_mod.F90 b/src/io/model_interface/surf_couple_explicit_mod.F90 new file mode 100644 index 00000000..6f922834 --- /dev/null +++ b/src/io/model_interface/surf_couple_explicit_mod.F90 @@ -0,0 +1,848 @@ +! *****************************COPYRIGHT**************************************** +! (c) Crown copyright, Met Office. All rights reserved. +! +! This routine has been licensed to the other JULES partners for use and +! distribution under the JULES collaboration agreement, subject to the terms and +! conditions set out therein. +! +! [Met Office Ref SC0237] +! *****************************COPYRIGHT**************************************** + +MODULE surf_couple_explicit_mod + +USE jules_gridinit_sf_explicit_mod, ONLY: jules_gridinit_sf_explicit +USE jules_land_sf_explicit_mod, ONLY: jules_land_sf_explicit +USE jules_ssi_sf_explicit_mod, ONLY: jules_ssi_sf_explicit +USE jules_griddiag_sf_explicit_mod, ONLY: jules_griddiag_sf_explicit + +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +PRIVATE +PUBLIC :: surf_couple_explicit + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='SURF_COUPLE_EXPLICIT_MOD' + +CONTAINS + +!=============================================================================== +! Public subroutine +!=============================================================================== +SUBROUTINE surf_couple_explicit( & + !Arguments used by JULES-standalone + !Misc INTENT(IN) DONE + bq_1, bt_1, photosynth_act_rad, & + curr_day_number, curr_time, & + !Diagnostics, INTENT(INOUT) + sf_diag, & + !Fluxes INTENT(OUT) DONE + fqw_1,ftl_1, & + !Misc INTENT(OUT) + radnet_sice, rhokm_1, rhokm_land, rhokm_ssi, & + !Out of explicit and into implicit only INTENT(OUT) + cdr10m, & + alpha1, alpha1_sea, alpha1_sice, ashtf_prime, ashtf_prime_sea, & + ashtf_prime_surft, epot_surft, & + fracaero_t, fracaero_s, resfs, resft, & + rhokh, rhokh_surft, rhokh_sice, rhokh_sea, & + dtstar_surft, dtstar_sea, dtstar_sice, & + z0hssi, z0mssi, chr1p5m, chr1p5m_sice, canhc_surft, & + wt_ext_surft, flake, & + !Out of explicit and into extra only INTENT(OUT) + hcons_soilt, & + !Out of explicit and into implicit and extra INTENT(OUT) + tile_frac, & + !Additional arguments for the UM----------------------------------------- + !JULES prognostics module + ti_cat_sicat, & + !JULES ancil_info module + land_pts, ssi_pts, sea_pts, nsurft, surft_pts, & + ! IN input data from the wave model + charnock_w, & + !JULES coastal module + flandg, & + !JULES aero module + co2_3d_ij, & + !JULES trifctl module + asteps_since_triffid, & + !JULES p_s_parms module + soil_clay_ij, & + !JULES switches module + l_spec_z0, & + !Not in a JULES module + numcycles, cycleno, z1_uv_top, z1_tq_top, sky, ddmfx, & + l_aero_classic, z0m_scm, z0h_scm, recip_l_mo_sea, rib, & + flandfac, fseafac, fb_surf, u_s, t1_sd, q1_sd, rhostar, & + vshr, resp_s_tot_soilt, emis_soil, & + !TYPES containing field data (IN OUT) + crop_vars,psparms,ainfo,trif_vars,aerotype,urban_param,progs, & + trifctltype,coast, jules_vars, & + fluxes, & + lake_vars, & + forcing, & + !rivers, & + !veg3_parm, & + !veg3_field, & + chemvars, & + wtrac_jls, & + progs_cbl, work_cbl & + ) + +!Module imports + +!TYPE definitions +USE crop_vars_mod, ONLY: crop_vars_type +USE p_s_parms, ONLY: psparms_type +USE ancil_info, ONLY: ainfo_type +USE trif_vars_mod, ONLY: trif_vars_type +USE aero, ONLY: aero_type +USE urban_param_mod, ONLY: urban_param_type +USE prognostics, ONLY: progs_type +USE trifctl, ONLY: trifctl_type +USE coastal, ONLY: coastal_type +USE jules_vars_mod, ONLY: jules_vars_type +USE fluxes_mod, ONLY: fluxes_type +USE lake_mod, ONLY: lake_type +USE jules_forcing_mod, ONLY: forcing_type +! USE jules_rivers_mod, ONLY: rivers_type +! USE veg3_parm_mod, ONLY: in_dev +! USE veg3_field_mod, ONLY: in_dev +USE jules_chemvars_mod, ONLY: chemvars_type +USE jules_wtrac_type_mod, ONLY: jls_wtrac_type + +! In general CABLE utilizes a required subset of tbe JULES types, however; +USE progs_cbl_vars_mod, ONLY: progs_cbl_vars_type ! CABLE requires extra progs +USE work_vars_mod_cbl, ONLY: work_vars_type ! and some kept thru timestep + +!Common modules +USE ereport_mod, ONLY: & + ereport +USE jules_sea_seaice_mod, ONLY: & + nice, nice_use +USE jules_soil_mod, ONLY: & + sm_levels +USE jules_vegetation_mod, ONLY: & + l_phenol +USE jules_surface_types_mod, ONLY: & + npft +USE sf_diags_mod, ONLY: strnewsfdiag +USE atm_fields_bounds_mod, ONLY: & + pdims_s, pdims, tdims +USE ancil_info, ONLY: & + nsoilt +USE jules_water_tracers_mod, ONLY: & + n_wtrac_jls, n_evap_srce + + !New arguments replacing USE statements +USE veg3_field_mod, ONLY: & + veg_state + +! Modules used by UM only +#if defined(UM_JULES) +USE jules_soil_biogeochem_mod, ONLY: & +! imported scalar parameters + soil_model_4pool, & +! imported scalar variables (IN) + soil_bgc_model +#endif + +!Modules that change name between JULES and UM +#if defined(UM_JULES) +USE atm_step_local, ONLY: & + co2_dim_len,co2_dim_row +USE carbon_options_mod, ONLY: & + l_co2_interactive +USE dust_parameters_mod, ONLY: & + l_dust, l_dust_diag +! USE s_main_force, ONLY: & !A SCM module currently doesn't work +! l_spec_z0 +USE rad_input_mod, ONLY: & + co2_mmr +USE gen_phys_inputs_mod, ONLY: & + l_mr_physics +USE ancil_info, ONLY: & + dim_cslayer +#else +USE ancil_info, ONLY: & + co2_dim_len,co2_dim_row, dim_cslayer +USE switches, ONLY: & + l_co2_interactive, l_dust, l_mr_physics +USE aero, ONLY: & + co2_mmr + +#endif + +! reads the lsm switch set in a namelist +USE jules_model_environment_mod, ONLY: lsm_id, jules, cable + +! for testing lsm switch +USE jules_print_mgr, ONLY: jules_message, jules_print + +!Dr Hook +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +IMPLICIT NONE + +!----------------------------------------------------------------------------- +! Description: +! Coupling routine between the UM or JULES system code and land surface +! explicit science routines. Calls the appropriate LSM-specific code. +! +! Code Owner: Please refer to ModuleLeaders.txt +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +! Subroutine arguments + +!UM-only arguments +!JULES ancil_info module +INTEGER, INTENT(IN) :: & + land_pts, & + nsurft, ssi_pts, sea_pts +REAL(KIND=real_jlslsm), INTENT(IN) :: & + charnock_w(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) +INTEGER, INTENT(OUT) :: & + surft_pts(nsurft) + +!JULES prognostics module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + ti_cat_sicat(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice) +!JULES coastal module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + flandg(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end) +!JULES aero module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + co2_3d_ij(co2_dim_len, co2_dim_row) +!JULES trifctl module +INTEGER, INTENT(IN):: & + asteps_since_triffid +!JULES p_s_parms module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + soil_clay_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + +!JULES switches module +LOGICAL, INTENT(IN) :: & + l_spec_z0 +!Not in a JULES module +INTEGER, INTENT(IN) :: numcycles, cycleno +REAL(KIND=real_jlslsm), INTENT(IN) :: & + !These variables are INTENT(IN) to sf_expl, but not used with the + !current configuration of standalone JULES (initialised to 0 below) + z1_uv_top(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + ! Height of top of lowest uv-layer + z1_tq_top(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + ! Height of top of lowest Tq-layer + sky(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + ! Skyview correction factor for surface LW + ddmfx(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + ! Convective downdraught mass-flux at cloud base +LOGICAL, INTENT(IN) :: l_aero_classic +REAL(KIND=real_jlslsm), INTENT(IN) :: & + !These variables are required for prescribed roughness lengths in + !SCM mode in UM - not used standalone + z0m_scm(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Fixed Sea-surface roughness length for momentum (m).(SCM) + z0h_scm(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + !Fixed Sea-surface roughness length for heat (m). (SCM) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + recip_l_mo_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Reciprocal of the surface Obukhov length at sea points. (m-1). + rib(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Mean bulk Richardson number for lowest layer. + flandfac(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & + fseafac(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & + fb_surf(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface flux buoyancy over density (m^2/s^3) + u_s(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface friction velocity (m/s) + t1_sd(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Standard deviation of turbulent fluctuations of layer 1 temp; used in + !initiating convection. + q1_sd(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Standard deviation of turbulent flux of layer 1 humidity; used in + !initiating convection. + rhostar(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface air density + vshr(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Magnitude of surface-to-lowest atm level wind shear (m per s). + resp_s_tot_soilt(land_pts,nsoilt), & + !Total soil respiration (kg C/m2/s). + emis_soil(land_pts) + !Emissivity of underlying soil + +!Misc INTENT(IN) +REAL(KIND=real_jlslsm), INTENT(IN) :: & + bq_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !A buoyancy parameter (beta q tilde). + bt_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !A buoyancy parameter (beta T tilde). + photosynth_act_rad(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + !Net downward shortwave radiation in band 1 (w/m2). +INTEGER, INTENT(IN) :: & + curr_day_number, curr_time + +!diagnostic array +TYPE (strnewsfdiag), INTENT(IN OUT) :: sf_diag + +!Fluxes INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + fqw_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Moisture flux between layers (kg per square metre per sec). + ftl_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + !FTL(,K) contains net turbulent sensible heat flux into layer K from below; + !so FTL(,1) is the surface sensible heat, H.(W/m2) + +!Misc INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + radnet_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Surface net radiation on sea-ice (W/m2) + rhokm_1(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & + !Exchange coefficients for momentum on P-grid + rhokm_land(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & + rhokm_ssi(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end) + +!Out of explicit and into implicit only INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + cdr10m(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & +! Interpolation coefficient for diagnosis of 10 m winds + alpha1(land_pts,nsurft), & + !Mean gradient of saturated specific humidity with respect to temperature + !between the bottom model layer and tile surfaces + alpha1_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !ALPHA1 for sea-ice. + alpha1_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !ALPHA1 for sea. + ashtf_prime(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Coefficient to calculate surface heat flux into sea-ice. + ashtf_prime_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Coefficient to calculate surface heat flux into sea. + ashtf_prime_surft(land_pts,nsurft), & + !Coefficient to calculate surface heat flux into land tiles. + epot_surft(land_pts,nsurft), & + !Local EPOT for land tiles. + fracaero_t(land_pts,nsurft), & + !Total fraction of surface moisture flux with only aerodynamic resistance + fracaero_s(land_pts,nsurft), & + !Fraction of surface moisture flux with only aerodynamic resistance + !from the frozen part of the tile only + resfs(land_pts,nsurft), & + !Combined soil, stomatal and aerodynamic resistance factor for fraction + !(1-fracaero_t) of snow-free land tiles. + resft(land_pts,nsurft), & + !Total resistance factor. fracaero_t+(1-fracaero_t)*resfs & + !for snow-free land, 1 for snow. + rhokh(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Grid-box surface exchange coefficients + rhokh_surft(land_pts,nsurft), & + !Surface exchange coefficients for land tiles + rhokh_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Surface exchange coefficients for sea-ice + rhokh_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface exchange coefficients for sea + dtstar_surft(land_pts,nsurft), & + !Change in TSTAR over timestep for land tiles + dtstar_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Change is TSTAR over timestep for open sea + dtstar_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Change is TSTAR over timestep for sea-ice + z0hssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Roughness length for heat and moisture over sea (m). + z0mssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Roughness length for momentum over sea (m). + chr1p5m(land_pts,nsurft), & + !Ratio of coefffs for calculation of 1.5m temp for land tiles. + chr1p5m_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !CHR1P5M for sea and sea-ice (leads ignored). + canhc_surft(land_pts,nsurft), & + !Areal heat capacity of canopy for land tiles (J/K/m2). + wt_ext_surft(land_pts,sm_levels,nsurft), & + !Fraction of evapotranspiration which is extracted from each soil layer + !by each tile. + flake(land_pts,nsurft) + !Lake fraction. + +!Out of explicit and into extra only INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + hcons_soilt(land_pts,nsoilt) + !Soil thermal conductivity including water and ice + +!Out of explicit and into implicit and extra INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + tile_frac(land_pts,nsurft) + !Tile fractions including snow cover in the ice tile. + +#if !defined(UM_JULES) +LOGICAL :: l_dust_diag + !In standalone, this switch essentially does the same job as l_dust +#endif + +!TYPES containing field data (IN OUT) +TYPE(crop_vars_type), INTENT(IN OUT) :: crop_vars +TYPE(psparms_type), INTENT(IN OUT) :: psparms +TYPE(ainfo_type), INTENT(IN OUT) :: ainfo +TYPE(trif_vars_type), INTENT(IN OUT) :: trif_vars +TYPE(aero_type), INTENT(IN OUT) :: aerotype +TYPE(urban_param_type), INTENT(IN OUT) :: urban_param +TYPE(progs_type), INTENT(IN OUT) :: progs +TYPE(trifctl_type), INTENT(IN OUT) :: trifctltype +TYPE(coastal_type), INTENT(IN OUT) :: coast +TYPE(jules_vars_type), INTENT(IN OUT) :: jules_vars +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes +TYPE(lake_type), INTENT(IN OUT) :: lake_vars +TYPE(forcing_type), INTENT(IN OUT) :: forcing +!TYPE(rivers_type), INTENT(IN OUT) :: rivers +!TYPE(in_dev), INTENT(IN OUT) :: veg3_parm +!TYPE(in_dev), INTENT(IN OUT) :: veg3_field +TYPE(chemvars_type), INTENT(IN OUT) :: chemvars +TYPE(jls_wtrac_type), INTENT(IN OUT) :: wtrac_jls + +!CABLE TYPES containing field data (IN OUT) +TYPE(progs_cbl_vars_type), INTENT(IN OUT) :: progs_cbl +TYPE(work_vars_type), INTENT(IN OUT) :: work_cbl + +INTEGER :: & + i,j,l,n !Various counters +INTEGER :: errorstatus + +!----------------------------------------------------------------------------- +! Local variables +!----------------------------------------------------------------------------- + +! Workspace for sea and sea-ice leads +REAL(KIND=real_jlslsm) :: & + cd_ssi_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Bulk transfer coefficient for + ! momentum over sea mean. +,ch_ssi_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Bulk transfer coefficient for heat + ! and/or moisture over sea mean. +,rhokh_1_sice_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Surface exchange coefficient for + ! sea-ice. +,rib_sea_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Bulk Richardson number +,z0h_sea_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + ! Roughness length for heat and + ! moisture transport + +! Workspace for sea-ice and marginal ice zone +REAL(KIND=real_jlslsm) :: & + cd_land_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Bulk transfer coefficient for + ! momentum over land. +,rib_ice_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end, & + nice_use) & + ! Bulk Richardson number +,rib_surft(land_pts,nsurft) & + ! RIB for land tiles. +,z0m_ice_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end, & + nice_use) & + ! Momentum Roughness length. +,z0h_ice_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end, & + nice_use) & + ! Thermal Roughness length. +,ice_fract_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + ! Sea ice fraction summed over all cats + +! Temporary array for flandg variable +REAL(KIND=real_jlslsm) :: & + flandg_tmp(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end) + +! Workspace for land tiles +REAL(KIND=real_jlslsm) :: & + ch_surft_classic(land_pts,nsurft) & + ! Bulk transfer coefficient for + ! heat for aerosol deposition. +,cd_std_classic(land_pts,nsurft) + ! Bulk transfer coefficient for + ! momentum for aerosol deposition. +! Workspace for air density calculation +REAL(KIND=real_jlslsm) :: & + rhostar_land(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + rhostar_ssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + +! Workspace for variables passed between land and sea and sea-ice routines +LOGICAL :: & + l_cdr10m_snow + ! IN Flag indicating if cdr10m + ! (an interpolation coefficient) is + ! to be calculated for use with + ! snow unloading. + +! Temp until the model switching is implemented for coupled mode +! Having a parameter until then should hopefully help the compiler eliminate +! dead code + +!Dr Hook variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='SURF_COUPLE_EXPLICIT' + +!----------------------------------------------------------------------------- +!End of header +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +SELECT CASE( lsm_id ) +CASE ( jules ) + + ! Change 2d UM clay to 1d jules clay content for soil respiration + ! Soil tiling not currently in the UM, so broadcast ij value to all tiles. + ! Multi-layer clay not currently in UM so set all layers to same value. +#if defined(UM_JULES) + IF ( soil_bgc_model == soil_model_4pool ) THEN + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / pdims%i_end + 1 + i = ainfo%land_index(l) - (j-1) * pdims%i_end + DO n = 1, dim_cslayer + psparms%clay_soilt(l,1,n) = soil_clay_ij(i,j) + END DO + END DO + END IF +#else + l_dust_diag = l_dust +#endif + + CALL jules_gridinit_sf_explicit ( & + !IN soil/vegetation/land surface data : + flandg, & + !IN everything not covered so far : + forcing%pstar_ij,fluxes%tstar_ij, & + !IN variables for message passing + jules_vars%u_1_p_ij, jules_vars%v_1_p_ij, & + jules_vars%u_0_p_ij, jules_vars%v_0_p_ij, & + ! INOUT + sf_diag, & + !OUT Diagnostic not requiring STASH flags : + fqw_1,ftl_1, & + !OUT variables for message passing + flandfac, fseafac, cdr10m, & + !OUT data required for mineral dust scheme + t1_sd,q1_sd, & + !OUT data required elsewhere in boundary layer or surface code + rhostar,vshr,coast%vshr_land_ij,coast%vshr_ssi_ij,coast%mapi, coast%mapj & + ) + + DO i = tdims%i_start,tdims%i_end + DO j = tdims%j_start,tdims%j_end + rhostar_land(i,j) = rhostar(i,j) + rhostar_ssi(i,j) = rhostar(i,j) + END DO + END DO + + ! Populate temporary flandg array + flandg_tmp(:,:) = flandg(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end) + + CALL jules_land_sf_explicit ( & + !IN date- and time-related values + curr_day_number, curr_time, & + !IN values defining field dimensions and subset to be processed : + land_pts, & + !IN values defining water tracer field dimensions + n_wtrac_jls, n_evap_srce, & + !IN parameters for iterative SISL scheme + numcycles, cycleno, & + !IN parameters required from boundary-layer scheme : + bq_1,bt_1,ainfo%z1_uv_ij,z1_uv_top,ainfo%z1_tq_ij,z1_tq_top, & + forcing%qw_1_ij,forcing%tl_1_ij, & + !IN soil/vegetation/land surface data : + ainfo%land_index,nsurft,sm_levels,progs%canopy_surft,psparms%catch_surft, & + psparms%catch_snow_surft, psparms%hcon_soilt,jules_vars%ho2r2_orog_gb, & + flandg_tmp, fluxes%fsnow, & + jules_vars%sil_orog_land_gb,psparms%smvccl_soilt, & + psparms%smvcst_soilt,psparms%smvcwt_soilt, & + psparms%sthf_soilt, psparms%sthu_soilt,psparms%z0_surft, & + psparms%z0h_bare_surft, psparms%z0m_soil_gb, & + !IN input data from the wave model + charnock_w, & + !IN everything not covered so far : + forcing%pstar_ij,forcing%lw_down_ij,fluxes%sw_surft,sky,jules_vars%zh, & + ddmfx, co2_mmr,co2_3d_ij,l_co2_interactive,l_phenol, & + asteps_since_triffid,progs%cs_pool_soilt,veg_state,ainfo%frac_surft, & + progs%canht_pft, & + photosynth_act_rad, progs%lai_pft, & + l_mr_physics,progs%t_soil_soilt,progs%tsurf_elev_surft, & + progs%tstar_surft,jules_vars%z_land_ij, & + psparms%albsoil_soilt, & + psparms%cosz_ij, & + l_aero_classic,l_dust,l_dust_diag,psparms%clay_soilt,chemvars%o3_gb, & + fluxes%l_emis_surft_set, & + ainfo%latitude,ainfo%longitude, & + !INOUT diagnostics + sf_diag, & + !INOUT data : + fluxes%emis_surft,progs%gs_gb,chemvars%gc_corr,trifctltype%g_leaf_acc_pft, & + trifctltype%npp_acc_pft,trifctltype%resp_w_acc_pft, & + trifctltype%resp_s_acc_soilt,rhostar_land, & + !INOUT Diagnostic not requiring STASH flags : + fqw_1,ftl_1,t1_sd,q1_sd,vshr,coast%vshr_land_ij, & + !OUT Diagnostic not requiring STASH flags : + fluxes%ftl_surft, & + !OUT variables for message passing + rhokm_land, cdr10m, & + !OUT data required for mineral dust scheme + aerotype%u_s_std_surft, & + !OUT data required elsewhere in boundary layer or surface code + alpha1,ashtf_prime_surft,fluxes%fqw_surft,epot_surft,fracaero_t, & + fracaero_s(:,:), & + resfs,resft,rhokh_surft,dtstar_surft,fluxes%z0h_surft,fluxes%z0m_surft, & + chr1p5m,progs%smc_soilt,hcons_soilt,trifctltype%gpp_gb,trifctltype%npp_gb, & + trifctltype%resp_p_gb,trifctltype%g_leaf_pft, & + trifctltype%gpp_pft,trifctltype%npp_pft,trifctltype%resp_p_pft, & + trifctltype%resp_s_soilt,resp_s_tot_soilt, & + trif_vars%resp_l_pft,trif_vars%resp_r_pft,trifctltype%resp_w_pft, & + trif_vars%n_leaf_pft,trif_vars%n_root_pft,trif_vars%n_stem_pft, & + trif_vars%lai_bal_pft, & + progs%gc_surft,canhc_surft,wt_ext_surft,flake, & + ainfo%surft_index,surft_pts,tile_frac,fluxes%fsmc_pft,emis_soil, & + trifctltype%growth_sug_pft,trifctltype%growth_sug_gb,progs%f_nsc_pft, & + trifctltype%lwp_c_pft, trifctltype%psi_root_zone_pft, & + ! OUT required for classic aerosols + cd_land_ij,rib_surft,ch_surft_classic,cd_std_classic, & + ! OUT required for sea and sea-ice calculations + l_cdr10m_snow, & + !New arguments replacing USE statements + ! Prognostics (IN) + progs%t_home_gb, progs%t_growth_gb, & + !urban_param (IN) + urban_param%emisr_gb, urban_param%emisw_gb, urban_param%hwr_gb, & + !jules_vars_mod (IN OUT) + jules_vars%albobs_scaling_surft, & + !jules_chemvars_mod (OUT) + chemvars%isoprene_gb, chemvars%isoprene_pft, chemvars%terpene_gb , & + chemvars%terpene_pft, chemvars%methanol_gb, chemvars%methanol_pft, & + chemvars%acetone_gb, chemvars%acetone_pft, & + !trif_vars_mod (OUT) + trif_vars%fapar_diag_pft, trif_vars%apar_diag_pft, trif_vars%apar_diag_gb, & + trif_vars%gpp_gb_acc, trif_vars%gpp_pft_acc, & + !crop_vars_mod (IN) + crop_vars%rootc_cpft, crop_vars%sthu_irr_soilt, & + crop_vars%frac_irr_soilt, crop_vars%frac_irr_surft, crop_vars%dvi_cpft, & + !crop_vars_mod (IN OUT) + crop_vars%resfs_irr_surft, & + !crop_vars_mod (OUT) + crop_vars%gs_irr_surft, crop_vars%smc_irr_soilt, & + crop_vars%wt_ext_irr_surft, crop_vars%gc_irr_surft, & + !p_s_parms (IN) + psparms%bexp_soilt, psparms%sathh_soilt, psparms%v_close_pft, & + psparms%v_open_pft, & + !urban_param (IN) + urban_param%wrr_gb, & + !Fluxes (IN OUT) + fluxes%anthrop_heat_surft, & + !snow prognostics (IN) + progs%nsnow_surft, progs%sice_surft, progs%sliq_surft, & + progs%tsnow_surft, progs%ds_surft, & + !snow prognostics (IN OUT) + progs%snowdepth_surft, progs%snow_surft, & + !c_elevate (OUT) + jules_vars%surf_hgt_surft, jules_vars%lw_down_elevcorr_surft, & + !jules_vars_mod (OUT) + jules_vars%snowdep_surft, & + !urban_param (IN) + urban_param%hgt_gb, urban_param%disp_gb, & + !lake_mod (IN) + lake_vars%lake_t_ice_gb, lake_vars%lake_t_mxl_gb, lake_vars%lake_h_ice_gb, & + lake_vars%lake_depth_gb, lake_vars%g_dt_gb, lake_vars%non_lake_frac, & + !lake_mod (OUT) + lake_vars%nusselt_gb, lake_vars%ts1_lake_gb, lake_vars%hcon_lake, & + !ancil_info (IN) + ainfo%l_lice_point, ainfo%l_soil_point, ainfo%l_lice_surft, & + !jules_surface_types (IN) + jules_vars%diff_frac, & + !chemvars (OUT) + chemvars%flux_o3_pft, chemvars%fo3_pft, & + !Water tracers (IN) + wtrac_jls%snow_surft, wtrac_jls%canopy_surft, & + wtrac_jls%sthu_soilt, wtrac_jls%qw_1_ij, & + !Water tracers (OUT) + wtrac_jls%fqw_1, wtrac_jls%fqw_surft, wtrac_jls%fqw_evapsrce, & + wtrac_jls%smc_soilt & + ) + + CALL jules_ssi_sf_explicit ( & + !IN values defining field dimensions and subset to be processed : + nice, nice_use, n_wtrac_jls, & + !IN parameters required from boundary-layer scheme : + bq_1,bt_1,ainfo%z1_uv_ij,z1_uv_top,ainfo%z1_tq_ij,z1_tq_top, & + forcing%qw_1_ij,forcing%tl_1_ij, & + !IN soil/vegetation/land surface data : + flandg_tmp, & + !IN sea/sea-ice data : + ainfo%ice_fract_ncat_sicat,progs%k_sice_sicat, & + progs%z0m_sice_fmd, progs%z0m_sice_skin, & + !IN input data from the wave model + charnock_w, & + !IN everything not covered so far : + forcing%pstar_ij,forcing%lw_down_ij,jules_vars%zh,ddmfx, & + l_mr_physics,progs%ti_sicat, & + ainfo%ti_cat_sicat, & + coast%tstar_sea_ij,coast%tstar_sice_sicat,jules_vars%z_land_ij, & + !IN idealised and SCM things + l_spec_z0, z0m_scm, z0h_scm, & + ! IN calculated in land code + l_cdr10m_snow, & + !INOUT diagnostics + sf_diag, & + !INOUT data : + progs%z0msea_ij,rhostar_ssi,fqw_1,ftl_1,t1_sd,q1_sd,cdr10m, & + !OUT Diagnostic not requiring STASH flags : + recip_l_mo_sea,radnet_sice, & + !OUT variables for message passing + rhokm_ssi, & + !OUT data required elsewhere in boundary layer or surface code + alpha1_sea,alpha1_sice,ashtf_prime,ashtf_prime_sea,fluxes%fqw_sicat, & + fluxes%ftl_sicat, & + rhokh_sice,rhokh_sea,dtstar_sea,dtstar_sice,chr1p5m_sice, & + coast%vshr_ssi_ij, & + ! OUT required for classic aerosols + cd_ssi_ij,ch_ssi_ij,rhokh_1_sice_ij,rib_sea_ij,z0h_sea_ij, & + rib_ice_ij,z0m_ice_ij,z0h_ice_ij,ice_fract_ij, & + !New arguments replacing USE statements + !ancil_info (IN) + ainfo%ssi_index, ainfo%sea_index, ainfo%fssi_ij, ainfo%sea_frac, & + ainfo%sice_index_ncat, ainfo%sice_frac_ncat, ainfo%sice_pts_ncat, & + ssi_pts, sea_pts, & + !fluxes (IN) + fluxes%sw_sicat, fluxes%sw_sea, & + !Water tracers (IN) + wtrac_jls%qw_1_ij, wtrac_jls%r_sea, wtrac_jls%r_sice, & + !Water tracers (INOUT) + wtrac_jls%fqw_1, & + !Water tracers (OUT) + wtrac_jls%fqw_sea, wtrac_jls%fqw_sicat & + ) + + CALL jules_griddiag_sf_explicit ( & + !IN values defining field dimensions and subset to be processed : + land_pts, nice_use, & + !IN parameters required from boundary-layer scheme : + bq_1,bt_1,ainfo%z1_uv_ij, & + !IN soil/vegetation/land surface data : + ainfo%land_index,nsurft,jules_vars%ho2r2_orog_gb, & + flandg_tmp, & + jules_vars%sil_orog_land_gb, & + !IN everything not covered so far : + forcing%pstar_ij,rhostar_land,rhostar_ssi,jules_vars%zh, & + fluxes%tstar_ij,l_aero_classic, & + !INOUT data : + progs%z0msea_ij,cdr10m, & + !IN Diagnostic not requiring STASH flags : + fqw_1,ftl_1, & + !IN variables for message passing + rhokm_land, rhokm_ssi, & + !IN data required elsewhere in boundary layer or surface code + rhokh_surft,rhokh_sice,rhokh_sea,fluxes%z0h_surft,fluxes%z0m_surft, & + coast%vshr_land_ij,coast%vshr_ssi_ij,ainfo%surft_index,surft_pts,tile_frac,& + ! IN required for classic aerosols + cd_ssi_ij,ch_ssi_ij,rhokh_1_sice_ij,rib_sea_ij,z0h_sea_ij, & + cd_land_ij,rib_ice_ij,rib_surft,z0m_ice_ij,z0h_ice_ij,ice_fract_ij, & + ch_surft_classic,cd_std_classic, & + ! INOUT data + rhostar, & + !INOUT diagnostics + sf_diag, & + !OUT Diagnostic not requiring STASH flags : + rhokm_1,rib, & + !OUT data required for tracer mixing : + aerotype%rho_aresist_ij,aerotype%aresist_ij,aerotype%resist_b_ij, & + aerotype%rho_aresist_surft,aerotype%aresist_surft,aerotype%resist_b_surft, & + !OUT data required for mineral dust scheme + aerotype%r_b_dust_ij,aerotype%cd_std_dust_ij, & + !OUT data required elsewhere in UM system : + fb_surf,u_s,jules_vars%blendht, & + !OUT data required elsewhere in boundary layer or surface code + rhokh,jules_vars%h_blend_orog_ij,z0hssi,z0mssi,jules_vars%z0m_eff_ij, & + !ancil_info (IN) + ainfo%ssi_index, ainfo%sice_index_ncat, ainfo%sice_frac_ncat, & + ainfo%sice_pts_ncat, & + !jules_internal + jules_vars%unload_backgrnd_pft & + ) + +CASE ( cable ) + ! for testing LSM switch + WRITE(jules_message,'(A)') "CABLE not yet implemented" + CALL jules_print(RoutineName, jules_message) + + ! initialise all INTENT(OUT) for now until CABLE is implemented + fqw_1(:,:) = 0.0 + ftl_1(:,:) = 0.0 + fluxes%ftl_surft(:,:) = 0.0 + fluxes%fqw_surft(:,:) = 0.0 + fluxes%fqw_sicat(:,:,:) = 0.0 + fluxes%ftl_sicat(:,:,:) = 0.0 + fluxes%fsmc_pft(:,:) = 0.0 + fluxes%emis_surft(:,:) = 0.0 + radnet_sice(:,:,:) = 0.0 + rhokm_1(:,:) = 0.0 + rhokm_land(:,:) = 0.0 + rhokm_ssi(:,:) = 0.0 + cdr10m(:,:) = 0.0 + alpha1(:,:) = 0.0 + alpha1_sea(:,:) = 0.0 + alpha1_sice(:,:,:) = 0.0 + ashtf_prime(:,:,:) = 0.0 + ashtf_prime_sea(:,:) = 0.0 + ashtf_prime_surft(:,:) = 0.0 + epot_surft(:,:) = 0.0 + fracaero_t(:,:) = 0.0 + fracaero_s(:,:) = 0.0 + resfs(:,:) = 0.0 + resft(:,:) = 0.0 + rhokh(:,:) = 0.0 + rhokh_surft(:,:) = 0.0 + rhokh_sice(:,:,:) = 0.0 + rhokh_sea(:,:) = 0.0 + dtstar_surft(:,:) = 0.0 + dtstar_sea(:,:) = 0.0 + dtstar_sice(:,:,:) = 0.0 + z0hssi(:,:) = 0.0 + fluxes%z0h_surft(:,:) = 0.0 + z0mssi(:,:) = 0.0 + fluxes%z0m_surft(:,:) = 0.0 + chr1p5m(:,:) = 0.0 + chr1p5m_sice(:,:) = 0.0 + canhc_surft(:,:) = 0.0 + wt_ext_surft(:,:,:) = 0.0 + flake(:,:) = 0.0 + hcons_soilt(:,:) = 0.0 + tile_frac(:,:) = 0.0 + + progs_cbl%SoilTemp_CABLE(:,:,:) = 0.0 + progs_cbl%SoilMoisture_CABLE(:,:,:) = 0.0 + progs_cbl%FrozenSoilFrac_CABLE(:,:,:) = 0.0 + progs_cbl%SnowDepth_CABLE(:,:,:) = 0.0 + progs_cbl%SnowMass_CABLE(:,:,:) = 0.0 + progs_cbl%SnowTemp_CABLE(:,:,:) = 0.0 + progs_cbl%SnowDensity_CABLE(:,:,:) = 0.0 + progs_cbl%ThreeLayerSnowFlag_CABLE(:,:) = 0.0 + progs_cbl%OneLyrSnowDensity_CABLE(:,:) = 0.0 + progs_cbl%SnowAge_CABLE(:,:) = 0.0 + progs_cbl%snowOsurft(:,:) = 0.0 + + work_cbl%snow_tile(:,:) = 0.0 + work_cbl%lying_snow(:) = 0.0 + work_cbl%surf_roff(:) = 0.0 + work_cbl%sub_surf_roff(:) = 0.0 + work_cbl%tot_tfall(:) = 0.0 + +CASE DEFAULT + errorstatus = 101 + WRITE(jules_message,'(A,I0)') 'Unrecognised surface scheme. lsm_id = ', & + lsm_id + CALL ereport(RoutineName, errorstatus, jules_message) + +END SELECT + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE surf_couple_explicit + +END MODULE surf_couple_explicit_mod diff --git a/src/io/model_interface/surf_couple_extra_mod.F90 b/src/io/model_interface/surf_couple_extra_mod.F90 new file mode 100644 index 00000000..2fb70b3a --- /dev/null +++ b/src/io/model_interface/surf_couple_extra_mod.F90 @@ -0,0 +1,1422 @@ +MODULE surf_couple_extra_mod +! *****************************COPYRIGHT**************************************** +! (c) Crown copyright, Met Office. All rights reserved. +! +! This routine has been licensed to the other JULES partners for use and +! distribution under the JULES collaboration agreement, subject to the terms and +! conditions set out therein. +! +! [Met Office Ref SC0237] +! *****************************COPYRIGHT**************************************** +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +PRIVATE + +PUBLIC :: surf_couple_extra + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='SURF_COUPLE_EXTRA_MOD' + +CONTAINS + +!=============================================================================== +!Note that the commented intents may not be correct. The declarations section +!will be more correct + +!Routine Layout **PLEASE READ** + +!-Header (please follow ordering guidance) +!-Type correction for timestep +!-SELECT for LSM +!--JULES +!---Redimensioning of arrays (Note 1) +!---Compression to land points (Note 2) +!---Science calls (Note 3) +!---Expansion to full grid +!---UM Diagnostics calls +!--CABLE +!---Details to follow +!-END + +!Note 1- This is for any variable that might have differnt numbers of +! dimensions in standalone vs UM mode. These generally arise due to the +! differing science payloads. + +!Note 2- This is where all gridded variables (with _ij suffix) are compressed +! onto land points. This allows efficient use of OpenMP + +!Note 3- New additions here should ideally consist ONLY of calls to other +! routines. Please encapsulate longer code within a subroutine. + +!=============================================================================== +SUBROUTINE surf_couple_extra( & + ! Arguments used by JULES-standalone + u_1_ij, v_1_ij, & + !Misc INTENT(IN) + a_step, smlt, tile_frac, hcons_soilt, rhostar, & + !Arguments for the UM----------------------------------------- + !IN + land_pts, row_length, rows, river_row_length, river_rows, & + ls_graup_ij, & + cca_2d, nsurft, surft_pts, & + lice_pts, soil_pts, & + stf_sub_surf_roff, & + fexp_soilt, gamtot_soilt, ti_mean_soilt, ti_sig_soilt, & + flash_rate_ancil, pop_den_ancil, wealth_index_ancil, & + a_fsat_soilt, c_fsat_soilt, a_fwet_soilt, c_fwet_soilt, & + ntype, & + delta_lambda, delta_phi, xx_cos_theta_latitude, & + aocpl_row_length, aocpl_p_rows, xpa, xua, xva, ypa, yua, yva, & + g_p_field, g_r_field, n_proc, global_row_length, global_rows, & + global_river_row_length, global_river_rows, flandg, & + trivdir, trivseq, r_area, slope, flowobs1, r_inext, r_jnext, r_land, & + frac_agr_gb, soil_clay_ij,npp_gb, & + u_s_std_surft, & + !INOUT + a_steps_since_riv, & + fsat_soilt, fwetl_soilt, & + zw_soilt, sthzw_soilt, & + ls_rainfrac_gb, & + substore, surfstore, flowin, bflowin, & + tot_surf_runoff_gb, tot_sub_runoff_gb, acc_lake_evap_gb, twatstor, & + asteps_since_triffid, & + inlandout_atm_gb, & + !OUT + dhf_surf_minus_soil, & + land_sea_mask, & + !TYPES containing field data (IN OUT) + crop_vars,psparms,toppdm,fire_vars,ainfo,trif_vars,soilecosse,urban_param, & + progs,trifctltype,coast,jules_vars, & + fluxes, & + lake_vars, & + forcing, & + rivers, & + chemvars, water_resources, & + wtrac_jls, & + work_cbl & + ) + +!Module imports + +!TYPE definitions +USE crop_vars_mod, ONLY: crop_vars_type +USE p_s_parms, ONLY: psparms_type +USE top_pdm, ONLY: top_pdm_type +USE fire_vars_mod, ONLY: fire_vars_type +USE ancil_info, ONLY: ainfo_type +USE trif_vars_mod, ONLY: trif_vars_type +USE soil_ecosse_vars_mod, ONLY: soil_ecosse_vars_type +USE urban_param_mod, ONLY:urban_param_type +USE prognostics, ONLY: progs_type +USE coastal, ONLY: coastal_type +USE trifctl, ONLY: trifctl_type +USE jules_vars_mod, ONLY: jules_vars_type +USE fluxes_mod, ONLY: fluxes_type +USE lake_mod, ONLY: lake_type +USE jules_forcing_mod, ONLY: forcing_type +USE jules_rivers_mod, ONLY: rivers_type +! USE veg3_parm_mod, ONLY: in_dev +! USE veg3_field_mod, ONLY: in_dev +USE jules_chemvars_mod, ONLY: chemvars_type +USE water_resources_vars_mod, ONLY: water_resources_type +USE jules_wtrac_type_mod, ONLY: jls_wtrac_type + +! In general CABLE utilizes a required subset of tbe JULES types, however; +USE work_vars_mod_cbl, ONLY: work_vars_type ! and some kept thru timestep + +!Import interfaces to subroutines called +USE hydrol_mod, ONLY: hydrol +USE snow_mod, ONLY: snow +USE jules_rivers_mod, ONLY: l_rivers, l_inland, rivers_call + +! Code which isn't currently suitable for building into LFRic +#if !defined(LFRIC) +USE surf_couple_rivers_mod, ONLY: surf_couple_rivers +USE inferno_mod, ONLY: calc_soil_carbon_pools +USE inferno_io_mod, ONLY: inferno_io +USE irrigation_mod, ONLY: irrigation_control +USE veg_control_mod, ONLY: veg_control +USE next_gen_biogeochem_mod, ONLY: next_gen_biogeochem +USE sparm_mod, ONLY: sparm +USE infiltration_rate_mod, ONLY: infiltration_rate +USE flake_interface_mod, ONLY: flake_interface + +!ifdef required to manage the different science payloads available for +!coupled and standalone use +#if defined(UM_JULES) +USE diagnostics_riv_mod, ONLY: diagnostics_riv +USE diagnostics_hyd_mod, ONLY: diagnostics_hyd +USE diagnostics_veg_mod, ONLY: diagnostics_veg +#else +USE crop_mod, ONLY: crop +USE fire_timestep_mod, ONLY: fire_timestep +USE metstats_timestep_mod, ONLY: metstats_timestep +USE gridbox_mean_mod, ONLY: soiltiles_to_gbm, surftiles_to_gbm +USE fao_evapotranspiration, ONLY: fao_ref_evapotranspiration +USE water_resources_control_mod, ONLY: water_resources_control +USE zenith_mod, ONLY: photoperiod +! End of standalone Jules code +#endif + +! End of code excluded from LFRic builds +#else +! Code only used for LFRic +USE missing_data_mod, ONLY: rmdi +#endif + +!Variables- The rules here are: +!Each module is only USED once, using an ifdef to control the variables if req. +!Whole modules are grouped into a single ifdef/else +!Alphabetical order + +!Modules common to JULES and UM + +USE ancil_info, ONLY: & +#if !defined(UM_JULES) + dim_cs1, & +#endif + dim_cslayer, nsoilt, nmasst, & + ! Replaces USE statement + dim_soil_n_pool + +USE atm_fields_bounds_mod, ONLY: tdims, tdims_s, pdims_s + +USE conversions_mod, ONLY: rsec_per_day + +! Used for atmospheric deposition +USE jules_deposition_mod, ONLY: l_deposition, l_deposition_from_ukca + +USE deposition_from_surf_couple_extra_mod, & + ONLY: deposition_from_surf_couple_extra + +USE jules_hydrology_mod, ONLY: l_hydrology, l_pdm, l_top, l_var_rainfrac + +USE lake_mod, ONLY: h_snow_sw_att & + , trap_frozen & + , trap_unfrozen + +USE jules_snow_mod, ONLY: nsmax + +USE jules_soil_biogeochem_mod, ONLY: & + soil_model_4pool, soil_bgc_model + +USE jules_soil_mod, ONLY: sm_levels, l_soil_sat_down, confrac + +USE jules_surface_mod, ONLY: l_aggregate, l_flake_model + +USE jules_surface_types_mod, ONLY: npft, ncpft, nnpft, lake + +USE jules_vegetation_mod, ONLY: & +#if !defined(UM_JULES) + l_crop, l_fao_ref_evapotranspiration, & +#endif + l_triffid, l_trif_eq, l_phenol, phenol_period, triffid_period, & + l_inferno, ignition_method, & + ignition_vary_natural, ignition_vary_natural_human, & + l_acclim, n_day_photo_acclim, l_red + +USE jules_irrig_mod, ONLY: l_irrig_dmd + +USE jules_water_tracers_mod, ONLY: n_wtrac_jls + +USE jules_water_resources_mod, ONLY: & + l_water_resources, sw_river_source + +USE jules_water_tracers_mod, ONLY: l_wtrac_jls +USE wtrac_extra_mod, ONLY: wtrac_alloc_extra, wtrac_dealloc_extra, & + wtrac_ex_type +USE wtrac_checks_jls_mod, ONLY: wtrac_checks_sno, wtrac_checks_hyd, & + wtrac_checks_grid +USE wtrac_correct_jls_mod, ONLY: wtrac_correct_sno, wtrac_correct_hyd, & + wtrac_correct_grid + +USE sf_diags_mod, ONLY: sf_diag + +USE timestep_mod, ONLY: timestep + +USE veg3_parm_mod, ONLY: veg3_ctrl,litter_parms,red_parms + +USE veg3_field_mod, ONLY: veg_state,red_state + +USE water_constants_mod, ONLY: rho_water + +USE ereport_mod, ONLY: ereport + +!Modules specific to the UM and not needed in LFRic +#if defined(UM_JULES) && !defined(LFRIC) +USE atm_fields_mod, ONLY: disturb_veg_prev + +USE atm_step_local, ONLY: STASHwork19, STASHwork8, STASHwork26 + +USE model_domain_mod, ONLY: model_type, mt_single_column + +USE stash_array_mod, ONLY: sf +#endif +! Modules used by both UM and LFRic +#if defined(UM_JULES) +USE atm_step_local, ONLY: dim_cs1 +#else +!Modules specific to JULES +USE fire_mod, ONLY: fire_prog, fire_diag, l_fire + +USE metstats_mod, ONLY: metstats_prog, l_metstats + +USE model_grid_mod, ONLY: grid_area_ij + +USE model_time_mod, ONLY: current_time +#endif + +#if !defined(UM_JULES) +USE theta_field_sizes, ONLY: t_i_length, t_j_length +#endif + +!Technical modules +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook +USE ereport_mod, ONLY: ereport +USE jules_print_mgr, ONLY: jules_message, jules_print +USE jules_model_environment_mod, ONLY: lsm_id, jules, cable + + +IMPLICIT NONE + +!----------------------------------------------------------------------------- +! Description: +! Private subroutine for accessing the JULES science routines that are called +! after the implicit code +! +! Code Owner: Please refer to ModuleLeaders.txt +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +!Subroutine Arguments, ordered by intent and type + +INTEGER, INTENT(IN) :: & + land_pts, & + row_length, & + rows, & + river_row_length, & + river_rows, & + nsurft, & + a_step, & + surft_pts(nsurft), & + lice_pts, & + soil_pts, & + ntype, & + aocpl_row_length, & + aocpl_p_rows, & + g_p_field, & + g_r_field, & + n_proc, & + global_row_length, & + global_rows, & + global_river_row_length, & + global_river_rows + +LOGICAL, INTENT(IN) :: & + smlt, & + stf_sub_surf_roff, & + land_sea_mask(row_length, rows) + +REAL(KIND=real_jlslsm), INTENT(IN) :: & + ls_graup_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + u_1_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + v_1_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + tile_frac(land_pts,nsurft), & + rhostar(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface air density + cca_2d(row_length,rows), & + fexp_soilt(land_pts,nsoilt), & + gamtot_soilt(land_pts,nsoilt), & + ti_mean_soilt(land_pts,nsoilt), & + ti_sig_soilt(land_pts,nsoilt), & + npp_gb(land_pts), & + a_fsat_soilt(land_pts,nsoilt), & + c_fsat_soilt(land_pts,nsoilt), & + a_fwet_soilt(land_pts,nsoilt), & + c_fwet_soilt(land_pts,nsoilt), & + frac_agr_gb(land_pts), & + soil_clay_ij(row_length,rows), & + flash_rate_ancil(row_length,rows), & + pop_den_ancil(row_length,rows), & + wealth_index_ancil(row_length,rows), & + u_s_std_surft(land_pts, nsurft), & + + !River routing + delta_lambda, & + delta_phi, & + xx_cos_theta_latitude(tdims_s%i_start:tdims_s%i_end, & + tdims_s%j_start:tdims_s%j_end), & + xpa(aocpl_row_length+1), & + xua(0:aocpl_row_length), & + xva(aocpl_row_length+1), & + ypa(aocpl_p_rows), & + yua(aocpl_p_rows), & + yva(0:aocpl_p_rows), & + flandg(pdims_s%i_start:pdims_s%i_end,pdims_s%j_start:pdims_s%j_end), & + trivdir(river_row_length, river_rows), & + trivseq(river_row_length, river_rows), & + r_area(row_length, rows), & + slope(row_length, rows), & + flowobs1(row_length, rows), & + r_inext(row_length, rows), & + r_jnext(row_length, rows), & + r_land(row_length, rows) + +INTEGER, INTENT(IN OUT) :: & + a_steps_since_riv, & + asteps_since_triffid + +REAL(KIND=real_jlslsm), INTENT(IN OUT) :: & + hcons_soilt(land_pts,nsoilt), & + dhf_surf_minus_soil(land_pts), & + ! Heat flux difference across the FLake snowpack (W/m2) + ls_rainfrac_gb(land_pts), & + fsat_soilt(land_pts,nsoilt), & + fwetl_soilt(land_pts,nsoilt), & + zw_soilt(land_pts,nsoilt), & + sthzw_soilt(land_pts,nsoilt), & + substore(row_length, rows), & + surfstore(row_length, rows), & + flowin(row_length, rows), & + bflowin(row_length, rows), & + acc_lake_evap_gb(row_length,rows), & + tot_surf_runoff_gb(land_pts), & + tot_sub_runoff_gb(land_pts), & + twatstor(river_row_length, river_rows), & + inlandout_atm_gb(land_pts) + +!TYPES containing field data (IN OUT) +TYPE(crop_vars_type), INTENT(IN OUT) :: crop_vars +TYPE(psparms_type), INTENT(IN OUT) :: psparms +TYPE(top_pdm_type), INTENT(IN OUT) :: toppdm +TYPE(fire_vars_type), INTENT(IN OUT) :: fire_vars +TYPE(ainfo_type), INTENT(IN OUT) :: ainfo +TYPE(trif_vars_type), INTENT(IN OUT) :: trif_vars +TYPE(soil_ecosse_vars_type), INTENT(IN OUT) :: soilecosse +TYPE(urban_param_type), INTENT(IN OUT) :: urban_param +TYPE(progs_type), INTENT(IN OUT) :: progs +TYPE(trifctl_type), INTENT(IN OUT) :: trifctltype +TYPE(coastal_type), INTENT(IN OUT) :: coast +TYPE(jules_vars_type), INTENT(IN OUT) :: jules_vars +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes +TYPE(lake_type), INTENT(IN OUT) :: lake_vars +TYPE(forcing_type), INTENT(IN OUT) :: forcing +TYPE(rivers_type), INTENT(IN OUT) :: rivers +!TYPE(in_dev), INTENT(IN OUT) :: veg3_parm +!TYPE(in_dev), INTENT(IN OUT) :: veg3_field +TYPE(chemvars_type), INTENT(IN OUT) :: chemvars +TYPE(water_resources_type), INTENT(IN OUT) :: water_resources +TYPE(jls_wtrac_type), INTENT(IN OUT) :: wtrac_jls + +!CABLE TYPES containing field data (IN OUT) +TYPE(work_vars_type), INTENT(OUT) :: work_cbl + +!============================================================================== +!Local variables + +INTEGER :: & + i,j,k,l,n,m,i_wt, & + !Various counters + errcode + ! error code to pass to ereport. + +#if !defined(UM_JULES) +INTEGER :: crop_call + !indicates whether crop model is to be called +#endif + +INTEGER, PARAMETER :: crop_period = 1 + ! Crop code hard wired to run daily : crop_period = 1 + +REAL(KIND=real_jlslsm) :: & + !Gridbox versions of forcing data (ls_rainfrac_gb is an argument) + con_rainfrac_gb(land_pts), & + ls_rain_gb(land_pts), & + con_rain_gb(land_pts), & + ls_snow_gb(land_pts), & + ls_graup_gb(land_pts), & + con_snow_gb(land_pts), & + pstar_gb(land_pts), & + tl_1_gb(land_pts), & + qw_1_gb(land_pts), & + u_1_gb(land_pts), & + v_1_gb(land_pts), & + + !Passed between different science schemes (various) + qbase_l_soilt(land_pts,nsoilt,sm_levels+1), & + ! Base flow from each soil layer (kg m-2 s-1). + w_flux_soilt(land_pts,nsoilt,0:sm_levels), & + ! Fluxes of water between layers (kg m-2 s-1). + surf_ht_flux_ld(land_pts), & + dwsw_sub_snow, & + !remaining SW flux under the snowpack (W/m2), for FLake with ML-snow. + + !Passed between fire & INFERNO routines only + c_soil_dpm_gb(land_pts), & + ! Gridbox soil C in the Decomposable Plant Material pool (kg m-2). + c_soil_rpm_gb(land_pts), & + ! Gridbox soil C in the Resistant Plant Material pool (kg m-2). +#if !defined(UM_JULES) + smc_gb(land_pts), & + !To allow GBM soil moisture to be passed down to fire. + + !Passed between evapotranspiration routines only + trad(land_pts), & + ! gridbox effective radiative temperature (assuming emissivity=1) +#endif + + !Passed from snow or hydrol to diagnostics_hyd only + snow_mass_gb(land_pts), & + + !Passed between river routing and diagnostics_riv only + riverout(row_length, rows), & + riverout_rgrid(river_row_length, river_rows), & + box_outflow(river_row_length, river_rows), & + box_inflow(river_row_length, river_rows), & + inlandout_riv(river_row_length,river_rows) + +! Smoothing factor of exponential filter used with temperature acclimation. +REAL(KIND=real_jlslsm) :: alpha_acclim + +! Water tracer local arrays +TYPE(wtrac_ex_type) :: wtrac_ex + +! Compressed water tracer rain and snow fields on land points only +REAL(KIND=real_jlslsm), ALLOCATABLE :: & + con_snow_wtrac(:,:), & + con_rain_wtrac(:,:), & + ls_rain_wtrac(:,:), & + ls_snow_wtrac(:,:), & + ls_graup_wtrac(:,:) + +INTEGER :: errorstatus +CHARACTER(LEN=*), PARAMETER :: RoutineName = 'SURF_COUPLE_EXTRA' + +!Dr Hook variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +!=============================================================================== +!End of header +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +!CABLE_LSM: implement switching based on lsm_id +SELECT CASE( lsm_id ) + + !============================================================================= +CASE ( jules ) + + !============================================================================= + ! Redimensioning of arrays + +#if defined(UM_JULES) + !Dimensionality of variables differ betweeen UM and standalone, so copy across + + + ! Change 2d to 1d soil clay content for soil respiration- dimensionalities + ! differ between UM and standalone. + ! Soil tiling not currently in the UM, so broadcast ij value to all tiles. + ! Multi-layer clay not currently in UM so set all layers to same value. + IF ( soil_bgc_model == soil_model_4pool ) THEN + m = 1 + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + DO n = 1, dim_cslayer + psparms%clay_soilt(l,m,n) = soil_clay_ij(i,j) + END DO + END DO + END IF + + ! Compress pop_den and flash rate fields to land points if required. + ! In the UM the ancils are 2D fields. + IF (l_inferno) THEN + IF (ignition_method == ignition_vary_natural) THEN + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + fire_vars%flash_rate(l) = flash_rate_ancil(i,j) + END DO + END IF + IF (ignition_method == ignition_vary_natural_human) THEN + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + fire_vars%flash_rate(l) = flash_rate_ancil(i,j) + fire_vars%pop_den(l) = pop_den_ancil(i,j) + fire_vars%wealth_index(l) = wealth_index_ancil(i,j) + END DO + END IF + END IF +#endif + + ! Set up water tracer working arrays + CALL wtrac_alloc_extra(land_pts, nsurft, n_wtrac_jls, wtrac_ex) + + !Encompassing IF statement for land_pts > 0. We exit/re-enter this IF about + !half way down to allow for river routing + + IF (l_hydrology .AND. land_pts > 0 ) THEN + + !=========================================================================== + ! Compression to land points + + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + ls_rain_gb(l) = forcing%ls_rain_ij(i,j) + con_rain_gb(l) = forcing%con_rain_ij(i,j) + con_snow_gb(l) = forcing%con_snow_ij(i,j) + ls_snow_gb(l) = forcing%ls_snow_ij(i,j) + ls_graup_gb(l) = ls_graup_ij(i,j) + pstar_gb(l) = forcing%pstar_ij(i,j) + tl_1_gb(l) = forcing%tl_1_ij(i,j) + qw_1_gb(l) = forcing%qw_1_ij(i,j) + u_1_gb(l) = u_1_ij(i,j) + v_1_gb(l) = v_1_ij(i,j) + END DO + + IF (l_wtrac_jls) THEN + ! Repeat for water tracers + ALLOCATE(con_snow_wtrac(land_pts,n_wtrac_jls)) + ALLOCATE(con_rain_wtrac(land_pts,n_wtrac_jls)) + ALLOCATE(ls_snow_wtrac(land_pts,n_wtrac_jls)) + ALLOCATE(ls_rain_wtrac(land_pts,n_wtrac_jls)) + ALLOCATE(ls_graup_wtrac(land_pts,n_wtrac_jls)) + DO i_wt = 1, n_wtrac_jls + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + con_snow_wtrac(l,i_wt) = wtrac_jls%con_snow_ij(i,j,i_wt) + con_rain_wtrac(l,i_wt) = wtrac_jls%con_rain_ij(i,j,i_wt) + ls_rain_wtrac(l,i_wt) = wtrac_jls%ls_rain_ij(i,j,i_wt) + ls_snow_wtrac(l,i_wt) = wtrac_jls%ls_snow_ij(i,j,i_wt) + ls_graup_wtrac(l,i_wt) = 0.0 + END DO + END DO + ELSE + ALLOCATE(con_snow_wtrac(1,1)) + ALLOCATE(con_rain_wtrac(1,1)) + ALLOCATE(ls_snow_wtrac(1,1)) + ALLOCATE(ls_rain_wtrac(1,1)) + ALLOCATE(ls_graup_wtrac(1,1)) + END IF + + !Pass jules the modelled rain fractions + !In standalone mode, both rainfracs are passed in as zeroed arrays + IF (l_var_rainfrac) THEN + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + + con_rainfrac_gb(l) = MIN(cca_2d(i,j),0.5) !As in diagnostics_conv + + !provide some safety checking for convective rain with no CCA + IF (con_rain_gb(l) > 0.0 ) THEN + IF (con_rainfrac_gb(l) == 0.0) THEN + con_rainfrac_gb(l) = confrac + !and for very small CCA amounts + ELSE IF (con_rainfrac_gb(l) < 0.01) THEN + con_rainfrac_gb(l) = 0.01 + END IF + END IF + + !provide some safety checking for ls rain with no rainfrac + IF (ls_rain_gb(l) > 0.0) THEN + IF (ls_rainfrac_gb(l) == 0.0) THEN + ls_rainfrac_gb(l) = 0.5 + !and for very small rainfrac amounts + ELSE IF (ls_rainfrac_gb(l) < 0.01) THEN + ls_rainfrac_gb(l) = 0.01 + END IF + END IF + END DO + + ELSE !use original default values + DO l = 1, land_pts + con_rainfrac_gb(l) = confrac + ls_rainfrac_gb(l) = 1.0 + END DO + END IF !l_var_rainfrac + + !=========================================================================== + ! Science calls + + !--------------------------------------------------------------------------- + ! Calculation of atmospheric deposition parameters and fluxes + !--------------------------------------------------------------------------- + + ! For JULES standalone applications, the JULES-based atmospheric deposition + ! routines are called here, by setting the switches l_deposition to true + ! and l_deposition_from_ukca to false. + + ! For coupled UM_JULES applications, the JULES-based atmospheric deposition + ! routines can only be called from the UKCA routine ukca_chemistry_ctl + ! (or its equivalents, ukca_chemistry_ctl_BE and ukca_chemistry_ctl_col), + ! setting the switch l_deposition_from_ukca to true. This is to demonstrate + ! that the JULES-based routines give outputs that bit compare with + ! UKCA known good outputs. + + ! The longer-term plan is for the UKCA to provide surface atmospheric + ! concentrations of the deposited species with JULES returning deposition + ! fluxes, with the call routed through here. + + IF (l_deposition .AND. (.NOT. l_deposition_from_ukca)) THEN + + CALL deposition_from_surf_couple_extra( & + a_step, timestep, & + ! passed to calling routine, defined in module ainfo + row_length, rows, land_pts, ainfo%land_index, & + ! From TYPE ainfo, coast + surft_pts, ainfo%surft_index, ainfo%ice_fract_ij, coast%fland, & + ! From TYPE forcing: pressure and surface water mixing ratio + forcing%pstar_ij, forcing%qw_1_ij, & + ! From TYPE fluxes: temperature, SW radiation & surface roughness + fluxes%tstar_ij, fluxes%sw_surft, fluxes%z0h_surft, & + ! From TYPE progs + progs%gc_surft, progs%lai_pft, progs%canht_pft, progs%tstar_surft, & + progs%canopy_surft, progs%smc_soilt, & + ! From TYPES jules_vars + jules_vars%zh, jules_vars%dzl, & + ! passed From calling routine + tile_frac, & + ! From TYPE chemvars: atmospheric variables & diagnostics + chemvars, & + ! Parameters for Paulot et al H2 deposition scheme + nsurft, nsoilt, sm_levels, dim_cslayer, dim_cs1, & + psparms%smvcst_soilt, progs%snowdepth_surft, progs%cs_pool_soilt, & + progs%t_soil_soilt ) + + END IF + + !Snow (standalone and UM) + CALL snow ( land_pts,timestep,smlt,nsurft,n_wtrac_jls,surft_pts, & + ainfo%surft_index,psparms%catch_snow_surft,con_snow_gb, & + con_rain_gb, tile_frac,ls_snow_gb,ls_graup_gb,ls_rain_gb, & + fluxes%ei_surft,psparms%hcap_soilt(:,:,1),hcons_soilt, & + fluxes%melt_surft,fluxes%snowinc_surft, & + progs%smcl_soilt(:,:,1),psparms%sthf_soilt(:,:,1), & + fluxes%surf_htf_surft, progs%t_soil_soilt(:,:,1), & + progs%tsurf_elev_surft, progs%tstar_surft, psparms%z0_surft, & + psparms%smvcst_soilt(:,:,1), con_snow_wtrac, & + wtrac_jls%ei_surft, fluxes%fsnow, progs%rgrain_surft, & + progs%rgrainl_surft, progs%rho_snow_grnd_surft, & + progs%sice_surft,progs%sliq_surft,progs%snow_grnd_surft, & + progs%snow_surft, progs%snowdepth_surft, progs%tsnow_surft, & + progs%nsnow_surft, con_rain_wtrac, ls_rain_wtrac, & + ls_snow_wtrac, ls_graup_wtrac, wtrac_ex%melt_surft, & + wtrac_jls%snow_grnd_surft, wtrac_jls%snow_surft, & + wtrac_jls%sice_surft, wtrac_jls%sliq_surft, progs%ds_surft, & + fluxes%hf_snow_melt_gb,snow_mass_gb,progs%rho_snow_surft, & + fluxes%snomlt_sub_htf_gb, fluxes%snow_melt_gb, & + fluxes%snow_soil_htf,surf_ht_flux_ld, & + sf_diag, dhf_surf_minus_soil, wtrac_ex%snow_melt, & + wtrac_ex%lake_snow_melt, & + ! New Arguments to replace USE statements + ! jules_internal + jules_vars%unload_backgrnd_pft, npft, l_flake_model, & + !Ancil info (IN) + ainfo%l_lice_point, ainfo%l_lice_surft, & + ! Types Variables + lake_vars%lake_h_ice_gb, lake_vars%hcon_lake, & + lake_vars%ts1_lake_gb, lake_vars%lake_snow_melt, & + lake_vars%non_lake_frac, lake_vars%lake_h_mxl_gb, & + lake_vars%lake_depth_gb) + + IF (l_wtrac_jls) THEN + + ! Check that the snow normal water tracer stores are not significantly + ! diverging from the water fields + CALL wtrac_checks_sno(land_pts, nsurft, nsmax, surft_pts, & + ainfo%surft_index, & + progs%snow_surft, progs%snow_grnd_surft, & + progs%sice_surft,progs%sliq_surft, & + wtrac_jls%snow_surft(:,:,1), & + wtrac_jls%snow_grnd_surft(:,:,1), & + wtrac_jls%sice_surft(:,:,:,1), & + wtrac_jls%sliq_surft(:,:,:,1)) + + ! Correct snow water tracer stores to prevent build up of numerical error + CALL wtrac_correct_sno(land_pts, nsurft, nsmax, n_wtrac_jls, & + surft_pts, ainfo%surft_index, & + progs%snow_surft, progs%snow_grnd_surft, & + progs%sice_surft,progs%sliq_surft, & + wtrac_jls%snow_surft, wtrac_jls%snow_grnd_surft, & + wtrac_jls%sice_surft, wtrac_jls%sliq_surft) + + END IF ! l_wtrac_jls + + !Hydrology (standalone and UM) + CALL hydrol ( & + land_pts, soil_pts, lice_pts, sm_levels, npft, nsurft, dim_cs1, & + asteps_since_triffid, n_wtrac_jls, timestep, ainfo%l_lice_surft, & + l_inland, l_pdm, l_soil_sat_down, l_top, stf_sub_surf_roff, & + ainfo%soil_index, ainfo%lice_index, & + surft_pts, ainfo%surft_index, progs%nsnow_surft, & + psparms%bexp_soilt, psparms%hcap_soilt, psparms%hcon_soilt, & + psparms%satcon_soilt, psparms%sathh_soilt, & + psparms%smvcst_soilt, psparms%smvcwt_soilt, & + psparms%catch_surft, crop_vars%frac_irr_soilt, psparms%infil_surft, & + lake_vars%non_lake_frac, toppdm%slope_gb, tile_frac, & + con_rainfrac_gb, con_rain_gb, ls_rainfrac_gb, ls_rain_gb, & + con_rain_wtrac, ls_rain_wtrac, & + surf_ht_flux_ld, & + fluxes%ecan_surft, crop_vars%ext_irr_soilt, fluxes%ext_soilt, & + wtrac_jls%ecan_surft, wtrac_jls%ext_soilt, & + progs%snowdepth_surft, fluxes%melt_surft, fluxes%snow_melt_gb, & + wtrac_ex%melt_surft, wtrac_ex%snow_melt, fluxes%snow_soil_htf, & + a_fsat_soilt,c_fsat_soilt,a_fwet_soilt, c_fwet_soilt, & + fexp_soilt, ti_mean_soilt, & + npp_gb, inlandout_atm_gb, wtrac_jls%inlandout_atm_gb, & + progs%canopy_surft, wtrac_jls%canopy_surft, & + progs%smcl_soilt, psparms%sthf_soilt, & + psparms%sthu_soilt, crop_vars%sthu_irr_soilt, progs%tsoil_deep_gb, & + progs%t_soil_soilt, progs%t_soil_soilt_acc, progs%tsurf_elev_surft, & + wtrac_jls%smcl_soilt, wtrac_jls%sthf_soilt, wtrac_jls%sthu_soilt, & + fsat_soilt, fwetl_soilt, sthzw_soilt, zw_soilt, wtrac_jls%sthzw_soilt, & + progs%cs_pool_soilt, trifctltype%resp_s_soilt, & + toppdm%cs_ch4_soilt, toppdm%fch4_wetl_acc_soilt, & + progs%substr_ch4, progs%mic_ch4, progs%mic_act_ch4, progs%acclim_ch4, & + progs%n_inorg_avail_pft, progs%n_inorg_soilt_lyrs, & + trif_vars%n_leach_gb_acc, & + progs%canopy_gb, progs%smc_soilt, & + toppdm%drain_soilt, toppdm%dun_roff_soilt, fluxes%sub_surf_roff_gb, & + fluxes%surf_roff_gb, fluxes%tot_tfall_gb, fluxes%tot_tfall_surft, & + w_flux_soilt, wtrac_jls%surf_roff_gb, wtrac_jls%sub_surf_roff_gb, & + toppdm%qbase_soilt, qbase_l_soilt, toppdm%qbase_zw_soilt, & + toppdm%fch4_wetl_soilt, toppdm%fch4_wetl_cs_soilt, & + toppdm%fch4_wetl_npp_soilt, toppdm%fch4_wetl_resps_soilt, & + trif_vars%n_leach_soilt) + + + DEALLOCATE(ls_graup_wtrac) + DEALLOCATE(ls_rain_wtrac) + DEALLOCATE(ls_snow_wtrac) + DEALLOCATE(con_rain_wtrac) + DEALLOCATE(con_snow_wtrac) + + IF (l_wtrac_jls) THEN + + ! Check that the hydrological normal water tracer stores are not + ! significantly diverging from the water fields + CALL wtrac_checks_hyd(land_pts, nsurft, nsoilt, sm_levels, soil_pts, & + ainfo%soil_index, surft_pts, ainfo%surft_index, & + progs%canopy_surft, psparms%sthu_soilt, & + psparms%sthf_soilt, progs%smcl_soilt, sthzw_soilt, & + wtrac_jls%canopy_surft(:,:,1), & + wtrac_jls%sthu_soilt(:,:,:,1), & + wtrac_jls%sthf_soilt(:,:,:,1), & + wtrac_jls%smcl_soilt(:,:,:,1), & + wtrac_jls%sthzw_soilt(:,:,1)) + + ! Correct hydrological water tracer stores to prevent build up of + ! numerical errors + CALL wtrac_correct_hyd(land_pts, nsurft, nsoilt, sm_levels, n_wtrac_jls, & + soil_pts, ainfo%soil_index, & + surft_pts, ainfo%surft_index, & + progs%canopy_surft, psparms%sthu_soilt, & + psparms%sthf_soilt, progs%smcl_soilt, sthzw_soilt, & + wtrac_jls%canopy_surft, wtrac_jls%sthu_soilt, & + wtrac_jls%sthf_soilt, wtrac_jls%smcl_soilt, & + wtrac_jls%sthzw_soilt) + + END IF ! l_wtrac_jls + + END IF ! ( l_hydrology .AND. land_pts /= 0 ) + + ! Code not yet ported to LFRic +#if !defined(LFRIC) + + !Here we need to exit the land_pts IF to allow river routing to be called + !on all MPI ranks. This is because it is possible that a rank with + !land_pts = 0 still has a river routing point. + +#if !defined(UM_JULES) + ! Water resources (standalone; not yet allowed in UM). + IF ( l_water_resources ) THEN + CALL water_resources_control( & + rivers%global_land_index, ainfo%land_index, & + rivers%map_river_to_land_points, rivers%rivers_index_rp, & + forcing%con_rain_ij, forcing%con_snow_ij, & + water_resources%conv_loss_frac, & + water_resources%demand_rate_domestic, & + water_resources%demand_rate_industry, & + water_resources%demand_rate_livestock, & + water_resources%demand_rate_transfers, crop_vars%dvi_cpft, & + flandg, crop_vars%frac_irr_soilt, ainfo%frac_soilt, & + ainfo%frac_surft, grid_area_ij, & + forcing%ls_rain_ij, forcing%ls_snow_ij, forcing%lw_down_ij, & + psparms%smvccl_soilt, & + psparms%smvcst_soilt, psparms%smvcwt_soilt, psparms%sthf_soilt, & + fluxes%sw_surft, forcing%tl_1_ij, progs%tstar_surft, & + crop_vars%icntmax_gb, crop_vars%plant_n_gb, & + water_resources%demand_accum, & + crop_vars%prec_1_day_av_gb, crop_vars%prec_1_day_av_use_gb, & + rivers%rfm_surfstore_rp, rivers%rivers_sto_rp, & + crop_vars%rn_1_day_av_gb, & + crop_vars%rn_1_day_av_use_gb, water_resources%sfc_water_frac, & + progs%smcl_soilt, & + crop_vars%sthu_irr_soilt, psparms%sthu_soilt, & + sthzw_soilt, fluxes%sub_surf_roff_gb, crop_vars%tl_1_day_av_gb, & + crop_vars%tl_1_day_av_use_gb, water_resources%priority_order, & + water_resources%demand_unmet, water_resources%gw_abstracted, & + water_resources%gw_avail, water_resources%gw_nr_abstracted, & + crop_vars%irrig_water_gb, water_resources%net_abstracted_river, & + water_resources%sw_abstracted, water_resources%sw_avail_total, & + water_resources%water_removed ) + END IF +#endif + + IF ( l_rivers ) THEN + CALL surf_couple_rivers( & + !INTEGER, INTENT(IN) + land_pts, n_wtrac_jls, & + !REAL, INTENT(IN) + water_resources%net_abstracted_river, & + fluxes%sub_surf_roff_gb, fluxes%surf_roff_gb, & + wtrac_jls%sub_surf_roff_gb, wtrac_jls%surf_roff_gb, & + !INTEGER, INTENT(INOUT) + a_steps_since_riv, & + !REAL, INTENT (INOUT) + tot_surf_runoff_gb, tot_sub_runoff_gb, acc_lake_evap_gb, & + wtrac_jls%tot_surf_runoff_gb, wtrac_jls%tot_sub_runoff_gb, & + wtrac_jls%acc_lake_evap_gb, & + !REAL, INTENT (OUT) + rivers%rivers_sto_per_m2_on_landpts, fluxes%rflow_gb, fluxes%rrun_gb, & + !Arguments for the UM----------------------------------------- + !INTEGER, INTENT(IN) + n_proc, row_length, rows, river_row_length, river_rows, & + ainfo%land_index, aocpl_row_length, aocpl_p_rows, g_p_field, & + g_r_field, global_row_length, global_rows, global_river_row_length, & + global_river_rows, & + !REAL, INTENT(IN) + fluxes%lake_evap, delta_lambda, delta_phi, xx_cos_theta_latitude, & + xpa, xua, xva, ypa, yua, yva, flandg, trivdir, & + trivseq, r_area, slope, flowobs1, r_inext, r_jnext, r_land, & + psparms%smvcst_soilt, psparms%smvcwt_soilt, ainfo%frac_surft, & + wtrac_jls%lake_evap, & + !REAL, INTENT(INOUT) + substore, surfstore, flowin, bflowin, twatstor, & + progs%smcl_soilt, psparms%sthu_soilt, & + wtrac_jls%twatstor, wtrac_jls%smcl_soilt, wtrac_jls%sthu_soilt, & + !REAL, INTENT(OUT) + inlandout_atm_gb, inlandout_riv, riverout, box_outflow, box_inflow, & + riverout_rgrid, wtrac_jls%inlandout_atm_gb, & + ! imported rivers arrays + rivers) + END IF ! l_rivers (ATMOS) + + IF (l_rivers .AND. l_wtrac_jls) THEN + + ! Check that the river normal water tracer store is not diverging + ! signifantly from the water field. + CALL wtrac_checks_grid(river_row_length, river_rows, & + twatstor, wtrac_jls%twatstor(:,:,1),'twatstor') + + ! Correct river water tracer store to prevent build up of numerical errors + CALL wtrac_correct_grid(river_row_length, river_rows, n_wtrac_jls, & + twatstor, wtrac_jls%twatstor) + + ! Check that the hydrological normal water tracer stores are not + ! significantly diverging from the water fields + CALL wtrac_checks_hyd(land_pts, nsurft, nsoilt, sm_levels, soil_pts, & + ainfo%soil_index, surft_pts, ainfo%surft_index, & + progs%canopy_surft, psparms%sthu_soilt, & + psparms%sthf_soilt, progs%smcl_soilt, sthzw_soilt, & + wtrac_jls%canopy_surft(:,:,1), & + wtrac_jls%sthu_soilt(:,:,:,1), & + wtrac_jls%sthf_soilt(:,:,:,1), & + wtrac_jls%smcl_soilt(:,:,:,1), & + wtrac_jls%sthzw_soilt(:,:,1)) + + END IF ! (l_rivers .AND. l_wtrac_jls) + + !Irrigation (standalone and UM for some options). + IF ( l_irrig_dmd .AND. .NOT. l_water_resources ) THEN + CALL irrigation_control( a_step, land_pts, & + ainfo%land_index, & + forcing%con_rain_ij, forcing%con_snow_ij, & + crop_vars%dvi_cpft, & + crop_vars%frac_irr_soilt, ainfo%frac_surft, & + forcing%ls_rain_ij, forcing%ls_snow_ij, & + forcing%lw_down_ij, & + psparms%smvccl_soilt, psparms%smvcst_soilt, & + psparms%smvcwt_soilt, psparms%sthf_soilt, & + fluxes%sw_surft, forcing%tl_1_ij, progs%tstar_surft, & + crop_vars%icntmax_gb, crop_vars%plant_n_gb, & + crop_vars%irrDaysDiag_gb, crop_vars%prec_1_day_av_gb, & + crop_vars%prec_1_day_av_use_gb, & + crop_vars%rn_1_day_av_gb, crop_vars%rn_1_day_av_use_gb, & + crop_vars%tl_1_day_av_gb, crop_vars%tl_1_day_av_use_gb, & + progs%smcl_soilt, crop_vars%sthu_irr_soilt, & + psparms%sthu_soilt, sthzw_soilt, & + crop_vars%irrig_water_gb, & + !New arguments replacing USE statements + !jules_rivers_mod + rivers%rivers_sto_per_m2_on_landpts, & + rivers%rivers_adj_on_landpts, & + rivers ) + ELSE IF ( .NOT. l_water_resources ) THEN + ! Set sthu_irr_soilt to 0.0 in case it is still reported + ! It would be better to not allocate this array when it is not being used. + ! This led to a memory leak in the D1 array (UM). + crop_vars%sthu_irr_soilt(:,:,:) = 0.0 + END IF ! l_irrig_dmd + + !Restart the IF for land_pts > 0 now that river routing is done with + IF (land_pts > 0) THEN + + !Crops (standalone only) +#if !defined(UM_JULES) + IF ( l_crop ) THEN + crop_call = MOD ( REAL(a_step), & + REAL(crop_period) * rsec_per_day / timestep ) + + DO n = 1,ncpft + DO l = 1,land_pts + trifctltype%npp_acc_pft(l,nnpft + n) = & + trifctltype%npp_acc_pft(l,nnpft + n) & + + (trifctltype%npp_pft(l,nnpft + n) * timestep) + END DO + END DO + + CALL photoperiod(t_i_length * t_j_length, crop_vars%phot, & + crop_vars%dphotdt) + + CALL crop(t_i_length * t_j_length, land_pts, ainfo%land_index, a_step, & + crop_call, sm_levels, ainfo%frac_surft, crop_vars%phot, & + crop_vars%dphotdt, & + sf_diag%t1p5m_surft, progs%t_soil_soilt, psparms%sthu_soilt, & + psparms%smvccl_soilt, & + psparms%smvcst_soilt, trifctltype%npp_acc_pft, & + progs%canht_pft, progs%lai_pft, crop_vars%dvi_cpft, & + crop_vars%rootc_cpft, & + crop_vars%harvc_cpft, & + crop_vars%reservec_cpft, crop_vars%croplai_cpft, & + crop_vars%cropcanht_cpft, & + psparms%catch_surft, psparms%z0_surft, & + !New arguments replacing USE statements + !crop_vars_mod + crop_vars%sow_date_cpft, crop_vars%tt_veg_cpft, & + crop_vars%tt_rep_cpft, & + crop_vars%latestharv_date_cpft, & + crop_vars%yield_diag_cpft, crop_vars%stemc_diag_cpft, & + crop_vars%leafc_diag_cpft, crop_vars%nonyield_diag_cpft, & + crop_vars%harvest_trigger_cpft, crop_vars%harvest_counter_cpft,& + !ancil_info (IN) + ainfo%l_lice_point, ainfo%l_lice_surft) + + END IF ! l_crop + + !Metstats (standalone only) + !Beware- does not currently account for graupel + IF ( l_metstats ) THEN + CALL metstats_timestep(tl_1_gb, qw_1_gb, u_1_gb, v_1_gb, ls_rain_gb, & + con_rain_gb, ls_snow_gb, con_snow_gb, pstar_gb, & + metstats_prog, & + !Vars that should be USED but can't due to + !UM/standalone differences + current_time%TIME, timestep,land_pts) + END IF + + !Fire (standalone only) + IF ( l_fire ) THEN + !Calculate the gridbox mean soil moisture + smc_gb = soiltiles_to_gbm(progs%smc_soilt, ainfo) + CALL fire_timestep(metstats_prog, smc_gb, fire_prog, fire_diag, & + !Vars that should be USED but can't due to + !UM/standalone differences + current_time%TIME, current_time%month, timestep, & + land_pts) + END IF +#endif + + !INFERNO (standalone and UM) + IF ( l_inferno ) THEN + CALL calc_soil_carbon_pools(land_pts, soil_pts, ainfo%soil_index, & + dim_cs1, progs%cs_pool_soilt, & + c_soil_dpm_gb, c_soil_rpm_gb) + + CALL inferno_io( sf_diag%t1p5m_surft, sf_diag%q1p5m_surft, pstar_gb, & + psparms%sthu_soilt, sm_levels, & + ainfo%frac_surft, c_soil_dpm_gb, c_soil_rpm_gb, & + progs%canht_pft, & + ls_rain_gb, con_rain_gb, & + fire_vars, & + land_pts, ignition_method, & + nsurft, asteps_since_triffid, & + ! New Arguments to replace USE statements + ! TRIF_VARS_MOD + trif_vars%g_burn_pft_acc) + END IF ! l_inferno + + IF (l_phenol .OR. l_triffid) THEN + IF (l_red) THEN + + !Update veg3 structures + !Structures can then be passed straight through + ! Copy prognostics and forcings in and then out + ! Same calling when coulpled to atmos + DO n = 1, nnpft + DO k = 1, land_pts + l = ainfo%land_index(k) + veg_state%npp_acc(l,n) = trifctltype%npp_acc_pft(l,n) + END DO + END DO + + CALL next_gen_biogeochem( & + !IN control vars + asteps_since_triffid,land_pts,nnpft,nmasst,veg3_ctrl, & + ainfo, & + !IN parms + litter_parms,red_parms, & + !INOUT data structures + veg_state,red_state & + !OUT diagnostics + ) + + ! Update the physical state of the land + ! Outside of main call as + ! - Uses non-veg tilepts + ! - Accesses fields via USE statements + + ! Note use lai_bal for now as we don't yet have phenology + CALL sparm (land_pts,nsurft,surft_pts,ainfo%surft_index, & + veg_state%frac,veg_state%canht, & + veg_state%lai_bal,psparms%z0m_soil_gb, & + psparms%catch_snow_surft,psparms%catch_surft, & + psparms%z0_surft,psparms%z0h_bare_surft,urban_param%ztm_gb) + + CALL infiltration_rate(land_pts,nsurft,surft_pts,ainfo%surft_index, & + psparms%satcon_soilt,veg_state%frac, & + psparms%infil_surft) + + DO n = 1, nnpft + DO k = 1, land_pts + l = ainfo%land_index(k) + trifctltype%npp_acc_pft(l,n) = veg_state%npp_acc(l,n) + END DO + END DO + + ELSE ! Triffid based model + + !Vegetation (standalone and UM- differing functionality) + + CALL veg_control( & + land_pts, nsurft, dim_cs1, & + a_step, asteps_since_triffid, & + phenol_period, triffid_period, & + l_phenol, l_triffid, l_trif_eq, & + timestep, frac_agr_gb, trif_vars%frac_past_gb, & + trif_vars%frac_biocrop_gb, psparms%satcon_soilt, & + trifctltype%g_leaf_acc_pft, trifctltype%g_leaf_phen_acc_pft, & + trifctltype%npp_acc_pft, & + trifctltype%resp_w_acc_pft, & + ainfo%frac_surft, progs%lai_pft, & + psparms%clay_soilt, psparms%z0m_soil_gb, & + progs%canht_pft, & + psparms%catch_snow_surft, psparms%catch_surft, & + psparms%infil_surft, psparms%z0_surft, psparms%z0h_bare_surft, & + trifctltype%c_veg_pft, trifctltype%cv_gb, trifctltype%lit_c_pft, & + trifctltype%lit_c_mn_gb, trifctltype%g_leaf_day_pft, & + trifctltype%g_leaf_phen_pft,trifctltype%lai_phen_pft, & + trifctltype%g_leaf_dr_out_pft, trifctltype%npp_dr_out_pft, & + trifctltype%resp_w_dr_out_pft, qbase_l_soilt, & + psparms%sthf_soilt, psparms%sthu_soilt, & + w_flux_soilt, progs%t_soil_soilt,progs%cs_pool_soilt, & + progs%frac_c_label_pool_soilt, & + !New arguments replacing USE statements + !trif_vars_mod (IN OUT) + trif_vars, & + !crop_vars_mod (IN) + crop_vars%rootc_cpft, crop_vars%harvc_cpft, crop_vars%reservec_cpft, & + crop_vars%stemc_diag_cpft, & + crop_vars%leafc_diag_cpft, crop_vars%dvi_cpft, & + ! prognostics (IN) + progs%wood_prod_fast_gb, progs%wood_prod_med_gb, & + progs%wood_prod_slow_gb, progs%frac_agr_prev_gb, & + progs%frac_past_prev_gb, progs%frac_biocrop_prev_gb, & + progs%n_inorg_gb, progs%n_inorg_soilt_lyrs, & + progs%n_inorg_avail_pft, progs%ns_pool_gb, & + progs%triffid_co2_gb, progs%t_soil_soilt_acc, & + progs%years_since_harvest, & + ! p_s_parms (IN) + psparms%bexp_soilt, psparms%sathh_soilt, & + psparms%smvcst_soilt, psparms%smvcwt_soilt, & + psparms%soil_ph_soilt, & + ! soil_ecosse_vars_mod (OUT) + soilecosse%n_soil_pool_soilt, dim_soil_n_pool, & + soilecosse%co2_soil_gb, soilecosse%n2o_soil_gb, & + soilecosse%n2o_denitrif_gb, soilecosse%n2o_nitrif_gb, & + soilecosse%n2o_partial_nitrif_gb, soilecosse%n2_denitrif_gb, & + soilecosse%n_denitrification_gb, soilecosse%n_leach_amm_gb, & + soilecosse%n_leach_nit_gb, soilecosse%n_nitrification_gb, & + soilecosse%no_soil_gb, soilecosse%soil_c_add, soilecosse%soil_n_add, & + ! soil_ecosse_vars_mod (IN) + soilecosse%qbase_l_driver, soilecosse%sthf_driver, & + soilecosse%sthu_driver, soilecosse%tsoil_driver, & + soilecosse%wflux_driver, & + !ancil_info (IN) + ainfo%l_lice_point, ainfo%l_lice_surft, ainfo%soil_index, & + !TYPES + soilecosse, urban_param,trifctltype) + + END IF + END IF + + IF ( l_acclim ) THEN + ! Update the growth temperature for the thermal acclimation of + ! photosynthesis. + alpha_acclim = 1.0 - EXP( -1.0 * timestep & + / ( n_day_photo_acclim * rsec_per_day ) ) + progs%t_growth_gb(:) = progs%t_growth_gb(:) & + + alpha_acclim * ( tl_1_gb(:) - progs%t_growth_gb(:) ) + END IF + +#if !defined(UM_JULES) + !Reference evapotranspiration (standalone only) + IF (l_fao_ref_evapotranspiration) THEN + trad = ( surftiles_to_gbm(progs%tstar_surft**4, ainfo) )**0.25 + CALL fao_ref_evapotranspiration(soil_pts, ainfo%soil_index, & + land_pts, ainfo%land_index, sf_diag%t1p5m, & + forcing%sw_down_ij, forcing%lw_down_ij, fluxes%surf_ht_flux_ij, & + sf_diag%u10m, sf_diag%v10m, sf_diag%q1p5m, forcing%pstar_ij, & + trad, trif_vars%fao_et0) + END IF + ! End of standalone Jules code +#endif + + !----------------------------------------------------------------------- + ! call to the FLake interface + !----------------------------------------------------------------------- + IF ( l_flake_model & + .AND. ( .NOT. l_aggregate) & + .AND. (land_pts > 0 ) ) THEN + + DO k = 1,surft_pts(lake) + l = ainfo%surft_index(k,lake) + j=(ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + + ! U* of lake obtained by matching surface stress + lake_vars%u_s_lake_gb(l) = u_s_std_surft(l, lake) & + * SQRT( rhostar(i,j) / rho_water ) + + ! Downwelling SW on lake tile + lake_vars%sw_down_gb(l) = fluxes%sw_surft(l,lake) / & + (1.0 - lake_vars%lake_albedo_gb(l)) + + ! Take the net SW flux out of the surface heat flux + ! since this is done separately within FLake. + lake_vars%surf_ht_flux_lk_gb(l) = lake_vars%surf_ht_flux_lake_ij(i,j) -& + fluxes%sw_surft(l,lake) + + IF ( (nsmax > 0) .AND. & + (lake_vars%lake_h_snow_gb(l) > EPSILON(1.0)) ) THEN + ! For the new snow scheme, FLake is forced with zero snow + ! and the forcing fluxes are those at the bottom of the snowpack. + ! The order of the following calculations is important. + + ! attenuated DWSW below snow + dwsw_sub_snow = lake_vars%sw_down_gb(l) * & + EXP( -lake_vars%lake_h_snow_gb(l) / h_snow_sw_att ) + + ! following last calculation, set snow depth to zero + lake_vars%lake_h_snow_gb(l) = 0.0 + + ! heat flux into the lake becomes + ! the sub-snow value minus remaining DWSW + + lake_vars%surf_ht_flux_lk_gb(l) = lake_vars%surf_ht_flux_lake_ij(i,j)& + - dhf_surf_minus_soil(l) & + - dwsw_sub_snow + + ! now overwrite the DWSW with the remaining sub-snow amount + lake_vars%sw_down_gb(l) = dwsw_sub_snow + END IF + + END DO + + trap_frozen = 0 + trap_unfrozen = 0 + + CALL flake_interface( land_pts & + ,surft_pts(lake) & + ,ainfo%surft_index(:,lake) & + ,lake_vars%u_s_lake_gb & + ,lake_vars%surf_ht_flux_lk_gb & + ,lake_vars%sw_down_gb & + ,lake_vars%lake_depth_gb & + ,lake_vars%coriolis_param_gb & + ,timestep & + ,lake_vars%lake_albedo_gb & + ,lake_vars%lake_t_snow_gb & + ,lake_vars%lake_t_ice_gb & + ,lake_vars%lake_t_mean_gb & + ,lake_vars%lake_t_mxl_gb & + ,lake_vars%lake_shape_factor_gb & + ,lake_vars%lake_h_snow_gb & + ,lake_vars%lake_h_ice_gb & + ,lake_vars%lake_h_mxl_gb & + ,lake_vars%lake_t_sfc_gb & + ,lake_vars%ts1_lake_gb & + ,lake_vars%g_dt_gb & + ,trap_frozen & + ,trap_unfrozen ) + + IF ( trap_frozen > 0 ) THEN + errcode = -1 + CALL ereport('control', errcode, & + 'surf_couple_extra-FLake: # zero-divide (frozen) avoided =') + END IF + IF ( trap_unfrozen > 0 ) THEN + errcode = -1 + CALL ereport('control', errcode, & + 'surf_couple_extra-FLake: # zero-divide (unfrozen) avoided =') + END IF + + END IF ! Flake + + ! End of code excluded from LFRic builds +#else + ! These variables are intent out and so must be initialised for LFRic to + ! compile. Set to rmdi in the hope this will show errors if they are used! + fluxes%rrun_gb = rmdi + fluxes%rflow_gb = rmdi + psparms%z0_surft = rmdi + psparms%z0h_bare_surft = rmdi + + IF (land_pts > 0) THEN + + ! End of bespoke LFRic code +#endif + + !End of science calls + !=========================================================================== + !Expansion to full grid + + DO l = 1, land_pts + j=(ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + progs%snow_mass_ij(i,j) = snow_mass_gb(l) + END DO + + END IF ! land_pts > 0 + + !============================================================================= + !UM Diagnostics calls + + !For the UM, soil tiling has not been implemented, ie nsoilt = 1, so we can + !hard-code _soilt variables with index 1 using m = 1 +#if defined(UM_JULES) && !defined(LFRIC) + IF (model_type /= mt_single_column) THEN + m = 1 + IF (l_hydrology .AND. sf(0,8) ) THEN + CALL diagnostics_hyd( & + row_length, rows, & + land_pts, sm_levels, & + ainfo%land_index,inlandout_atm_gb, & + progs%smc_soilt(:,m), fluxes%surf_roff_gb, fluxes%sub_surf_roff_gb, & + snow_mass_gb, fluxes%snow_melt_gb, & + progs%t_soil_soilt(:,m,:), & + fluxes%snow_soil_htf, & + progs%smcl_soilt(:,m,:), & + nsurft, fluxes%hf_snow_melt_gb, psparms%sthu_soilt(:,m,:), & + psparms%sthf_soilt(:,m,:), & + fluxes%tot_tfall_gb, fluxes%melt_surft, & + land_sea_mask, & + toppdm%dun_roff_soilt(:,m), toppdm%drain_soilt(:,m), & + toppdm%qbase_soilt(:,m), & + toppdm%qbase_zw_soilt(:,m), toppdm%fch4_wetl_soilt(:,m), & + toppdm%fch4_wetl_cs_soilt(:,m), & + toppdm%fch4_wetl_npp_soilt(:,m),toppdm%fch4_wetl_resps_soilt(:,m), & + fexp_soilt(:,m),gamtot_soilt(:,m),ti_mean_soilt(:,m), & + ti_sig_soilt(:,m), & + fsat_soilt(:,m),fwetl_soilt(:,m),zw_soilt(:,m),sthzw_soilt(:,m), & + timestep, & + STASHwork8, & + sf_diag, & + !JULES TYPEs + fire_vars, progs) + END IF + + IF ( l_rivers .AND. rivers_call .AND. sf(0,26) ) THEN + CALL diagnostics_riv( & + row_length, rows, & + river_row_length, river_rows, & + riverout, & + riverout_rgrid, & + box_outflow, box_inflow, & + twatstor,inlandout_riv, & + STASHwork26 & + ) + END IF + + IF (sf(0,19)) THEN + CALL diagnostics_veg( & + row_length, rows, & + land_pts, sm_levels, & + ainfo%land_index, & + ntype,npft, & + trifctltype%c_veg_pft,trifctltype%cv_gb,trifctltype%g_leaf_phen_pft, & + trifctltype%lit_c_pft,trifctltype%lit_c_mn_gb, & + trifctltype%g_leaf_day_pft, trifctltype%lai_phen_pft, & + trifctltype%g_leaf_dr_out_pft,trifctltype%npp_dr_out_pft, & + trifctltype%resp_w_dr_out_pft,frac_agr_gb, & + disturb_veg_prev, & + ainfo%frac_surft, & + STASHwork19, & + !JULES TYPEs + trif_vars, progs, trifctltype & + ) + END IF + + END IF +#endif + + ! Deallocate water tracer working arrays + CALL wtrac_dealloc_extra(wtrac_ex) + + !============================================================================= +CASE ( cable ) + ! for testing LSM switch + WRITE(jules_message,'(A)') "CABLE not yet implemented" + CALL jules_print(RoutineName, jules_message) + + ! initialise all INTENT(OUT) for now until CABLE is implemented + fluxes%melt_surft(:,:) = 0.0 + fluxes%snowinc_surft(:,:) = 0.0 + fluxes%hf_snow_melt_gb(:) = 0.0 + fluxes%snomlt_sub_htf_gb(:) = 0.0 + fluxes%sub_surf_roff_gb(:) = 0.0 + fluxes%surf_roff_gb(:) = 0.0 + fluxes%tot_tfall_gb(:) = 0.0 + fluxes%snow_melt_gb(:) = 0.0 + fluxes%rrun_gb(:) = 0.0 + fluxes%rflow_gb(:) = 0.0 + fluxes%snow_soil_htf(:,:) = 0.0 + + work_cbl%snow_tile(:,:) = 0.0 + work_cbl%lying_snow(:) = 0.0 + work_cbl%surf_roff(:) = 0.0 + work_cbl%sub_surf_roff(:) = 0.0 + work_cbl%tot_tfall(:) = 0.0 + +CASE DEFAULT + errorstatus = 101 + WRITE(jules_message,'(A,I0)') 'Unrecognised surface scheme. lsm_id = ', & + lsm_id + CALL ereport(RoutineName, errorstatus, jules_message) + +END SELECT + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE surf_couple_extra +END MODULE surf_couple_extra_mod diff --git a/src/io/model_interface/surf_couple_implicit_mod.F90 b/src/io/model_interface/surf_couple_implicit_mod.F90 new file mode 100644 index 00000000..b917440d --- /dev/null +++ b/src/io/model_interface/surf_couple_implicit_mod.F90 @@ -0,0 +1,622 @@ +! *****************************COPYRIGHT**************************************** +! (c) Crown copyright, Met Office. All rights reserved. +! +! This routine has been licensed to the other JULES partners for use and +! distribution under the JULES collaboration agreement, subject to the terms and +! conditions set out therein. +! +! [Met Office Ref SC0237] +! *****************************COPYRIGHT**************************************** + +MODULE surf_couple_implicit_mod + +USE jules_grid_update_implicit_mod, ONLY: jules_grid_update_implicit +USE jules_land_sf_implicit_mod, ONLY: jules_land_sf_implicit +USE jules_ssi_sf_implicit_mod, ONLY: jules_ssi_sf_implicit +USE jules_griddiag_sf_implicit_mod, ONLY: jules_griddiag_sf_implicit + +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +PRIVATE +PUBLIC :: surf_couple_implicit + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='SURF_COUPLE_IMPLICIT_MOD' + +CONTAINS + +!=============================================================================== +! Public subroutine +!=============================================================================== +SUBROUTINE surf_couple_implicit( & + !Important switch + l_correct, & + !Forcing INTENT(IN) + u_1, v_1, & + !Misc INTENT(IN) Many of these simply come out of explicit and into here. + rhokm_u, rhokm_v, gamma1, gamma2, alpha1, alpha1_sea, alpha1_sice, & + ashtf, ashtf_sea, ashtf_surft, du, dv, fracaero_t, fracaero_s, resfs, resft, & + rhokh, rhokh_surft, rhokh_sice, rhokh_sea, & + z0hssi, z0mssi, chr1p5m, & + chr1p5m_sice, canhc_surft, flake, tile_frac, wt_ext_surft, & + cdr10m_u, cdr10m_v, r_gamma, & + !Diagnostics, INTENT(INOUT) + sf_diag, & + !Fluxes INTENT(INOUT) + fqw_1, ftl_1, & + !Misc INTENT(INOUT) + epot_surft, dtstar_surft, dtstar_sea, dtstar_sice, radnet_sice, olr, & + !Fluxes INTENT(OUT) + taux_1, tauy_1, & + !Misc INTENT(OUT) + ERROR, & + !UM-only arguments + !JULES ancil_info module + nsurft, land_pts, ssi_pts, sice_pts, sea_pts, surft_pts, & + !JULES coastal module + flandg, & + ! Coastal OUT - do this here as needs to be OUT + tstar_land_ij, tstar_sice_ij, & + ! Water tracer switch (IN) + l_wtrac_bl, & + !JULES switches module + l_co2_interactive, l_mr_physics, & + co2_3d_ij, & + !Arguments without a JULES module + ctctq1,dqw1_1,dtl1_1,du_star1,dv_star1,cq_cm_u_1,cq_cm_v_1,flandg_u, & + flandg_v, rho1, f3_at_p, uStarGBM,tscrndcl_ssi,tscrndcl_surft,tStbTrans, & + rhokh_mix, ti_gb, sky, & + !TYPES containing field data (IN OUT) + crop_vars,ainfo,aerotype,progs,coast, jules_vars, & + fluxes, & + lake_vars, & + forcing, & + !rivers, & + !veg3_parm, & + !veg3_field, & + !chemvars, & + wtrac_jls, & + progs_cbl, & + work_cbl & + ) + +!Module Imports + +!TYPE definitions +USE crop_vars_mod, ONLY: crop_vars_type +USE ancil_info, ONLY: ainfo_type +USE aero, ONLY: aero_type +USE prognostics, ONLY: progs_type +USE coastal, ONLY: coastal_type +USE jules_vars_mod, ONLY: jules_vars_type +USE fluxes_mod, ONLY: fluxes_type +USE lake_mod, ONLY: lake_type +USE jules_forcing_mod, ONLY: forcing_type +! USE jules_rivers_mod, ONLY: rivers_type +! USE veg3_parm_mod, ONLY: in_dev +! USE veg3_field_mod, ONLY: in_dev +! USE jules_chemvars_mod, ONLY: chemvars_type +USE jules_wtrac_type_mod, ONLY: jls_wtrac_type + + + +! In general CABLE utilizes a required subset of tbe JULES types, however; +USE progs_cbl_vars_mod, ONLY: progs_cbl_vars_type ! CABLE requires extra progs +USE work_vars_mod_cbl, ONLY: work_vars_type ! and some kept thru timestep + +!Common modules +USE ereport_mod, ONLY: ereport + +USE jules_soil_mod, ONLY: sm_levels +USE jules_sea_seaice_mod, ONLY: nice_use, nice +USE sf_diags_mod, ONLY: strnewsfdiag +USE jules_surface_types_mod, ONLY: ntype +USE jules_water_tracers_mod, ONLY: l_wtrac_imp_jls, l_wtrac_jls, & + n_wtrac_jls, n_evap_srce + +!Potential troublemakers +USE atm_fields_bounds_mod, ONLY: tdims, udims, vdims, pdims, & + tdims_s, udims_s, vdims_s + +! Module switches name between UM and JULES-standalone +#if defined(UM_JULES) +USE rad_input_mod, ONLY: co2_mmr +#else +USE aero, ONLY: & + co2_mmr +#endif + +USE jules_model_environment_mod, ONLY: & + lsm_id, jules, cable + +USE parkind1, ONLY: jprb, jpim +USE yomhook, ONLY: lhook, dr_hook + +! for testing LSM switch +USE jules_print_mgr, ONLY: jules_message, jules_print + +IMPLICIT NONE + + +!----------------------------------------------------------------------------- +! Description: +! Coupling routine between the UM or JULES system code and land surface +! implicit science routines. Calls the appropriate LSM-specific code. +! +! Code Owner: Please refer to ModuleLeaders.txt +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + +!Arguments + +!UM-only arguments +!JULES ancil_info module +INTEGER, INTENT(IN) :: & + nsurft, land_pts, ssi_pts, sice_pts, sea_pts, & + ! ainfo%surft_index(nsurft), & + surft_pts(ntype) +!JULES water_tracer module +LOGICAL, INTENT(IN) :: & + l_wtrac_bl ! Logical to control the water tracer calculations in + ! the implicit code +!JULES coastal module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + flandg(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) +!JULES switches module +LOGICAL, INTENT(IN) :: & + l_co2_interactive, l_mr_physics + +REAL(KIND=real_jlslsm), INTENT(IN) :: & + co2_3d_ij(tdims_s%i_start:tdims_s%i_end,tdims_s%j_start:tdims_s%j_end) + +!Arguments without a JULES module +REAL(KIND=real_jlslsm), INTENT(IN) :: & + ctctq1(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end), & + dqw1_1(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end), & + dtl1_1(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end), & + du_star1(udims_s%i_start:udims_s%i_end,udims_s%j_start:udims_s%j_end), & + dv_star1(vdims_s%i_start:vdims_s%i_end,vdims_s%j_start:vdims_s%j_end), & + cq_cm_u_1(udims%i_start:udims%i_end,udims%j_start:udims%j_end), & + ! Coefficient in U tri-diagonal implicit matrix + cq_cm_v_1(vdims%i_start:vdims%i_end,vdims%j_start:vdims%j_end), & + ! Coefficient in V tri-diagonal implicit matrix + flandg_u(udims%i_start:udims%i_end,udims%j_start:udims%j_end), & + !Land frac (on U-grid, with 1st and last rows undefined or, at present, + !set to "missing data") + flandg_v(vdims%i_start:vdims%i_end,vdims%j_start:vdims%j_end), & + !Land frac (on V-grid, with 1st and last rows undefined or, at present, + !set to "missing data") + rho1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Density on lowest level + f3_at_p(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Coriolis parameter + sky(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + ! Skyview correction factor for surface LW + ustargbm(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + +REAL(KIND=real_jlslsm), INTENT(IN OUT) :: & + ! BM surface friction velocity + tscrndcl_ssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Decoupled screen-level temperature over sea or sea-ice + tscrndcl_surft(land_pts,nsurft), & + !Decoupled screen-level temperature over land tiles + tstbtrans(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + tstar_land_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + tstar_sice_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + rhokh_mix(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Exchange coeffs for moisture. + ti_gb(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + +! Important switch +LOGICAL, INTENT(IN) :: & + l_correct ! flag used by the new BL solver + +!Forcing INTENT(IN) +REAL(KIND=real_jlslsm), INTENT(IN) :: & + u_1(udims_s%i_start:udims_s%i_end,udims_s%j_start:udims_s%j_end), & + !W'ly wind component (m/s) + v_1(vdims_s%i_start:vdims_s%i_end,vdims_s%j_start:vdims_s%j_end) + !S'ly wind component (m/s) + +!Misc INTENT(IN) Many of these simply come out of explicit and come +!back into here. Need a module for them. + +!Arrays +REAL(KIND=real_jlslsm), INTENT(IN) :: & + rhokm_u(udims%i_start:udims%i_end,udims%j_start:udims%j_end), & + !Exchange coefficients for momentum (on U-grid, with 1st and last rows + !undefined or, at present, set to "missing data") + rhokm_v(vdims%i_start:vdims%i_end,vdims%j_start:vdims%j_end), & + !Exchange coefficients for momentum (on V-grid, with 1st and last rows + !undefined or, at present, set to "missing data") + gamma1(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end), & + !weights for new BL solver + gamma2(pdims%i_start:pdims%i_end,pdims%j_start:pdims%j_end), & + alpha1(land_pts,nsurft), & + !Mean gradient of saturated specific humidity with respect to temperature + !between the bottom model layer and tile surfaces + alpha1_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !ALPHA1 for sea. + alpha1_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !ALPHA1 for sea-ice. + ashtf(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Adjusted SEB coefficient for sea-ice + ashtf_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Adjusted SEB coefficient for sea + ashtf_surft(land_pts,nsurft), & + !Adjusted SEB coefficient for land tiles. + du(udims_s%i_start:udims_s%i_end,udims_s%j_start:udims_s%j_end), & + !Level 1 increment to u wind field + dv(vdims_s%i_start:vdims_s%i_end,vdims_s%j_start:vdims_s%j_end), & + !Level 1 increment to v wind field + fracaero_t(land_pts,nsurft), & + !Total fraction of surface moisture flux with only aerodynamic resistance + fracaero_s(land_pts,nsurft), & + !Fraction of surface moisture flux with only aerodynamic resistance + !from the frozen portion of the tile + resfs(land_pts,nsurft), & + !Combined soil, stomatal and aerodynamic resistance factor for fraction + !(1-fracaero_t) of snow-free land tiles. + resft(land_pts,nsurft), & + !Total resistance factor. fracaero_t+(1-fracaero_t)*resfs + !for snow-free land, 1 for snow. + rhokh(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Grid-box surface exchange coefficients (not used for JULES) + rhokh_surft(land_pts,nsurft), & + !Surface exchange coefficients for land tiles + rhokh_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Surface exchange coefficients for sea sea-ice + rhokh_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Surface exchange coefficients for sea + z0hssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + z0mssi(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Roughness lengths over sea (m) + chr1p5m(land_pts,nsurft), & + !Ratio of coefffs for calculation of 1.5m temp for land tiles. + chr1p5m_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !CHR1P5M for sea and sea-ice (leads ignored). + canhc_surft(land_pts,nsurft), & + !Areal heat capacity of canopy for land tiles (J/K/m2). + flake(land_pts,nsurft), & + !Lake fraction. + tile_frac(land_pts,nsurft), & + !Tile fractions including snow cover in the ice tile. + wt_ext_surft(land_pts,sm_levels,nsurft), & + !Fraction of evapotranspiration extracted from each soil layer by each tile. + cdr10m_u(udims%i_start:udims%i_end,udims%j_start:udims%j_end), & + !Ratio of CD's reqd for calculation of 10 m wind. On U-grid; comments as + !per RHOKM. + cdr10m_v(vdims%i_start:vdims%i_end,vdims%j_start:vdims%j_end) + !Ratio of CD's reqd for calculation of 10 m wind. On V-grid; comments as + +! Implicit weighting for B.L. +REAL(KIND=real_jlslsm), INTENT(IN) :: r_gamma + +!diagnostic array +TYPE (strnewsfdiag), INTENT(IN OUT) :: sf_diag + +!Fluxes INTENT(INOUT) +REAL(KIND=real_jlslsm), INTENT(IN OUT) :: & + fqw_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Moisture flux between layers (kg per square metre per sec) FQW(,1) is + !total water flux from surface, 'E'. + ftl_1(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + !FTL(,K) contains net turbulent sensible heat flux into layer K from below; + !so FTL(,1) is the surface sensible heat, H.(W/m2) + +!Misc INTENT(INOUT) +REAL(KIND=real_jlslsm), INTENT(IN OUT) :: & + epot_surft(land_pts,nsurft), & + !surface tile potential evaporation + dtstar_surft(land_pts,nsurft), & + !Change in TSTAR over timestep for land tiles + dtstar_sea(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end), & + !Change is TSTAR over timestep for open sea + dtstar_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Change is TSTAR over timestep for sea-ice + radnet_sice(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end,nice_use), & + !Surface net radiation on sea-ice (W/m2) + olr(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) + ! IN TOA - surface upward LW on last radiation timestep + ! OUT Corrected TOA outward LW + +!Fluxes INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + taux_1(udims%i_start:udims%i_end,udims%j_start:udims%j_end), & + !W'ly component of surface wind stress (N/sq m). (On UV-grid with first + !and last rows undefined or, at present, set to missing data + tauy_1(vdims%i_start:vdims%i_end,vdims%j_start:vdims%j_end) + !S'ly component of surface wind stress (N/sq m). On UV-grid; comments as + !per TAUX + +!Misc INTENT(OUT) +INTEGER, INTENT(OUT) :: & + ERROR !0 - AOK; 1 to 7 - bad grid definition detected + +TYPE(jls_wtrac_type), INTENT(IN OUT) :: wtrac_jls +!TYPES containing field data (IN OUT) +TYPE(crop_vars_type), INTENT(IN OUT) :: crop_vars +TYPE(ainfo_type), INTENT(IN OUT) :: ainfo +TYPE(aero_type), INTENT(IN OUT) :: aerotype +TYPE(progs_type), INTENT(IN OUT) :: progs +TYPE(coastal_type), INTENT(IN OUT) :: coast +TYPE(jules_vars_type), INTENT(IN OUT) :: jules_vars +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes +TYPE(lake_type), INTENT(IN OUT) :: lake_vars +TYPE(forcing_type), INTENT(IN OUT) :: forcing +!TYPE(rivers_type), INTENT(IN OUT) :: rivers +!TYPE(in_dev), INTENT(IN OUT) :: veg3_parm +!TYPE(in_dev), INTENT(IN OUT) :: veg3_field +!TYPE(chemvars_type), INTENT(IN OUT) :: chemvars + +!CABLE TYPES containing field data (IN OUT) +TYPE(progs_cbl_vars_type), INTENT(IN OUT) :: progs_cbl +TYPE(work_vars_type), INTENT(IN OUT) :: work_cbl + +!----------------------------------------------------------------------------- +! Workspace +!----------------------------------------------------------------------------- +REAL(KIND=real_jlslsm) :: & + ice_fract_cat_use(tdims%i_start:tdims%i_end, & + tdims%j_start:tdims%j_end,nice_use) & + ! Sea ice category fractions + ! If nice_use=1, this is the total ice + ! fraction +,surf_ht_flux_sice_sm_ij(tdims%i_start:tdims%i_end, & + tdims%j_start:tdims%j_end) & + ! Sea area mean seaice surface heat flux +,ei_sice_sm_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Sea area mean sea ice sublimation +,ei_land_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Sublimation from lying snow + ! (kg/m2/s). +,tstar_ssi_old_ij(tdims%i_start:tdims%i_end,tdims%j_start:tdims%j_end) & + ! Sea and sea-ice surface temperature + ! at beginning of timestep -- + ! Required only for decoupled diagnosis, + ! so allocatable, and local since it is + ! used only on the predictor step +,tstar_surft_old(land_pts,nsurft) + ! Tile surface temperatures at + ! beginning of timestep. + +!----------------------------------------------------------------------------- +! Local variables +!----------------------------------------------------------------------------- + +!Dr Hook variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='SURF_COUPLE_IMPLICIT' + +!----------------------------------------------------------------------------- +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +! Set JULES logical to control water tracer calculations in the implicit +! scheme (they only need to be done on the final dynamics loop in UMJULES) +l_wtrac_imp_jls = (l_wtrac_bl .AND. l_wtrac_jls) + +SELECT CASE( lsm_id ) +CASE ( jules ) + + CALL jules_grid_update_implicit( & + ! IN values defining field dimensions and subset to be processed : + land_pts,ainfo%land_index,nice,nice_use,nsurft,ainfo%surft_index, & + surft_pts, tile_frac,flandg, & + ! IN sea/sea-ice data : + ainfo%ice_fract_ij,ainfo%ice_fract_ncat_sicat, & + ! IN everything not covered so far : + rhokm_u,rhokm_v,r_gamma, & + gamma1,gamma2,alpha1,alpha1_sea,alpha1_sice, & + ashtf,ashtf_sea,ashtf_surft, & + du,dv,resft,rhokh_surft,rhokh_sice,rhokh_sea,ctctq1, & + dqw1_1,dtl1_1,du_star1,dv_star1,cq_cm_u_1,cq_cm_v_1, & + l_correct,flandg_u,flandg_v,progs%snow_surft, & + ! INOUT data : + epot_surft,fluxes%fqw_sicat,fluxes%ftl_sicat,dtstar_surft, & + dtstar_sea,dtstar_sice, & + fluxes%fqw_surft, fqw_1,ftl_1,fluxes%ftl_surft, & + coast%taux_land_ij,coast%taux_ssi_ij,coast%tauy_land_ij, & + coast%tauy_ssi_ij, coast%taux_land_star,coast%tauy_land_star, & + coast%taux_ssi_star, coast%tauy_ssi_star, & + ! OUT data required elsewhere in UM system : + fluxes%e_sea_ij,fluxes%h_sea_ij,taux_1,tauy_1,ice_fract_cat_use & + ) + + IF ( .NOT. l_correct ) THEN + CALL jules_land_sf_implicit ( & + ! IN values defining field dimensions and subset to be processed : + land_pts,ainfo%land_index,nsurft,ainfo%surft_index,surft_pts, & + sm_levels, & + canhc_surft,progs%canopy_surft,flake,progs%smc_soilt,tile_frac, & + wt_ext_surft,coast%fland,flandg, & + !IN values defining water tracer field dimensions + n_wtrac_jls, n_evap_srce, & + ! IN everything not covered so far : + forcing%lw_down_ij,fluxes%sw_surft,sky,progs%t_soil_soilt,r_gamma, & + alpha1,ashtf_surft, & + jules_vars%dtrdz_charney_grid_1_ij,fracaero_t,fracaero_s, & + resfs,resft, rhokh_surft, & + fluxes%emis_surft,fluxes%fsnow,progs%snow_surft,dtstar_surft, & + ! INOUT data : + progs%tstar_surft,fluxes%fqw_surft,fqw_1,ftl_1,fluxes%ftl_surft, & + sf_diag, & + ! OUT Diagnostic not requiring STASH flags : + fluxes%ecan_ij,fluxes%ei_surft,fluxes%esoil_surft, & + coast%surf_ht_flux_land_ij, & + ei_land_ij,fluxes%surf_htf_surft, & + ! OUT data required elsewhere in UM system : + tstar_land_ij,fluxes%le_surft,fluxes%radnet_surft, & + fluxes%ecan_surft,fluxes%esoil_ij_soilt, & + fluxes%ext_soilt,fluxes%melt_surft,fluxes%snowinc_surft, & + tstar_surft_old,ERROR, & + !New arguments replacing USE statements + ! lake_mod (IN) + lake_vars%lake_h_ice_gb, & + ! lake_mod (OUT) + lake_vars%surf_ht_flux_lake_ij, lake_vars%non_lake_frac, & + ! fluxes (IN) + fluxes%anthrop_heat_surft, & + ! fluxes (OUT) + fluxes%surf_ht_store_surft, & + fluxes%lake_evap, & + ! c_elevate (IN) + jules_vars%lw_down_elevcorr_surft, & + ! prognostics (IN) + progs%nsnow_surft, & + ! jules_vars_mod (IN) + jules_vars%snowdep_surft, & + !TYPES containing field data (IN OUT) + crop_vars, & + !Water tracers (IN) + wtrac_jls%snow_surft, wtrac_jls%smc_soilt, & + wtrac_jls%canopy_surft, wtrac_jls%fqw_evapsrce, & + !Water tracers (INOUT) + wtrac_jls%fqw_surft, & + !Water tracers (OUT) + wtrac_jls%ei_surft, wtrac_jls%ei_ij, & + wtrac_jls%esoil_surft, wtrac_jls%esoil_ij_soilt, & + wtrac_jls%ext_soilt, wtrac_jls%ecan_surft, & + wtrac_jls%ecan_ij, wtrac_jls%lake_evap, & + wtrac_jls%fqw_1, wtrac_jls%dfqw_imp) + + END IF ! IF .NOT. L_correct + + CALL jules_ssi_sf_implicit ( & + ! IN values defining field dimensions and subset to be processed : + nice,nice_use,flandg, & + !IN values defining water tracer field dimensions + n_wtrac_jls, & + ! IN sea/sea-ice data : + ainfo%ice_fract_ij,ainfo%ice_fract_ncat_sicat,ice_fract_cat_use, & + progs%k_sice_sicat,progs%di_ncat_sicat,ainfo%sstfrz_ij, & + ! IN everything not covered so far : + forcing%lw_down_ij,r_gamma,alpha1_sice,ashtf, & + jules_vars%dtrdz_charney_grid_1_ij,rhokh_sice,l_correct, & + ! INOUT data : + fluxes%fqw_sicat,fluxes%ftl_sicat,coast%tstar_sice_sicat, & + coast%tstar_ssi_ij, & + coast%tstar_sea_ij, & + radnet_sice,fqw_1,ftl_1,progs%ti_sicat,sf_diag, & + ! OUT Diagnostic not requiring STASH flags : + ti_gb,fluxes%sea_ice_htf_sicat, & + fluxes%surf_ht_flux_sice, & + ! OUT data required elsewhere in UM system : + tstar_sice_ij,fluxes%e_sea_ij,fluxes%h_sea_ij,fluxes%ei_sice, & + dtstar_sea,dtstar_sice, & + surf_ht_flux_sice_sm_ij,ei_sice_sm_ij,tstar_ssi_old_ij, & + fluxes%sice_melt, & + ! ancil_info (IN) + ainfo%ssi_index, ainfo%sice_index, ainfo%sice_index_ncat, & + ainfo%fssi_ij, ainfo%sice_frac, ainfo%sice_frac_ncat, & + ainfo%ocn_cpl_point, ainfo%sice_pts_ncat, ainfo%sea_index, & + ssi_pts, sice_pts, sea_pts, & + ! fluxes (IN) + fluxes%sw_sicat, & + ! water tracers (IN) + wtrac_jls%fqw_sea, wtrac_jls%fqw_sicat, & + ! water tracers (INOUT) + wtrac_jls%fqw_1, wtrac_jls%dfqw_imp, & + ! water tracers (OUT) + wtrac_jls%e_sea_ij, wtrac_jls%ei_sice_ij) + + CALL jules_griddiag_sf_implicit ( & + ! IN values defining field dimensions and subset to be processed : + land_pts,ainfo%land_index,nice_use,nsurft,ainfo%surft_index, & + surft_pts, tile_frac,flandg,l_mr_physics, & + ! IN sea/sea-ice data : + ainfo%ice_fract_ij,ice_fract_cat_use,forcing%u_0_ij, & + forcing%v_0_ij, & + ! IN everything not covered so far : + forcing%pstar_ij,forcing%qw_1_ij,forcing%tl_1_ij, & + u_1,v_1,du,dv, & + resft,rhokh,ainfo%z1_tq_ij, & + z0hssi,z0mssi,fluxes%z0h_surft,fluxes%z0m_surft, & + cdr10m_u,cdr10m_v, & + chr1p5m,chr1p5m_sice,ctctq1, & + dqw1_1,dtl1_1,du_star1,dv_star1,cq_cm_u_1,cq_cm_v_1, & + l_correct,fluxes%emis_surft, & + coast%tstar_sice_sicat,coast%tstar_ssi_ij,progs%tstar_surft, & + coast%tstar_sea_ij,fqw_1, & + coast%surf_ht_flux_land_ij,surf_ht_flux_sice_sm_ij,ei_land_ij, & + ei_sice_sm_ij,tstar_land_ij,tstar_ssi_old_ij,tstar_surft_old, & + taux_1,tauy_1, & + ! IN variables used to calculate cooling at the screen level + l_co2_interactive, co2_mmr, co2_3d_ij,rho1, f3_at_p, & + ustargbm, & + ! INOUT data : + ftl_1,olr, & + tscrndcl_ssi,tscrndcl_surft,tStbTrans,sf_diag, & + ! OUT Diagnostic not requiring STASH flags : + fluxes%surf_ht_flux_ij, & + ! OUT data required elsewhere in UM system : + fluxes%tstar_ij,fluxes%ei_ij,rhokh_mix & + ) + +CASE ( cable ) + ! for testing LSM switch + WRITE(jules_message,'(A)') "CABLE not yet implemented" + CALL jules_print(RoutineName, jules_message) + + ! initialise all INTENT(OUT) for now until CABLE is implemented + fluxes%tstar_ij(:,:) = 0.0 + fluxes%le_surft(:,:) = 0.0 + fluxes%radnet_surft(:,:) = 0.0 + fluxes%e_sea_ij(:,:) = 0.0 + fluxes%h_sea_ij(:,:) = 0.0 + taux_1(:,:) = 0.0 + tauy_1(:,:) = 0.0 + fluxes%ecan_surft(:,:) = 0.0 + fluxes%ei_ij(:,:) = 0.0 + fluxes%esoil_ij_soilt(:,:,:) = 0.0 + fluxes%ext_soilt(:,:,:) = 0.0 + fluxes%melt_surft(:,:) = 0.0 + fluxes%snowinc_surft(:,:) = 0.0 + fluxes%ecan_ij(:,:) = 0.0 + fluxes%ei_surft(:,:) = 0.0 + fluxes%esoil_surft(:,:) = 0.0 + fluxes%sea_ice_htf_sicat(:,:,:) = 0.0 + fluxes%surf_ht_flux_ij(:,:) = 0.0 + fluxes%surf_htf_surft(:,:) = 0.0 + ERROR = 0 + + progs_cbl%SoilTemp_CABLE(:,:,:) = 0.0 + progs_cbl%SoilMoisture_CABLE(:,:,:) = 0.0 + progs_cbl%FrozenSoilFrac_CABLE(:,:,:) = 0.0 + progs_cbl%SnowDepth_CABLE(:,:,:) = 0.0 + progs_cbl%SnowMass_CABLE(:,:,:) = 0.0 + progs_cbl%SnowTemp_CABLE(:,:,:) = 0.0 + progs_cbl%SnowDensity_CABLE(:,:,:) = 0.0 + progs_cbl%ThreeLayerSnowFlag_CABLE(:,:) = 0.0 + progs_cbl%OneLyrSnowDensity_CABLE(:,:) = 0.0 + progs_cbl%SnowAge_CABLE(:,:) = 0.0 + progs_cbl%snowOsurft(:,:) = 0.0 + + work_cbl%snow_tile(:,:) = 0.0 + work_cbl%lying_snow(:) = 0.0 + work_cbl%surf_roff(:) = 0.0 + work_cbl%sub_surf_roff(:) = 0.0 + work_cbl%tot_tfall(:) = 0.0 + +CASE DEFAULT + ERROR = 101 + WRITE(jules_message,'(A,I0)') 'Unrecognised surface scheme. lsm_id = ', & + lsm_id + CALL ereport(RoutineName, ERROR, jules_message) + +END SELECT + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE surf_couple_implicit + +END MODULE surf_couple_implicit_mod diff --git a/src/io/model_interface/surf_couple_radiation_mod.F90 b/src/io/model_interface/surf_couple_radiation_mod.F90 new file mode 100644 index 00000000..b6260df8 --- /dev/null +++ b/src/io/model_interface/surf_couple_radiation_mod.F90 @@ -0,0 +1,333 @@ +! *****************************COPYRIGHT**************************************** +! (c) Crown copyright, Met Office. All rights reserved. +! +! This routine has been licensed to the other JULES partners for use and +! distribution under the JULES collaboration agreement, subject to the terms and +! conditions set out therein. +! +! [Met Office Ref SC0237] +! *****************************COPYRIGHT**************************************** + +MODULE surf_couple_radiation_mod + +USE um_types, ONLY: real_jlslsm + +IMPLICIT NONE + +PRIVATE +PUBLIC :: surf_couple_radiation + +CHARACTER(LEN=*), PARAMETER, PRIVATE :: ModuleName='SURF_COUPLE_RADIATION_MOD' + +CONTAINS + +!=============================================================================== +! Public subroutine +!=============================================================================== +SUBROUTINE surf_couple_radiation( & + !Misc INTENT(IN) + ws10m, chloro, & + n_band, max_n_swbands, wavelength_short, wavelength_long, & + !Misc INTENT(OUT) + sea_ice_albedo, & + !(ancil_info mod) + nsurft, land_pts, sea_pts, surft_pts, row_length, rows, & + !(coastal mod) + flandg, & + !(prognostics mod) + ! Warning- snow_surft causes problems in the UM and passed as an array here + snow_surft, & + !INTENT(OUT) + albobs_sc_ij, open_sea_albedo, & + !TYPES containing field data (IN OUT) + psparms,ainfo,urban_param,progs,coast,jules_vars, & + fluxes, & + lake_vars, & + !forcing, & + !rivers, & + !veg3_parm, & + !veg3_field, & + !chemvars, & + progs_cbl & + ) + +!Module imports + +!TYPE definitions +USE p_s_parms, ONLY: psparms_type +USE ancil_info, ONLY: ainfo_type +USE urban_param_mod, ONLY: urban_param_type +USE prognostics, ONLY: progs_type +USE coastal, ONLY: coastal_type +USE jules_vars_mod, ONLY: jules_vars_type +USE fluxes_mod, ONLY: fluxes_type +USE lake_mod, ONLY: lake_type +! USE forcing, ONLY: forcing_type +! USE jules_rivers_mod, ONLY: rivers_type +! USE veg3_parm_mod, ONLY: in_dev +! USE veg3_field_mod, ONLY: in_dev +! USE jules_chemvars_mod, ONLY: chemvars_type + +USE jules_ssi_albedo_mod, ONLY: jules_ssi_albedo +USE jules_land_albedo_mod, ONLY: jules_land_albedo +USE cable_land_albedo_mod, ONLY: cable_land_albedo + +!Common modules +USE ereport_mod, ONLY: & + ereport + +USE jules_sea_seaice_mod, ONLY: & + nice, nice_use, & + alpham, alphac, alphab, dtice, dt_bare, dalb_bare_wet, & + pen_rad_frac, sw_beta, & + albicev_cice, albicei_cice, albsnowv_cice, albsnowi_cice, & + albpondv_cice, albpondi_cice, & + ahmax, dalb_mlt_cice, dalb_mlts_v_cice, dalb_mlts_i_cice, & + dt_bare_cice, dt_snow_cice, pen_rad_frac_cice, sw_beta_cice, & + snowpatch + +!Potential troublemaker +USE theta_field_sizes, ONLY: & + t_i_length, t_j_length + +!for testing LSM switch +USE jules_print_mgr, ONLY: jules_message, jules_print + +USE jules_model_environment_mod, ONLY: & + lsm_id, jules, cable + +! In general CABLE utilizes a required subset of tbe JULES types, however; +USE progs_cbl_vars_mod, ONLY: progs_cbl_vars_type ! CABLE requires extra progs +USE cable_fields_mod, ONLY: pars_io_cbl ! and veg/soil parameters +!data: constants +USE cable_surface_types_mod, ONLY: ICE_cable, lakes_cable + +USE jules_surface_types_mod, ONLY: npft + +!Dr Hook +USE parkind1, ONLY: & + jprb, jpim + +USE yomhook, ONLY: & + lhook, dr_hook + + +IMPLICIT NONE + +!----------------------------------------------------------------------------- +! Description: +! Coupling routine between the UM or JULES system code and land surface +! radiation science routines. Calls the appropriate LSM-specific code. +! +! Some variables exist in modules only in JULES, others only in the UM +! Options in order of preference +! -UM and JULES share the same module names. +! -UM and JULES have different module names and USE statements go on an ifdef +! -The UM flavour of the variable does not live in a module. Pass in using a +! an ifdef'ed argument list +! +! If there are lots of ifs, then we could cosider splitting lsm_couple into +! jules_couple and lsm_couple +! +! Code Owner: Please refer to ModuleLeaders.txt +! +! Code Description: +! Language: Fortran 90. +! This code is written to JULES coding standards v1. +!----------------------------------------------------------------------------- + + +! Subroutine arguments + +! Dimensioning variables +INTEGER, INTENT(IN) :: & + n_band, & + max_n_swbands + + !UM-only args: INTENT(IN) + !(ancil_info mod) +INTEGER, INTENT(IN):: & + nsurft, land_pts, sea_pts, surft_pts(nsurft), row_length, rows + +!(coastal mod) +REAL(KIND=real_jlslsm), INTENT(IN) :: & +flandg(row_length, rows) + +!(prognostics mod) +REAL(KIND=real_jlslsm), INTENT(IN) :: & + snow_surft(land_pts, nsurft) + +!UM-only args: INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + albobs_sc_ij(t_i_length,t_j_length,nsurft,2), & + !albedo scaling factors to obs + open_sea_albedo(row_length,rows,2,max_n_swbands) + !Surface albedo for Open Sea (direct and diffuse components, for each + !band, with zeros for safety where no value applies) + +!Misc INTENT(IN) +REAL(KIND=real_jlslsm), INTENT(IN) :: & + ws10m(row_length,rows), & + !10m wind speed + chloro(row_length,rows) !nr surface chlorophyll content + +REAL(KIND=real_jlslsm), INTENT(IN) :: & + wavelength_short(n_band), & + wavelength_long(n_band) + +!Misc INTENT(OUT) +REAL(KIND=real_jlslsm), INTENT(OUT) :: & + sea_ice_albedo(row_length,rows,4) !Surface Albedo for sea ice + ! (*,1) - direct beam visible + ! (*,2) - diffuse visible + ! (*,3) - direct beam near-ir + ! (*,4) - diffuse near-ir + +!TYPES containing field data (IN OUT) +TYPE(psparms_type), INTENT(IN OUT) :: psparms +TYPE(ainfo_type), INTENT(IN OUT) :: ainfo +TYPE(urban_param_type), INTENT(IN OUT) :: urban_param +TYPE(progs_type), INTENT(IN OUT) :: progs +TYPE(coastal_type), INTENT(IN OUT) :: coast +TYPE(jules_vars_type), INTENT(IN OUT) :: jules_vars +TYPE(fluxes_type), INTENT(IN OUT) :: fluxes +TYPE(lake_type), INTENT(IN OUT) :: lake_vars +!TYPE(forcing_type), INTENT(IN OUT) :: forcing +!TYPE(rivers_type), INTENT(IN OUT) :: rivers +!TYPE(in_dev), INTENT(IN OUT) :: veg3_parm +!TYPE(in_dev), INTENT(IN OUT) :: veg3_field +!TYPE(chemvars_type), INTENT(IN OUT) :: chemvars + +!CABLE TYPES containing field data +TYPE(progs_cbl_vars_type), INTENT(IN) :: progs_cbl + +!----------------------------------------------------------------------------- +! Local variables +!----------------------------------------------------------------------------- + +!Land point only versions of ij variables +REAL(KIND=real_jlslsm) :: soot_gb(land_pts) +REAL(KIND=real_jlslsm) :: cosz_gb(land_pts) + +!Counters +INTEGER :: i,j,l +INTEGER :: errorstatus + +!Dr Hook variables +INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 +INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 +REAL(KIND=jprb) :: zhook_handle + +CHARACTER(LEN=*), PARAMETER :: RoutineName='SURF_COUPLE_RADIATION' + +!----------------------------------------------------------------------------- +!End of header + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_in,zhook_handle) + +SELECT CASE( lsm_id ) +CASE ( jules ) + + !Compress gridded variables to land point only + DO l = 1, land_pts + j = (ainfo%land_index(l) - 1) / row_length + 1 + i = ainfo%land_index(l) - (j-1) * row_length + soot_gb(l) = progs%soot_ij(i,j) + cosz_gb(l) = psparms%cosz_ij(i,j) + END DO + + CALL jules_land_albedo ( & + !INTENT(IN) + t_i_length * t_j_length, & + land_pts, nsurft, & + ainfo%land_index, surft_pts, ainfo%surft_index, & + psparms%albsoil_soilt, psparms%albobs_sw_gb, psparms%albobs_vis_gb, & + psparms%albobs_nir_gb, & + cosz_gb, soot_gb, jules_vars%ho2r2_orog_gb, & + progs%lai_pft, progs%canht_pft, & + progs%rgrain_surft, snow_surft, progs%tstar_surft, psparms%z0_surft, & + ainfo%frac_surft, & + !INTENT(OUT) + fluxes%alb_surft,albobs_sc_ij,fluxes%fsnow,fluxes%land_albedo_ij, & + !New arguments replacing USE statements + !jules_vars_mod (IN OUT) + jules_vars%albobs_scaling_surft, & + !jules_vars_mod (OUT) + jules_vars%snowdep_surft, & + !urban_param (IN) + urban_param%albwl_gb, urban_param%albrd_gb, urban_param%hwr_gb, & + !lake_mod (IN) + lake_vars%lake_h_ice_gb, & + !ancil_info (IN) + ainfo%l_lice_point, ainfo%l_lice_surft, & + !prognostics (IN) + progs%snowdepth_surft, progs%rho_snow_grnd_surft, progs%nsnow_surft, & + progs%sice_surft, progs%sliq_surft, progs%ds_surft) + + +CASE ( cable ) + + CALL cable_land_albedo ( & + !OUT: (per rad band) albedos [GridBoxMean & per tile albedo] + fluxes%land_albedo_ij, fluxes%alb_surft, & + !IN: JULES dimensions and associated + row_length, rows, land_pts, nsurft, npft, & + surft_pts, ainfo%surft_index, ainfo%land_index, & + !IN: JULES Surface descriptions generally parametrized + ainfo%frac_surft, progs%LAI_pft, progs%canht_pft, & + psparms%albsoil_soilt(:,1), & + !IN: JULES timestep varying fields + psparms%cosz_ij, snow_surft, & + !IN: CABLE specific surface_type indexes + ICE_cable, lakes_cable, & + !IN: CABLE Veg parameters. ICE surface_type index + pars_io_cbl%vegin_xfang, pars_io_cbl%vegin_taul, pars_io_cbl%vegin_refl, & + !IN: CABLE prognostics. decl in progs_cbl_vars_mod.F90 [(:,:,1)= top layer] + progs_cbl%SoilTemp_CABLE(:,:,1), progs_cbl%OneLyrSnowDensity_CABLE, & + progs_cbl%SnowAge_CABLE ) + + +CASE DEFAULT + errorstatus = 101 + WRITE(jules_message,'(A,I0)') 'Unrecognised surface scheme. lsm_id = ', & + lsm_id + CALL ereport(RoutineName, errorstatus, jules_message) + +END SELECT + +CALL jules_ssi_albedo ( & + !INTENT(IN) + !input fields + flandg, ainfo%ice_fract_ij, fluxes%tstar_ij, coast%tstar_sice_sicat, & + psparms%cosz_ij, ws10m, chloro, & + progs%snow_mass_sea_sicat, progs%di_ncat_sicat, & + ainfo%pond_frac_cat_sicat, ainfo%pond_depth_cat_sicat, & + !max and min sea ice albedo specifications + alpham, alphac, alphab, dtice, & + dt_bare, dalb_bare_wet, pen_rad_frac, sw_beta, & + ! parameters for CICE multi-band albedo scheme: + albicev_cice, albicei_cice, albsnowv_cice, albsnowi_cice, & + albpondv_cice, albpondi_cice, & + ahmax, dalb_mlt_cice, dalb_mlts_v_cice, dalb_mlts_i_cice, & + dt_bare_cice, dt_snow_cice, & + pen_rad_frac_cice, sw_beta_cice, snowpatch, & + !size and control variables + row_length * rows, max_n_swbands, & + n_band, nice, nice_use, & + !spectral boundaries + wavelength_short, & + wavelength_long, & + !INTENT(OUT) + !output arguments + sea_ice_albedo, & + open_sea_albedo, & + !ancil_info (IN) + ainfo%sea_index, ainfo%ssi_index, ainfo%sice_index_ncat, & + ainfo%sice_frac_ncat, ainfo%sice_pts_ncat, sea_pts, & + !Fluxes (OUT) + fluxes%alb_sicat, fluxes%penabs_rad_frac) + +IF (lhook) CALL dr_hook(ModuleName//':'//RoutineName,zhook_out,zhook_handle) +RETURN +END SUBROUTINE surf_couple_radiation +END MODULE surf_couple_radiation_mod diff --git a/src/science/radiation/jules_land_albedo_jls_mod.F90 b/src/science/radiation/jules_land_albedo_jls_mod.F90 index 18386fcb..18050bd9 100644 --- a/src/science/radiation/jules_land_albedo_jls_mod.F90 +++ b/src/science/radiation/jules_land_albedo_jls_mod.F90 @@ -27,7 +27,7 @@ SUBROUTINE jules_land_albedo( & lai_pft, canht_pft, & rgrain_surft, snow_surft, tstar_surft, z0_surft, frac_surft, & !INTENT(OUT) - alb_surft,albobs_sc_ij,land_albedo_ij, & + alb_surft,albobs_sc_ij,fsnow,land_albedo_ij, & !New arguments replacing USE statements !jules_mod (IN OUT) albobs_scaling_surft, & @@ -60,10 +60,10 @@ SUBROUTINE jules_land_albedo( & USE jules_snow_mod, ONLY: & kland, maskd, tcland, rho_snow_const, rho_snow_fresh, cansnowtile, & l_snowdep_surf, can_clump, lai_alb_lim_sn, n_lai_exposed, amax, aicemax, & - rho_firn_albedo, nsmax + i_snow_tile, rho_firn_albedo, nsmax USE jules_surface_types_mod, ONLY: & - npft, ntype, lake, soil, urban_canyon + npft, ntype, ice, lake, soil, urban_canyon USE water_constants_mod, ONLY: & rho_ice, tm @@ -168,6 +168,8 @@ SUBROUTINE jules_land_albedo( & albobs_sc_ij(pfield,nsurft,2), & ! albedo scaling to obs in VIS and NIR ! for diagnostics output by the UM + fsnow(land_pts,nsurft), & + ! Snow cover fraction on tiles land_albedo_ij(pfield,4) ! GBM albedos. @@ -258,8 +260,6 @@ SUBROUTINE jules_land_albedo( & alb_snow_surft(land_pts,4), & !Tiled snow albedo: used with embedded !snow - fsnow(land_pts), & - !Weighting factor for albedo. lai(land_pts,npft), & !Adjusted leaf area index. snowd(land_pts), & @@ -1161,7 +1161,7 @@ SUBROUTINE jules_land_albedo( & END DO END IF ! Calculate snow albedo weighting factor. - fsnow(:) = 0.0 + fsnow(:,n) = 0.0 IF ( l_point_data .AND. .NOT. cansnowtile(n)) THEN IF (l_fix_snow_frac) THEN DO j = 1,surft_pts(n) @@ -1170,32 +1170,44 @@ SUBROUTINE jules_land_albedo( & ! Use linear expansion of exponential if non-linear term ! is of order EPSILON (i.e., x^2.0/2.0 ~ EPSILON) IF ( snowd(l) > SQRT(2.0*EPSILON(snowd))/maskd) THEN - fsnow(l) = 1.0 - EXP( -maskd * snowd(l) ) + fsnow(l,n) = 1.0 - EXP( -maskd * snowd(l) ) ELSE - fsnow(l) = maskd * snowd(l) + fsnow(l,n) = maskd * snowd(l) END IF END DO ELSE DO j = 1,surft_pts(n) l = surft_index(j,n) - IF ( snowd(l) > 0.0) fsnow(l) = 1.0 - EXP( -50.0 * snowd(l) ) + IF (snowd(l) > 0.0) fsnow(l,n) = 1.0 - EXP( -maskd * snowd(l) ) END DO END IF ELSE + IF (i_snow_tile(n) == 1) z0(:) = z0_surft(:,ice) DO j = 1,surft_pts(n) l = surft_index(j,n) - IF ( snowd(l) > 0.0) fsnow(l) = snowd(l) / & - ( snowd(l) + 10.0 * z0(l) ) + IF (snowd(l) > 0.0) fsnow(l,n) = snowd(l) / (snowd(l) + 10 * z0(l)) END DO END IF ! Calculate weighted tile albedo. - DO j = 1,surft_pts(n) - l = surft_index(j,n) - DO band = 1,4 - alb_type(l,n,band) = fsnow(l) * alb_snow(l,n,band) & - + (1.0 - fsnow(l)) * alb_type(l,n,band) + IF (i_snow_tile(n) == 1) THEN + ! Ice tile is used as a separate snow tile + ! Preserve snow-free albedo on snow-free fractions of selected tiles + DO j = 1,surft_pts(n) + l = surft_index(j,n) + fsnow(l,ice) = 1.0 + DO band = 1,4 + alb_type(l,ice,band) = alb_snow(l,n,band) + END DO END DO - END DO + ELSE + DO j = 1,surft_pts(n) + l = surft_index(j,n) + DO band = 1,4 + alb_type(l,n,band) = fsnow(l,n) * alb_snow(l,n,band) & + + (1.0 - fsnow(l,n)) * alb_type(l,n,band) + END DO + END DO + END IF END DO ! ntype ELSE ! l_snow_albedo @@ -1232,6 +1244,7 @@ SUBROUTINE jules_land_albedo( & snowd(l) = snowdep_surft(l,n) END DO END IF + fsnow(:,n) = 0.0 DO j = 1,surft_pts(n) l = surft_index(j,n) IF ( tstar(l) < tcland ) THEN @@ -1248,20 +1261,15 @@ SUBROUTINE jules_land_albedo( & ! Use linear expansion of exponential if non-linear term ! is of order EPSILON (i.e., x^2.0/2.0 ~ EPSILON) IF (snowd(l) > SQRT(2.0*EPSILON(snowd))/maskd) THEN - alb_type(l,n,band) = alb_type(l,n,band) + & - (dsa - alb_type(l,n,band)) * & - ( 1.0 - EXP(-maskd * snowd(l)) ) + fsnow(l,n) = 1.0 - EXP(-maskd * snowd(l)) ELSE - alb_type(l,n,band) = alb_type(l,n,band) + & - (dsa - alb_type(l,n,band)) * & - maskd * snowd(l) + fsnow(l,n) = maskd * snowd(l) END IF ELSE - alb_type(l,n,band) = alb_type(l,n,band) + & - (dsa - alb_type(l,n,band)) * & - ( 1.0 - EXP(-maskd * snowd(l)) ) - + fsnow(l,n) = 1.0 - EXP(-maskd * snowd(l)) END IF + alb_type(l,n,band) = alb_type(l,n,band) + & + (dsa - alb_type(l,n,band)) * fsnow(l,n) END DO END DO END DO @@ -1514,16 +1522,14 @@ SUBROUTINE jules_land_albedo( & ! Use linear expansion of exponential if non-linear term ! is of order EPSILON (i.e., x^2.0/2.0 ~ EPSILON) IF (snowd(l) > SQRT(2.0*EPSILON(snowd))/maskd) THEN - alb_type(l,n,1) = albsnf(l,n) + (dsa - albsnf(l,n)) * & - ( 1.0 - EXP(-maskd * snowd(l)) ) + fsnow(l,n) = 1.0 - EXP(-maskd * snowd(l)) ELSE - alb_type(l,n,1) = albsnf(l,n) + (dsa - albsnf(l,n)) * & - maskd * snowd(l) + fsnow(l,n) = maskd * snowd(l) END IF ELSE - alb_type(l,n,1) = albsnf(l,n) + (dsa - albsnf(l,n)) * & - ( 1.0 - EXP(-maskd * snowd(l)) ) + fsnow(l,n) = 1.0 - EXP(-maskd * snowd(l)) END IF + alb_type(l,n,1) = albsnf(l,n) + (dsa - albsnf(l,n)) * fsnow(l,n) END DO END DO @@ -1546,14 +1552,33 @@ SUBROUTINE jules_land_albedo( & DO i = 1,pfield land_albedo_ij(i,band) = 0.0 END DO - DO n = 1,ntype - DO j = 1,surft_pts(n) - l = surft_index(j,n) - i = land_index(l) - land_albedo_ij(i,band) = land_albedo_ij(i,band) + & - frac_surft(l,n) * alb_type(l,n,band) + IF (ANY(i_snow_tile == 1)) THEN + ! using a separate snow tile + DO n = 1,ntype + DO j = 1,surft_pts(n) + l = surft_index(j,n) + i = land_index(l) + IF (i_snow_tile(n) == 1) THEN + land_albedo_ij(i,band) = land_albedo_ij(i,band) + & + (1 - fsnow(l,n)) * frac_surft(l,n) * alb_type(l,n,band) & + + fsnow(l,n) * frac_surft(l,n) * alb_type(l,ice,band) + ELSE + land_albedo_ij(i,band) = land_albedo_ij(i,band) + & + frac_surft(l,n) * alb_type(l,n,band) + END IF + END DO END DO - END DO + ELSE + ! not using a separate snow tile + DO n = 1,ntype + DO j = 1,surft_pts(n) + l = surft_index(j,n) + i = land_index(l) + land_albedo_ij(i,band) = land_albedo_ij(i,band) + & + frac_surft(l,n) * alb_type(l,n,band) + END DO + END DO + END IF END DO !---------------------------------------------------------------------- diff --git a/src/science/snow/snow_mod.F90 b/src/science/snow/snow_mod.F90 index f76818e9..a373a92b 100644 --- a/src/science/snow/snow_mod.F90 +++ b/src/science/snow/snow_mod.F90 @@ -24,8 +24,8 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & hcaps1_soilt, hcons, melt_surft, snowinc_surft, & smcl1_soilt, sthf1_soilt, surf_htf_surft, & t_soil1_soilt, tsurf_elev_surft, tstar_surft, & - smvcst1_soilt, con_snow_wtrac, ei_surft_wtrac, & - rgrain, rgrainl, rho_snow_grnd, sice, & + z0_surft, smvcst1_soilt, con_snow_wtrac, ei_surft_wtrac, & + fsnow, rgrain, rgrainl, rho_snow_grnd, sice, & sliq, snow_grnd, snow_surft, snowdepth, & tsnow, nsnow, con_rain_wtrac, ls_rain_wtrac, & ls_snow_wtrac, ls_graup_wtrac, melt_surft_wtrac, & @@ -71,12 +71,20 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & ! Include graupel in the surface snowfall. ignore_graupel, & ! Ignore graupel in the surface snowfall. + i_snow_tile, & + ! Tiles without (0) or with (1) separate energy balance for snow. + maskd, & + ! Inverse masking depth for snow albedo (1/m). + rho_snow_fresh, & + ! Density of fresh snow (kg per m**3). r0, & ! Grain size for fresh snow (microns). cansnowtile ! Switch for canopy snow model. -USE jules_surface_types_mod, ONLY: lake, ntype +USE jules_surface_types_mod, ONLY: ice, lake, ntype + +USE jules_surface_mod, ONLY: l_point_data USE jules_radiation_mod, ONLY: l_snow_albedo, l_embedded_snow @@ -128,7 +136,7 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & con_snow(land_pts), & ! Convective snowfall rate (kg/m2/s). tile_frac(land_pts,nsurft), & - ! Tile fractions. + ! Tile fractions including snow cover in the ice tile. ei_surft(land_pts,nsurft), & ! Sublimation of snow (kg/m2/s). snowinc_surft(land_pts,nsurft), & @@ -157,6 +165,8 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & ! Surface heat flux (W/m2). tstar_surft(land_pts,nsurft), & ! Tile surface temperature (K). + z0_surft(land_pts,nsurft), & + ! Tile roughness lengths (m). smvcst1_soilt(land_pts,nsoilt), & ! Surface soil layer volumetric moisture concentration at saturation. con_snow_wtrac(land_pts,n_wtrac_jls), & @@ -181,6 +191,8 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & ! Large-scale graupel fall rate (kg/m2/s). ls_rain(land_pts), & ! Large-scale rainfall rate (kg/m2/s). + fsnow(land_pts,nsurft), & + ! Snow cover fractions on tiles. melt_surft(land_pts,nsurft), & ! Surface or canopy snowmelt rate (kg/m2/s). ! On output, this is the total melt rate for the tile @@ -315,6 +327,8 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & snowfall(land_pts), & ! Total frozen precip reaching the ground in timestep ! (kg/m2) - includes any canopy unloading. + snowonsnow(land_pts), & + ! Frozen precipitation on the snow tile (kg/m2). graupfall(land_pts), & ! Graupel reaching the ground in timestep(kg/m2). infiltration(land_pts), & @@ -359,6 +373,17 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & rgrainl_sl(land_pts,nsmax) ! Snow layer grain size (microns). +! Snow quantities on a separate snow tile +REAL(KIND=real_jlslsm) :: & + fsnow_old(land_pts), & + ! Snow cover fraction at start of timestep. + fsnow_new(land_pts), & + ! Snow cover fraction at end of timestep. + fsratio(land_pts), & + ! Ratio of snow cover fractions at start and end of timestep. + snowdepave(land_pts) + ! Average snow depth on surface types with a separate snow tile. + REAL(KIND=real_jlslsm), ALLOCATABLE :: & snow_surft_old(:,:), & snow_grnd_old(:,:), & @@ -805,6 +830,35 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & !$OMP END PARALLEL END IF + !--------------------------------------------------------------------------- + ! Scale snowfall on selected tiles to the snow tile (number nsurft) + !--------------------------------------------------------------------------- + IF (i_snow_tile(n) == 1) THEN + DO k = 1,surft_pts(n) + i = surft_index(k,n) + snowonsnow(i) = 0.0 + IF (fsnow(i,n) > EPSILON(0.0)) THEN + snowonsnow(i) = snowfall(i) / fsnow(i,n) + ELSE IF (snowfall(i) > 0.0) THEN + snowdepth(i,n) = snowfall(i) / rho_snow_fresh + IF (l_point_data) THEN + fsnow(i,n) = 1 - EXP(- maskd * snowdepth(i,n)) + ELSE + fsnow(i,n) = snowdepth(i,n) / & + (snowdepth(i,n) + 10*z0_surft(i,nsurft)) + END IF + snowonsnow(i) = snowfall(i) / fsnow(i,n) + END IF + fsnow_old(i) = fsnow(i,n) + snowfall(i) = 0.0 + END DO + END IF + IF ( ANY(i_snow_tile == 1) .AND. n == nsurft ) THEN + DO i = 1,land_pts + IF (.NOT. l_lice_point(i)) snowfall(i) = snowonsnow(i) + END DO + END IF + !========================================================================== ! *NOTICE REGARDING SOIL TILING** ! @@ -903,6 +957,35 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & surft_index(:,n), sice_sl, sliq_sl, & tsnow_sl, rho_snow_sl, ds_sl ) + !------------------------------------------------------------------------- + ! Update snow cover fraction and layers on the separate snow tile + !------------------------------------------------------------------------- + IF (ANY(i_snow_tile == 1) .AND. n == nsurft) THEN + DO k = 1,surft_pts(n) + i = surft_index(k,n) + IF (.NOT. l_lice_point(i)) THEN + snowdepth(i,n) = 0.0 + IF (nsnow(i,n) > 0) snowdepth(i,n) = SUM(ds_sl(i,1:nsnow(i,n))) + snowdepth(i,n) = snowdepth(i,n) + sice0(i) / rho0(i) + snowdepave(i) = fsnow_old(i) * snowdepth(i,n) + IF (l_point_data) THEN + fsnow_new(i) = 1 - EXP(- maskd * snowdepave(i)) + ELSE + fsnow_new(i) = snowdepave(i) / (snowdepave(i) + 10*z0_surft(i,n)) + END IF + fsratio(i) = 0.0 + IF (fsnow_new(i) > EPSILON(0.0)) & + fsratio(i) = fsnow_old(i) / fsnow_new(i) + sice0(i) = fsratio(i) * sice0(i) + snowfall(i) = fsratio(i) * snowfall(i) + snowmass(i) = fsratio(i) * snowmass(i) + ds_sl(i,:) = fsratio(i) * ds_sl(i,:) + sice_sl(i,:) = fsratio(i) * sice_sl(i,:) + sliq_sl(i,:) = fsratio(i) * sliq_sl(i,:) + END IF + END DO + END IF + !------------------------------------------------------------------------- ! Redivide snowpack after changes in depth, conserving mass and energy !------------------------------------------------------------------------- @@ -1058,6 +1141,26 @@ SUBROUTINE snow ( land_pts, timestep, stf_hf_snow_melt, nsurft, n_wtrac_jls, & END DO ! tiles +!----------------------------------------------------------------------- +! Store average snow masses on selected tiles for diagnostics +! Snow depth is required for snow fraction calculations +!----------------------------------------------------------------------- +IF (ANY(i_snow_tile == 1)) THEN + DO n = 1,nsurft-1 + IF (i_snow_tile(n) == 1) THEN + DO k = 1,surft_pts(n) + i = surft_index(k,n) + fsnow(i,n) = fsnow_new(i) + sice(i,n,:) = fsnow_new(i) * sice(i,nsurft,:) + sliq(i,n,:) = fsnow_new(i) * sliq(i,nsurft,:) + snow_surft(i,n) = fsnow_new(i) * snow_surft(i,nsurft) + rgrain(i,n) = rgrain(i,nsurft) + snowdepth(i,n) = snowdepave(i) + END DO + END IF + END DO +END IF + IF (l_flake_model) THEN !$OMP PARALLEL DEFAULT(NONE) & !$OMP PRIVATE(i,i_wt) & diff --git a/src/science/surface/jules_land_sf_explicit_jls.F90 b/src/science/surface/jules_land_sf_explicit_jls.F90 index f722d839..4d498753 100644 --- a/src/science/surface/jules_land_sf_explicit_jls.F90 +++ b/src/science/surface/jules_land_sf_explicit_jls.F90 @@ -39,8 +39,8 @@ SUBROUTINE jules_land_sf_explicit ( & bq_1,bt_1,z1_uv,z1_uv_top,z1_tq,z1_tq_top,qw_1,tl_1, & ! IN soil/vegetation/land surface data : land_index,nsurft,sm_levels,canopy,catch,catch_snow,hcon_soilt, & - ho2r2_orog, flandg, & - snow_surft,sil_orog_land,smvccl_soilt,smvcst_soilt,smvcwt_soilt,sthf_soilt, & + ho2r2_orog, flandg, fsnow, & + sil_orog_land,smvccl_soilt,smvcst_soilt,smvcwt_soilt,sthf_soilt, & sthu_soilt,z0_surft,z0h_surft_bare, z0m_soil_in, & ! IN input data from the wave model charnock_w, & @@ -99,9 +99,10 @@ SUBROUTINE jules_land_sf_explicit ( & wrr_gb, & !Fluxes (IN OUT) anthrop_heat_surft, & - !prognostics (IN) - nsnow_surft, sice_surft, sliq_surft, snowdepth_surft, & - tsnow_surft, ds_surft, & + !snow prognostics (IN) + nsnow_surft, sice_surft, sliq_surft, tsnow_surft, ds_surft, & + !snow prognostics (IN OUT) + snowdepth_surft, snow_surft, & !c_elevate (OUT) surf_hgt_surft, lw_down_elevcorr_surft, & !jules_mod (OUT) @@ -169,8 +170,8 @@ SUBROUTINE jules_land_sf_explicit ( & USE jules_soil_mod, ONLY: dzsoil, dzsoil_elev, hcice, hcwat, hcondeep -USE jules_surface_types_mod, ONLY: npft, nnpft, ntype, & - urban_canyon, urban_roof, soil, lake, ncpft +USE jules_surface_types_mod, ONLY: npft, ncpft, nnpft, ntype, & + urban_canyon, urban_roof, soil, lake, ice #if defined(UM_JULES) USE atm_step_local, ONLY: dim_cs1, co2_dim_len,co2_dim_row @@ -181,6 +182,7 @@ SUBROUTINE jules_land_sf_explicit ( & USE jules_snow_mod, ONLY: cansnowtile & ,rho_snow_const & ,snow_hcon & + ,i_snow_tile & ,l_snowdep_surf & ,l_snow_nocan_hc & ,nsmax & @@ -307,11 +309,11 @@ SUBROUTINE jules_land_sf_explicit ( & ,catch_snow(land_pts,nsurft) & ! IN Snow interception capacity of ! tiles (kg/m2). +,fsnow(land_pts,nsurft) & + ! IN Snow cover fractions on tiles ,hcon_soilt(land_pts,nsoilt) & ! IN Soil thermal conductivity ! (W/m/K). -,snow_surft(land_pts,nsurft) & - ! IN Lying snow on tiles (kg/m2) ,smvccl_soilt(land_pts,nsoilt,sm_levels) & ! IN Critical volumetric SMC ! (cubic m per cubic m of soil). @@ -707,14 +709,17 @@ SUBROUTINE jules_land_sf_explicit ( & wt_ext_irr_surft(land_pts,sm_levels,nsurft) REAL(KIND=real_jlslsm), INTENT(OUT) :: gc_irr_surft(land_pts,nsurft) -!prognostics (IN) +!snow prognostics (IN) INTEGER, INTENT(IN) :: nsnow_surft(land_pts,nsurft) REAL(KIND=real_jlslsm), INTENT(IN) :: sice_surft(land_pts,nsurft,nsmax), & sliq_surft(land_pts,nsurft,nsmax), & - snowdepth_surft(land_pts,nsurft), & tsnow_surft(land_pts,nsurft,nsmax), & ds_surft(land_pts,nsurft,nsmax) +!snow prognostics (IN OUT) +REAL(KIND=real_jlslsm), INTENT(IN OUT) :: snowdepth_surft(land_pts,nsurft), & + snow_surft(land_pts,nsurft) + !c_elevate (OUT) REAL(KIND=real_jlslsm), INTENT(OUT) :: surf_hgt_surft(land_pts,nsurft), & lw_down_elevcorr_surft(land_pts,nsurft) @@ -1433,6 +1438,28 @@ SUBROUTINE jules_land_sf_explicit ( & END DO END IF +!----------------------------------------------------------------------- +! Modify tile fractions and zero snow mass and depth on selected tiles +! if using the ice tile as a separate snow tile. +!----------------------------------------------------------------------- +IF (ANY(i_snow_tile == 1)) THEN + DO l = 1, land_pts + IF (.NOT. l_lice_point(l)) tile_frac(l,ice) = 0.0 + END DO + DO n = 1,nsurft-1 + DO l = 1, land_pts + tile_frac(l,n) = frac(l,n) + IF (i_snow_tile(n) == 1) THEN + snow_surft(l,n) = 0.0 + snowdep_surft(l,n) = 0.0 + snowdepth_surft(l,n) = 0.0 + tile_frac(l,n) = (1 - fsnow(l,n))*frac(l,n) + tile_frac(l,ice) = tile_frac(l,ice) + fsnow(l,n)*frac(l,n) + END IF + END DO + END DO +END IF + IF (land_pts > 0) THEN ! Omit if no land points !----------------------------------------------------------------------- diff --git a/src/science/surface/jules_land_sf_implicit.jls.F90 b/src/science/surface/jules_land_sf_implicit.jls.F90 index effb1549..167fde7f 100644 --- a/src/science/surface/jules_land_sf_implicit.jls.F90 +++ b/src/science/surface/jules_land_sf_implicit.jls.F90 @@ -38,7 +38,7 @@ SUBROUTINE jules_land_sf_implicit ( & ! IN everything not covered so far : lw_down,sw_surft,sky,t_soil_soilt,r_gamma,alpha1,ashtf_prime_surft, & dtrdz_charney_grid_1,fracaero_t,fracaero_s,resfs,resft,rhokh_surft, & - emis_surft,snow_surft,dtstar_surft, & + emis_surft,fsnow,snow_surft,dtstar_surft, & ! INOUT data : tstar_surft,fqw_surft,fqw_1,ftl_1,ftl_surft,sf_diag, & ! OUT Diagnostic not requiring STASH flags : @@ -88,7 +88,7 @@ SUBROUTINE jules_land_sf_implicit ( & USE jules_surface_mod, ONLY: l_aggregate, l_flake_model, ls USE jules_snow_mod, ONLY: & - nsmax, rho_snow_const, cansnowtile, l_snow_nocan_hc + nsmax, rho_snow_const, cansnowtile, i_snow_tile, l_snow_nocan_hc USE jules_surface_types_mod, ONLY: lake @@ -171,6 +171,8 @@ SUBROUTINE jules_land_sf_implicit ( & ! IN Land fraction on all pts. ,emis_surft(land_pts,nsurft) & ! IN Emissivity for land tiles +,fsnow(land_pts,nsurft) & + ! IN Snow cover fractions on tiles ,snow_surft(land_pts,nsurft) & ! IN Lying snow on tiles (kg/m2) ,dtstar_surft(land_pts,nsurft) @@ -567,7 +569,7 @@ SUBROUTINE jules_land_sf_implicit ( & land_pts,nsurft, & land_index,surft_index,surft_pts,sm_levels,fland, & ashtf_prime_surft,canopy,dtrdz_charney_grid_1,flake,fracaero_t,fracaero_s, & - snow_surft,resfs,resft,rhokh_surft,tile_frac,smc_soilt,wt_ext_surft, & + fsnow,snow_surft,resfs,resft,rhokh_surft,tile_frac,smc_soilt,wt_ext_surft, & timestep,r_gamma,fqw_1,fqw_surft,ftl_1,ftl_surft,tstar_surft, & ecan,ecan_surft,elake_surft,esoil_soilt,esoil_surft,ei_surft,ext_soilt, & sf_diag, non_lake_frac, & @@ -798,6 +800,19 @@ SUBROUTINE jules_land_sf_implicit ( & lw_down_elevcorr_surft(l,n) END DO !$OMP END DO + ! Tiles with partial snow cover + IF (i_snow_tile(n) == 1) THEN + DO k = 1,surft_pts(n) + l = surft_index(k,n) + j = (land_index(l) - 1) / tdims%i_end + 1 + i = land_index(l) - (j-1) * tdims%i_end + sf_diag%lw_up_surft(l,n) = & + (1.0 - fsnow(i,n))*sf_diag%lw_up_surft(l,n) + fsnow(i,n) * & + (emis_surft(l,nsurft) * sbcon * tstar_surft(l,nsurft)**4 & + + (1.0 - emis_surft(l,nsurft)) * & + (lw_down(i,j) + lw_down_elevcorr_surft(l,nsurft))) + END DO + END IF END DO END IF END IF diff --git a/src/science/surface/sf_evap_jls.F90 b/src/science/surface/sf_evap_jls.F90 index b5dfaeff..8d69783e 100644 --- a/src/science/surface/sf_evap_jls.F90 +++ b/src/science/surface/sf_evap_jls.F90 @@ -25,8 +25,8 @@ MODULE sf_evap_mod SUBROUTINE sf_evap ( & land_pts,nsurft, & land_index,surft_index,surft_pts,sm_levels,fland, & - ashtf_prime_surft,canopy,dtrdz_1,flake,fracaero_t,fracaero_s,snow_surft,resfs,& - resft,rhokh_1,tile_frac,smc_soilt,wt_ext_surft,timestep,r_gamma, & + ashtf_prime_surft,canopy,dtrdz_1,flake,fracaero_t,fracaero_s,fsnow,snow_surft,& + resfs,resft,rhokh_1,tile_frac,smc_soilt,wt_ext_surft,timestep,r_gamma, & fqw_1,fqw_surft,ftl_1,ftl_surft,tstar_surft, & ecan,ecan_surft,elake_surft,esoil_soilt,esoil_surft,ei_surft,ext_soilt, & sf_diag, non_lake_frac, & @@ -48,6 +48,8 @@ SUBROUTINE sf_evap ( & USE jules_irrig_mod, ONLY: l_irrig_dmd +USE jules_snow_mod, ONLY: i_snow_tile + USE jules_surface_mod, ONLY: l_flake_model USE jules_surface_types_mod, ONLY: lake @@ -94,6 +96,8 @@ SUBROUTINE sf_evap ( & ! ! with only aerodynamic resistance ! ! for land tiles from the frozen part of the ! ! surface alone +,fsnow(land_pts,nsurft) & + ! IN Snow cover fractions on tiles ,snow_surft(land_pts,nsurft) & ! ! IN Lying snow amount on tiles (kg/m2). ,resfs(land_pts,nsurft) & @@ -547,6 +551,22 @@ SUBROUTINE sf_evap ( & !$OMP END DO END DO !nsurft +!----------------------------------------------------------------------- +! Do not deposit frost on the snow-free fractions of tiles with a +! separate snow tile +!----------------------------------------------------------------------- +IF (ANY(i_snow_tile == 1)) THEN + DO n = 1,nsurft-1 + DO k = 1,surft_pts(n) + l = surft_index(k,n) + IF (ei_surft(l,n) < 0.0) THEN + fqw_surft(l,n) = ei_surft(l,n) + ei_surft(l,n) = 0.0 + END IF + END DO + END DO +END IF + IF ((nsoilt == 1) .AND. (l_flake_model)) THEN !$OMP DO SCHEDULE(STATIC) DO l = 1,land_pts diff --git a/src/science/vegetation/sparm_jls_mod.F90 b/src/science/vegetation/sparm_jls_mod.F90 index 5e71ff4c..75aa1a2e 100644 --- a/src/science/vegetation/sparm_jls_mod.F90 +++ b/src/science/vegetation/sparm_jls_mod.F90 @@ -31,12 +31,12 @@ SUBROUTINE sparm (land_pts, nsurft, surft_pts, surft_index, & USE nvg_sparm_mod, ONLY: nvg_sparm !Use in relevant variables -USE jules_surface_types_mod, ONLY: lake, npft, ntype, soil, urban_canyon, & - urban_roof +USE jules_surface_types_mod, ONLY: ice, lake, npft, ntype, soil, & + urban_canyon, urban_roof USE jules_vegetation_mod, ONLY: can_model, l_spec_veg_z0 USE blend_h, ONLY: lb USE jules_surface_mod, ONLY: i_aggregate_opt, l_aggregate -USE jules_snow_mod, ONLY: cansnowtile, snowloadlai +USE jules_snow_mod, ONLY: cansnowtile, i_snow_tile, snowloadlai USE c_z0h_z0m, ONLY: z0h_z0m USE jules_urban_mod, ONLY: l_moruses @@ -102,6 +102,8 @@ SUBROUTINE sparm (land_pts, nsurft, surft_pts, surft_index, & ! GBM canopy capacity (kg/m2). catch_t(land_pts,ntype), & ! Canopy capacities for types not tiles (kg/m2); required for l_aggregate. + fst(land_pts), & + ! Total fraction of tiles using a separate energy balance for snow. fz0(land_pts), & ! Aggregation function of Z0. fz0h(land_pts), & @@ -165,6 +167,34 @@ SUBROUTINE sparm (land_pts, nsurft, surft_pts, surft_index, & CALL nvg_sparm (land_pts, surft_pts, surft_index, z0m_soil_gb, ztm_gb, & catch_surft, z0_surft) + !----------------------------------------------------------------------------- + ! Set parameters for the snow tile, number ice = ntype + !----------------------------------------------------------------------------- + ! A separate snow tile can only be used if this is not a land ice point and + ! surface types are not aggregated + IF (ANY(i_snow_tile == 1)) THEN + fst(:) = 0 + fz0(:) = 0 + DO n = 1,ntype-1 + IF (i_snow_tile(n) == 1) THEN + DO j = 1,surft_pts(n) + l = surft_index(j,n) + IF (frac_surft(l,ice) == 0.0) THEN + fst(l) = fst(l) + frac_surft(l,n) + fz0(l) = fz0(l) + frac_surft(l,n) / (LOG(lb / z0_surft(l,n)))**2 + END IF + END DO + END IF + END DO + DO l = 1,land_pts + IF (frac_surft(l,ice) == 0.0) THEN + catch_surft(l,ice) = 0.0 + IF (fz0(l) > EPSILON(0.0)) & + z0_surft(l,ice) = lb * EXP(-SQRT(fst(l)/fz0(l))) + END IF + END DO + END IF + ELSE ! l_aggregate !--------------------------------------------------------------------------- ! Set parameters for vegetated surface types