Uvic ice jp

src/airsea/airsea.F90

@@ -65,6 +65,7 @@ module airsea_driver
 !
 !  Meteorological forcing variables
    integer,  public                    :: hum_method
+   integer,  public                    :: longwave_type
    character(len=PATH_MAX)   :: meteo_file
    type (type_scalar_input), public, target            :: u10,v10
    type (type_scalar_input), public, target            :: airp,airt
@@ -76,10 +77,10 @@ module airsea_driver
 !
 !  surface shortwave radiation
 !  and surface heat flux (W/m^2)
-   type (type_scalar_input), public, target :: I_0, ql
+   type (type_scalar_input), public, target :: I_0
    REALTYPE, public, target            :: albedo
-   type (type_scalar_input), public, target            :: heat
-   REALTYPE, public                    :: qe,qh
+   type (type_scalar_input), public, target :: heat, ql_
+   REALTYPE, public                    :: qe,qh,ql
 
 !  surface stress components (Pa)
    REALTYPE, public, target                 :: tx,ty
@@ -341,7 +342,7 @@ subroutine init_airsea_nml(namlst, fn)
    call tx_%configure(method=momentum_method, path=momentumflux_file, index=1, constant_value=const_tx)
    call ty_%configure(method=momentum_method, path=momentumflux_file, index=2, constant_value=const_ty)
    call heat%configure(method=heat_method, path=heatflux_file, index=1, scale_factor=shf_factor, constant_value=const_heat)
-   call ql%configure(method=back_radiation_method, path=back_radiation_file, index=1)
+   call ql_%configure(method=back_radiation_method, path=back_radiation_file, index=1)
    call sst_obs%configure(method=sst_method, path=sst_file, index=1)
    call sss%configure(method=sss_method, path=sss_file, index=1)
    call precip%configure(method=precip_method, path=precip_file, index=1, scale_factor=precip_factor, constant_value=const_precip)
@@ -491,12 +492,24 @@ subroutine init_airsea_yaml()
    call branch%get(calc_evaporation, 'calc_evaporation', 'calculate evaporation from meteorological conditions', &
                 default=.false.)
 
+ !jpnote: before cherry pick 
    call branch%get(I_0, 'swr', 'shortwave radiation', 'W/m^2', &
                 minimum=0._rk,default=0._rk, extra_options=(/option(3, 'from time, location and cloud cover', 'calculate')/))
-   call branch%get(ql, 'longwave_radiation', 'net longwave radiation', 'W/m^2', &
+   !call branch%get(ql, 'longwave_radiation', 'net longwave radiation', 'W/m^2', &
+               ! default=0._rk, method_file=0, method_constant=method_unsupported, &
+               !extra_options=(/option(CLARK, 'Clark et al. (1974)', 'Clark'), option(HASTENRATH_LAMB, 'Hastenrath and Lamb (1978)', 'Hastenrath_Lamb'), option(BIGNAMI, 'Bignami et al. (1995)', 'Bignami'), option(BERLIAND_BERLIAND, 'Berliand and Berliand (1952)', 'Berliand_Berliand'), option(JOSEY1, 'Josey et al. (2003) - 1', 'Josey1'), option(JOSEY2, 'Josey et al. (2003) - 2', 'Josey2')/), default_method=CLARK)
+!jpnote: after cherry pick 1
+   !call branch%get(ql, 'longwave_radiation', 'longwave radiation', 'W/m^2', &
+   !call twig%get(ssuv_method, 'ssuv_method', 'wind treatment', &
+                !options=(/option(0, 'use absolute wind speed'), option(1, 'use wind speed relative to current velocity')/), default=1, display=display_advanced)
+!jpnote: after cherry pick 2
+   call branch%get(ql_, 'longwave_radiation', 'longwave radiation', 'W/m^2', &
                 default=0._rk, method_file=0, method_constant=method_unsupported, &
-               extra_options=(/option(CLARK, 'Clark et al. (1974)', 'Clark'), option(HASTENRATH_LAMB, 'Hastenrath and Lamb (1978)', 'Hastenrath_Lamb'), option(BIGNAMI, 'Bignami et al. (1995)', 'Bignami'), option(BERLIAND_BERLIAND, 'Berliand and Berliand (1952)', 'Berliand_Berliand'), option(JOSEY1, 'Josey et al. (2003) - 1', 'Josey1'), option(JOSEY2, 'Josey et al. (2003) - 2', 'Josey2')/), default_method=CLARK)
+               extra_options=(/option(1, 'Clark'), option(2, 'Hastenrath'), option(3, 'Bignami'), option(4, 'Berliand'), option(5, 'Josey-1'), option(6, 'Josey-2')/), default_method=1, pchild=leaf)
+   call leaf%get(longwave_type, 'type', 'longwave type from file', &
+                 options=(/option(1, 'net longwave radiation'), option(2, 'incoming longwave radiation')/), default=1)
 
+
    twig => branch%get_typed_child('albedo')
    call twig%get(albedo_method, 'method', 'method to compute albedo', &
                 options=(/option(0, 'constant', 'constant'), option(PAYNE, 'Payne (1972)', 'Payne'), option(COGLEY, 'Cogley (1979)', 'Cogley')/), default=PAYNE)
@@ -687,11 +700,11 @@ subroutine post_init_airsea(lat,lon)
             LEVEL4 'using Fairall et. all formulation'
          case default
       end select
-      LEVEL3 'net longwave radiation:'
-      select case (ql%method)
+      LEVEL3 'longwave radiation:'
+      select case (ql_%method)
          case(0) ! Read from file instead of calculating
-            call register_input(ql)
-         case(CLARK)
+            call register_input(ql_) !jpnote commmit 2 change
+         case(1)
             LEVEL4 'using Clark formulation'
          case(HASTENRATH_LAMB)
             LEVEL4 'using Hastenrath formulation'
@@ -762,9 +775,11 @@ subroutine surface_fluxes_uvic(surface_temp,sensible,latent,longwave_radiation_v
    ! call humidity(hum_method,rh,airp,TTss-kelvin,airt) 
       call humidity(hum_method,hum%value,airp%value,tw,airt%value)
    ! call longwave_radiation(longwave_radiation_method, &
-                     !    lat,TTss,airt+kelvin,cloud,qb)    
-      call longwave_radiation(ql%method, &
-                        dlat,tw_k,ta_k,cloud%value,longwave_radiation_value)  
+                     !    lat,TTss,airt+kelvin,cloud,qb)   
+      call longwave_radiation(ql_%method,longwave_type, &
+                        dlat,tw_k,ta_k,cloud%value,ql_%value,longwave_radiation_value) 
+      !call longwave_radiation(ql%method, &
+                       ! dlat,tw_k,ta_k,cloud%value,longwave_radiation_value)  
       !call airsea_fluxes(fluxes_method, &
                      !TTss-kelvin,airt,u10,v10,precip,evap,tx,ty,qe,qh) 
       call airsea_fluxes(fluxes_method, &
@@ -792,10 +807,10 @@ subroutine surface_fluxes(surface_temp,sensible,latent,longwave_radiation)
    latent = qe
 #if 0
    if (qe .lt. _ZERO_) then
-      STDERR 'Stefan# ',qh/qe
+!KB      STDERR 'Stefan# ',qh/qe
    end if
 #endif
-   longwave_radiation = ql%value
+   longwave_radiation = ql
    return
    end subroutine surface_fluxes
 !EOC
@@ -1017,9 +1032,8 @@ subroutine flux_from_meteo(jul,secs)
       cloud1 = cloud2
 
       call humidity(hum_method,hum,airp,tw,ta)
-      if (ql%method .gt. 0) then
-         call longwave_radiation(ql%method, &
-                                 dlat,tw_k,ta_k,cloud,ql)
+      call longwave_radiation(ql_%method,longwave_type, &
+                              dlat,tw_k,ta_k,cloud,ql_%value,ql)
       end if
 #if 0
       call airsea_fluxes(fluxes_method,rain_impact,calc_evaporation, &
@@ -1028,7 +1042,7 @@ subroutine flux_from_meteo(jul,secs)
       call airsea_fluxes(fluxes_method, &
                          tw,ta,u10%value-ssu,v10%value-ssv,precip%value,evap,tx2,ty2,qe,qh)
 #endif
-      h2     = ql%value+qe+qh
+      h2     = ql+qe+qh
       cloud2 = cloud%value
 
       if (init_saved_vars) then
@@ -1069,13 +1083,11 @@ subroutine flux_from_meteo(jul,secs)
    end if
 
    call humidity(hum_method,hum%value,airp%value,tw,ta)
-   if (ql%method .gt. 0) then
-      call longwave_radiation(ql%method, &
-                          dlat,tw_k,ta_k,cloud%value,ql%value)
-   endif
+   call longwave_radiation(ql_%method,longwave_type, &
+                           dlat,tw_k,ta_k,cloud%value,ql_%value,ql)
    call airsea_fluxes(fluxes_method, &
                       tw,ta,u10%value-ssu,v10%value-ssv,precip%value,evap,tx_%value,ty_%value,qe,qh)
-   heat%value = (ql%value+qe+qh)
+   heat%value = (ql+qe+qh)
 #endif
 
    w = sqrt((u10%value-ssu)*(u10%value-ssu)+(v10%value-ssv)*(v10%value-ssv))

src/airsea/airsea_variables.F90

@@ -23,7 +23,8 @@ module airsea_variables
    REALTYPE, public, parameter         :: emiss=0.97
    REALTYPE, public, parameter         :: bolz=5.67e-8
    REALTYPE, public, parameter         :: kelvin=273.15
-   REALTYPE, public, parameter         :: const06=0.62198
+   REALTYPE, public, parameter         :: const06=0.62198 ! molecular weight ratio between water and dry air
+                                                          ! 18.01528 g/mol H2O, 28.97 g/mol dry air
    REALTYPE, public, parameter         :: rgas = 287.1    !
    REALTYPE, public, parameter         :: g = 9.81        ! [m/s2]
    REALTYPE, public, parameter         :: rho_0 = 1025.   ! [kg/m3]
@@ -46,6 +47,7 @@ module airsea_variables
    integer, public,  parameter         :: COGLEY=2
 
    ! Longwave radiation
+   integer, public,  parameter         :: FILE = 0                ! From file
    integer, public,  parameter         :: CLARK = 1               ! Clark et al, 1974
    integer, public,  parameter         :: HASTENRATH_LAMB = 2     ! Hastenrath and Lamb, 1978
    integer, public,  parameter         :: BIGNAMI = 3             ! Bignami et al., 1995 - Medsea

src/airsea/longwave_radiation.F90

@@ -5,13 +5,14 @@
 ! !ROUTINE: Calculate the net longwave radiation \label{sec:back-rad}
 !
 ! !INTERFACE:
-   subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
+   subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,qlobs,ql)
 !
 ! !DESCRIPTION:
 !
 ! Here, the net longwave radiation is calculated by means of one out
 ! of six methods, which depend on the value given to the parameter
 ! {\tt method}:
+! {\tt method}=0: read observations from a file
 ! {\tt method}=1: \cite{Clarketal74},
 ! {\tt method}=2: \cite{HastenrathLamb78},
 ! {\tt method}=3: \cite{Bignamietal95},
@@ -23,20 +24,29 @@ subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
 ! !USES:
    use airsea_variables, only: emiss,bolz
    use airsea_variables, only: ea,qa
-   use airsea_variables, only: CLARK, HASTENRATH_LAMB, BIGNAMI, BERLIAND_BERLIAND, JOSEY1, JOSEY2
+   use airsea_variables, only: FILE, CLARK, HASTENRATH_LAMB, BIGNAMI, BERLIAND_BERLIAND, JOSEY1, JOSEY2
    IMPLICIT NONE
 !
 ! !INPUT PARAMETERS:
-   integer, intent(in)                 :: method
-   REALTYPE, intent(in)                :: dlat,tw,ta,cloud
+   integer, intent(in)                 :: method,type
+   REALTYPE, intent(in)                :: dlat,tw,ta,cloud,qlobs
 !
 ! !OUTPUT PARAMETERS:
-   REALTYPE, intent(out)               :: ql
+   REALTYPE, intent(inout)               :: ql
 !
 ! !REVISION HISTORY:
 !  Original author(s): Adolf Stips, Hans Burchard & Karsten Bolding
 !
 ! !LOCAL VARIABLES:
+!cherry pick -- > jpnote to get rid of and add to airseavariables
+   !integer, parameter   :: from_file=0
+   !integer, parameter   :: clark=1      ! Clark et al, 1974
+   !integer, parameter   :: hastenrath=2 ! Hastenrath and Lamb, 1978
+   !integer, parameter   :: bignami=3    ! Bignami et al., 1995 - Medsea
+   !integer, parameter   :: berliand=4   ! Berliand and Berliand, 1952 - ROMS
+   !integer, parameter   :: josey1=5     ! Josey 2003, (J1,9)
+   !integer, parameter   :: josey2=6     ! Josey 2003, (J2,14)
+
    REALTYPE, parameter, dimension(91)  :: cloud_correction_factor = (/ &
      0.497202,     0.501885,     0.506568,     0.511250,     0.515933, &
      0.520616,     0.525299,     0.529982,     0.534665,     0.539348, &
@@ -60,14 +70,23 @@ subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
 
    REALTYPE                  :: ccf
    REALTYPE                  :: x1,x2,x3
-!
 !EOP
 !-----------------------------------------------------------------------
 !BOC
 !  calculate cloud correction factor,fortran counts from 1 !
    ccf= cloud_correction_factor(nint(abs(dlat))+1)
 
    select case(method)
+      !before
+      !case(CLARK)
+      !after
+      case(FILE)
+         select case(type)
+            case(1)
+               ql=qlobs
+            case(2)
+               ql = qlobs-bolz*emiss*(tw**4)
+         end select
       case(CLARK)
 !        Clark et al. (1974) formula.
 !        unit of ea is Pascal, must hPa

src/fabm/gotm_fabm.F90

@@ -105,6 +105,8 @@ module gotm_fabm
                                 bioshade_feedback,bioalbedo_feedback,biodrag_feedback, &
                                 freshwater_impact,salinity_relaxation_to_freshwater_flux, &
                                 save_inputs, no_surface
+
+
 
    ! Arrays for work, vertical movement, and cross-boundary fluxes
    REALTYPE,allocatable,dimension(:,:) :: ws
@@ -250,7 +252,7 @@ subroutine configure_gotm_fabm(cfg)
 !
 ! !INPUT PARAMETERS:
    class (type_settings), intent(inout) :: cfg
-   type (type_settings), pointer :: branch
+   type (type_settings), pointer :: branch, twig
 !
 ! !REVISION HISTORY:
 !  Original author(s): Jorn Bruggeman
@@ -266,6 +268,7 @@ subroutine configure_gotm_fabm(cfg)
    ! Initialize all namelist variables to reasonable default values.
    call cfg%get(fabm_calc, 'use', 'enable FABM', &
                 default=.false.)
+
    call cfg%get(freshwater_impact, 'freshwater_impact', 'enable dilution/concentration by precipitation/evaporation', &
                 default=.true.) ! disable to check mass conservation
    branch => cfg%get_child('feedbacks', 'feedbacks to physics')
@@ -277,6 +280,7 @@ subroutine configure_gotm_fabm(cfg)
                 default=.false.)
    call cfg%get(repair_state, 'repair_state', 'clip state to minimum/maximum boundaries', &
                 default=.false.)
+
    branch => cfg%get_child('numerics', display=display_advanced)
    call branch%get(ode_method, 'ode_method', 'time integration scheme applied to source terms', &
                 options=(/ option(1, 'Forward Euler', 'FE'), option(2, 'Runge-Kutta 2', 'RK2'), option(3, 'Runge-Kutta 4', 'RK4'), &
@@ -306,6 +310,8 @@ subroutine configure_gotm_fabm(cfg)
                 options=(/option(-1, 'auto-detect (prefer fabm.yaml)', 'auto'), option(0, 'fabm.nml', 'nml'), option(1, 'fabm.yaml', 'yaml')/), &
                 default=-1, display=display_advanced)
 
+
+
    LEVEL2 'done.'
 
    end subroutine configure_gotm_fabm
@@ -699,6 +705,8 @@ subroutine init_var_gotm_fabm(nlev)
    ! and link it to FABM.
    decimal_yearday = _ZERO_
    call model%link_scalar(standard_variables%number_of_days_since_start_of_the_year,decimal_yearday)
+   call model%link_scalar(standard_variables%timestep,dt)  !jpnote added 
+
 
    allocate(Qsour(0:nlev),stat=rc)
    if (rc /= 0) stop 'init_var_gotm_fabm(): Error allocating (Qsour)'

src/gotm/gotm.F90

@@ -65,7 +65,6 @@ module gotm
    use stim_variables,  only: ice_hs,ice_hi,ice_uvic_topmelt,ice_uvic_topgrowth,ice_uvic_termelt,ice_uvic_botmelt,ice_uvic_botgrowth,ice_uvic_tb,ice_uvic_ts
    use stim_variables,  only: ice_uvic_Fs,ice_uvic_Ff,ice_uvic_parb,ice_uvic_parui,ice_uvic_Amelt
 
-
    use turbulence,  only: turb_method
    use turbulence,  only: init_turbulence,post_init_turbulence,do_turbulence
    use turbulence,  only: num,nuh,nus
@@ -569,6 +568,27 @@ subroutine initialize_gotm()
       call model_fabm%link_horizontal_data(standard_variables_fabm%bottom_depth,depth)
       call model_fabm%link_horizontal_data(standard_variables_fabm%bottom_depth_below_geoid,depth0)
       call model_fabm%link_horizontal_data(standard_variables_fabm%bottom_roughness_length,z0b)
+
+      !ice vars--------- jpnote 
+
+      call model_fabm%link_horizontal_data(standard_variables_fabm%sea_ice_thickness,ice_hi)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%snow_thickness,ice_hs)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%topmelt,ice_uvic_topmelt)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%f_melt,ice_uvic_Amelt)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%topgrowth,ice_uvic_topgrowth)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%termelt,ice_uvic_termelt)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%tendency_of_sea_ice_thickness_due_to_thermodynamics_melt,ice_uvic_botmelt)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%tendency_of_sea_ice_thickness_due_to_thermodynamics_grow,ice_uvic_botgrowth)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%sea_ice_temperature,ice_uvic_tb)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%surface_ice_temperature,ice_uvic_ts)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%lowest_ice_layer_PAR,ice_uvic_parb)
+      call model_fabm%link_horizontal_data(standard_variables_fabm%under_ice_PAR,ice_uvic_parui)
+
+      call model_fabm%link_interior_data(standard_variables_fabm%zonal_current,u(1:nlev)) !jpnote : changed from bulk data to interior data
+      call model_fabm%link_interior_data(standard_variables_fabm%meridional_current,v(1:nlev)) !jpnote: changed from bulk data to interior data
+      !------------------
+
+
       if (fluxes_method /= 0) then
          call model_fabm%link_horizontal_data(standard_variables_fabm%surface_specific_humidity,qa)
          call model_fabm%link_horizontal_data(standard_variables_fabm%surface_air_pressure,airp%value)
@@ -596,6 +616,9 @@ subroutine initialize_gotm()
 !  Initialize FABM initial state (this is done after the first call to do_input,
 !  to allow user-specified observed values to be used as initial state)
    if (fabm_calc) call init_gotm_fabm_state(nlev)
+
+
+
 #endif
 
    if (restart) then
@@ -800,7 +823,9 @@ subroutine integrate_gotm()
 
 !     update temperature and salinity
       if (sprof%method .ne. 0) then
+
          call salinity(nlev,dt,cnpar,nus,gams,ice_uvic_Fs,ice_uvic_Ff) 
+
       endif
 
       if (tprof%method .ne. 0) then

src/gotm/register_all_variables.F90

@@ -162,7 +162,8 @@ subroutine register_airsea_variables(nlev)
    call fm%register('I_0', 'W/m2', 'incoming short wave radiation', standard_name='', data0d=I_0%value, category='surface/heat_fluxes')
    call fm%register('qh', 'W/m2', 'sensible heat flux', standard_name='', data0d=qh, category='surface/heat_fluxes')
    call fm%register('qe', 'W/m2', 'latent heat flux', standard_name='', data0d=qe, category='surface/heat_fluxes')
-   call fm%register('ql', 'W/m2', 'net longwave radiation', standard_name='', data0d=ql%value, category='surface/heat_fluxes')
+   call fm%register('ql', 'W/m2', 'net longwave radiation', standard_name='', data0d=ql, category='surface/heat_fluxes')
+   call fm%register('qlobs', 'W/m2', 'longwave radiation (obs)', standard_name='', data0d=ql_%value, category='surface/heat_fluxes')
    call fm%register('heat', 'W/m2', 'net surface heat flux', standard_name='', data0d=heat%value, category='surface/heat_fluxes')
    call fm%register('tx', 'm2/s2', 'wind stress (x)', standard_name='', data0d=tx, category='surface')
    call fm%register('ty', 'm2/s2', 'wind stress (y)', standard_name='', data0d=ty, category='surface')

src/meanflow/salinity.F90

@@ -91,7 +91,7 @@ subroutine salinity(nlev,dt,cnpar,nus,gams,Fs,Ff) !jpnote passing in Fs and Ff f
 !
 !
 !  ice-ocean interaction salt and freshwater fluxes
-   REALTYPE, intent(in)                :: Fs,Ff  !jpnote added 
+   REALTYPE, optional, intent(in)                :: Fs,Ff  
 
 !
 ! !REVISION HISTORY:
@@ -117,7 +117,9 @@ subroutine salinity(nlev,dt,cnpar,nus,gams,Fs,Ff) !jpnote passing in Fs and Ff f
    DiffBcup       = Neumann
    DiffBcdw       = Neumann
    DiffSup        = -S(nlev)*(precip%value+evap+Ff)-Fs !NSnote, check signs
-!  DiffSup        = -S(nlev)*(precip%value+evap) 
+
+  !DiffSup        = -S(nlev)*(precip%value+evap) 
+
    DiffSdw        = _ZERO_
 
    AdvBcup       = oneSided