From d2aec78be8631987b508ff66830dd46c6d723b9b Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Thu, 5 Aug 2021 17:22:30 -0700
Subject: [PATCH 01/20] adding vars from gotm_fabm.nml and fabm.nml

---
 src/fabm/gotm_fabm.F90    | 170 +++++++++++++++++++++++++++++++++++++-
 src/gotm/gotm.F90         |  34 +++++++-
 src/meanflow/salinity.F90 |   7 +-
 3 files changed, 207 insertions(+), 4 deletions(-)

diff --git a/src/fabm/gotm_fabm.F90 b/src/fabm/gotm_fabm.F90
index 40bc5f87..b85afa3d 100644
--- a/src/fabm/gotm_fabm.F90
+++ b/src/fabm/gotm_fabm.F90
@@ -105,6 +105,24 @@ module gotm_fabm
                                 bioshade_feedback,bioalbedo_feedback,biodrag_feedback, &
                                 freshwater_impact,salinity_relaxation_to_freshwater_flux, &
                                 save_inputs, no_surface
+   
+
+
+   !declare fabm variables here? jpnote 
+      !Yaml variables 
+      !
+   logical                    :: no_precipitation_dilution
+   REALTYPE                   :: r_pond,fmethod,fflush,drag,f_graze,zia,ac_ia,rnit,skno3_0,sknh4_0, &
+                                 sksil_0,ia_0,ia_b,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp, &
+                                 t_sens,nu,md_no3,md_sil,chl2n,sil2n
+
+   REALTYPE :: ac,f_seed,ph1_0,ph2_0,zo1_0,zo2_0,no3_0,nh4_0,de1_0,de2_0,bsi_0,sil_0,w1,w2,mu1,mu2,kn,rpp1,rpp2,mp1,mp2,gz1,kz1,az1,az2,mz1,rc,pp1,pp2,pd1,pd2,pz1,gz2,kz2,mz2,rd1,rd2,rd3,rpf,rn0,knt,qp,qz,qb,agg,rsin,ks,pmin
+   !real(rk) :: drag,f_graze,zia,ac_ia,ia_0,ia_b,rnit,skno3_0,sknh4_0,sksil_0,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp,rpi,t_sens,nu,md_no3,md_sil,chl2n,sil2n
+   REALTYPE:: dic_0, alk_0, dic_sw, alk_sw, dic_ice, alk_ice, ik_diff, ik_on, ice_on, IA_on, tplv, btlv, prop2sw, prop2sw_melt
+   
+
+
+
 
    ! Arrays for work, vertical movement, and cross-boundary fluxes
    REALTYPE,allocatable,dimension(:,:) :: ws
@@ -250,7 +268,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 +284,8 @@ subroutine configure_gotm_fabm(cfg)
    ! Initialize all namelist variables to reasonable default values.
    call cfg%get(fabm_calc, 'use', 'enable FABM', &
                 default=.false.)
+   !fabm_calc ??? same as fabm calc 
+   
    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')
@@ -306,6 +326,154 @@ 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)
 
+   
+!mortenson FABM VARS from gotm.yaml
+   
+!fabm.nml
+!uvic_eco
+   branch => cfg%get_child('uvic_eco', 'University of Victoria Ecosystem model')
+   call branch%get(ac, 'ac', 'light attenuation coefficient','m-1', default=0.03_rk)
+   call branch%get(f_seed, 'f_seed', 'fraction of ice algal fux as ph2 seeding','-', default=0.0_rk)
+   call branch%get(ph1_0, 'ph1_0', 'ph1 initial value','umol/L', default=1.0_rk )
+   call branch%get(ph2_0 , 'ph2_0 ', 'ph2 initial value','umol/L', default=0.5_rk)
+   call branch%get(zo1_0, 'zo1_0', 'zo1 initial value','umol/L', default=0.2_rk)
+   call branch%get(zo2_0, 'zo2_0', 'zo2 initial value','umol/L', default=0.1_rk)
+   call branch%get(nh4_0, 'nh4_0', 'nh4 initial value','umol/L', default=10.0_rk)
+   call branch%get(no3_0, 'no3_0', 'no3 initial value','umol/L', default=10.0_rk)
+   call branch%get(de1_0, 'de1_0', 'de1 initial value ','umol/L', default=1.0_rk)
+   call branch%get(de2_0, 'de2_0', 'de2 initial value','umol/L', default=1.0_rk)
+   call branch%get(bsi_0, 'bsi_0', 'bsi initial value','umol/L', default=1.0_rk)
+   call branch%get(sil_0 , 'sil_0 ', 'sil initial value','umol/L', default=5.0_rk)
+   call branch%get(w1 , 'w1 ', 'de1 sinking rate','m/d', default=6.0_rk)
+   call branch%get(w2 , 'w2 ', 'de2 sinking rate','m/d', default=6.0_rk)
+   call branch%get(mu1, 'mu1', 'ph1 maximum growth rate','1/d', default=2.0_rk)
+   call branch%get(mu2 , 'mu2 ', 'ph2 maximum growth rate','1/d', default=2.0_rk)
+   call branch%get(kn, 'kn', ' no3 & nh4 half saturation constant','umol/L', default=0.1_rk)
+   !call branch%get(al1, 'al1', 'initial slope of P-I curve ([time] is same unit as pmax)','m2/W/time', default=0)
+   !call branch%get(al2 , 'al2 ', 'initial slope of P-I curve ([time] is same unit as pmax)','m2/W/time', default=0)
+   call branch%get(rpp1, 'rpp1', 'maximum photosynthetic rate ([time] is same unit as alpha)','1/time', default=0.05_rk)
+   call branch%get(rpp2 , 'rpp2 ', 'maximum photosynthetic rate ([time] is same unit as alpha)','1/time', default=0.05_rk)
+   call branch%get(mp1, 'mp1', 'ph1 excretion rate ','1/d', default=0.05_rk)
+   call branch%get(mp2 , 'mp2 ', 'ph2 excretion rate','1/d', default=0.05_rk)
+   call branch%get(mz1, 'mz1', 'zo1 excretion rate ','1/d', default=0.1_rk)
+   call branch%get(mz2, 'mz2', 'zo2 excretion rate ','1/d', default=0.3_rk)
+   call branch%get(gz1, 'gz1', 'zo1 maximum grazing rate ','1/d', default=1.3_rk)
+   call branch%get(gz2, 'gz2', 'zo2 maximum grazing rate ','1/d', default=0.8_rk)
+   call branch%get(kz1, 'kz1', 'zo1 grazing half saturation constant ','umol/L', default=0.6_rk)
+   call branch%get(kz2, 'kz2', 'zo2 grazing half saturation constant ','umol/L', default=0.75_rk)
+   call branch%get(az1, 'az1', 'zo1 assimilation fraction ','-', default=0.7_rk)
+   call branch%get(az2, 'az2', 'zo2 assimilation fraction ','-', default=0.1_rk)
+   call branch%get(rc, 'rc', 'closure mortality rate','1/d', default=0.003_rk)
+  ! call branch%get(htlnh4, 'htlnh4', 'closure loss fraction to nh4','-', default=0)
+  ! call branch%get(htldet, 'htldet', 'closure loss fraction to det','-', default=0)
+   call branch%get(pp1, 'pp1', 'ph1 fraction as food for zo1','-', default=1.0_rk)
+   call branch%get(pp2, 'pp2', 'ph2 fraction as food for zo2','-', default=1.0_rk)
+   call branch%get(pd1, 'pd1', 'de1 fraction as food for zo1','-', default=0.5_rk)
+   call branch%get(pd2, 'pd2', 'de2 fraction as food for zo2','-', default=0.5_rk)
+   call branch%get(pz1, 'pz1', 'zo1 fraction as food for zo2','-', default=1.0_rk)
+   call branch%get(rd1, 'rd1', ' de1 remineralization rate','1/d', default=0.1_rk)
+   call branch%get(rd2, 'rd2', 'de2 remineralization rate','1/d', default=0.1_rk)
+   call branch%get(rd3, 'rd3', 'bsi remineralization rate','1/d', default=0.1_rk)
+   call branch%get(rpf, 'rpf', 'scale factor for nitrogen preference function','-', default=0.2_rk )
+   call branch%get(rn0, 'rn0', 'nitrification rate at 0 deg.C','1/d', default=0.03_rk)
+   call branch%get(knt, 'knt', 'nitrification temperature coefficient','1/deg.C', default=0.0693_rk)
+   call branch%get(qp, 'qp', 'ph1 & ph2 Q10 factor','-', default=2.0_rk)
+   call branch%get(qz, 'qz', 'zo1 & zo2 Q10 factor','-', default=3.0_rk )
+   call branch%get(qb, 'qb', 'bacteria Q10 factor','-', default=3.0_rk)
+   call branch%get(agg, 'agg', 'ph2 aggregation rate','1/d', default=0.07_rk)
+   call branch%get(rsin , 'rsin ', 'ph2 Si:N ratio','mol-Si/mol-N', default=2.0_rk)
+   call branch%get(ks, 'ks', 'si half saturation constant','umol/L', default=2.0_rk)
+   call branch%get(pmin, 'pmin', 'background plankton concentration','umol-N/L', default=0.01_rk)
+
+!ice_algea
+   branch => cfg%get_child('uvic_icealgae', 'University of Victoria Ice Algae model')
+   call branch%get(r_pond, 'r_pond', 'melt pond drainage rate','', default=0.0175_rk)
+  ! call branch%get(fmethod, 'fmethod', 'method for ice-ocean flux ', '', &
+   !            default=0._rk, pchild=twig)
+   !call twig%get(fmethod_method, 'type', 'method for ice-ocean flux ', &
+    !             options=(/option(0, 'for diffusion', 'diffusion'), option(1, 'for ice growth/melting', 'ice growth/melting'), &
+     !            option(2, 'for both diffusion and growth/melting', 'both diffusion and growth/melting')/), default=0)
+                  !(0 for diffusion; 1 for ice growth/melting; 2 for both diffusion and growth/melting)
+   call branch%get(fmethod, 'fmethod ', 'method for ice-ocean flux','', default=0.0_rk)
+  ! call branch%get(fflush, 'fflush', 'method for flushing', '', &
+   !              default=0._rk, pchild=twig)
+   !call twig%get(fflush , 'fflush ', 'method for flushing', &
+   !              options=(/option(0, 'for no flushing', 'no flushing'), option(1, 'for surface flushing', 'surface flushing'), &
+   !             option(2, 'for surface+internal flushing', 'surface+internal flushing')/), default=0)
+   call branch%get(fflush , 'fflush ', 'method for flushing','', default=0.0_rk)
+   call branch%get(drag , 'drag ', 'drag coefficient at the ice-water interface ','-', default=0.005_rk)
+   call branch%get(f_graze, 'f_graze', 'fraction of ice algal growth lost due to grazing ','-', default=0.1_rk)
+   call branch%get(zia, 'zia', 'ice algal layer thickness ','m', default=0.03_rk)
+   call branch%get(ac_ia, 'ac_ia', 'specific light attenuation coefficient for ice algae','', default=0.03_rk)
+   call branch%get(rnit , 'rnit ', 'nitrification rate ','per day', default=0.1_rk)
+   call branch%get(ia_0 , 'ia_0 ', 'ia initial value ','mmol-N/m3', default=0.16_rk)
+   call branch%get(ia_b , 'ia_b ', 'ia background value ','mmol-N/m3', default=0.01_rk)
+   call branch%get(skno3_0, 'skno3_0', 'no3 initial value ','mmol/m3', default=2.0_rk)
+   call branch%get(sknh4_0, 'sknh4_0', 'nh4 initial value ','mmol/m3', default=0.01_rk)
+   call branch%get(sksil_0, 'sksil_0', 'sil initial value ','mmol/m3', default=5.0_rk)
+   call branch%get(ks_no3, 'ks_no3', 'no3 half-saturation value ','mmol/m3', default=1.0_rk)
+   call branch%get(ks_sil, 'ks_sil', 'sil half-saturation value ','mmol/m3', default=4.0_rk)
+   call branch%get(maxg, 'maxg', 'maximum specific growth rate ','d-1', default=0.8511_rk)
+   call branch%get(mort , 'mort ', 'linear mortality rate','d-1', default=0.05_rk)
+   call branch%get(mort2, 'mort2', 'quadratic mortality rate ','d-1', default=0.05_rk)
+   call branch%get(crit_melt, 'crit_melt', 'critical melt rate [m d-1]','m d-1', default=0.015_rk)
+   call branch%get(lcompp, 'lcompp', '# compensation intensity ','umol m-2 s-1', default=0.4_rk)
+   call branch%get(rpp , 'rpp ', 'ratio of photosynthetic parameters (alpha and pbm) [W m-2]-1','[W m-2]-1', default=0.1_rk)
+   !call branch%get(rpi , 'rpi ', 'ratio of photoinhibition parameters (beta and pbm)', default=0)
+   call branch%get(t_sens , 't_sens ', 'temperature sensitivity ','deg.C-1', default=0.0633_rk)
+   call branch%get(nu , 'nu ', 'kinematic viscosity?','', default=1.86e-6_rk)
+   call branch%get(md_no3, 'md_no3', 'molecular diffusion coefficient for nitrate','', default=0.47e-9_rk)
+   call branch%get(md_sil , 'md_sil ', 'molecular diffusion coefficient for dissolved silica','', default=0.47e-9_rk)
+   call branch%get(chl2n , 'chl2n ', 'chl to nitrogen ratio','', default=2.8_rk)
+   call branch%get(sil2n , 'sil2n ', 'silicon to nitrogen ratio','', default=1.7_rk)
+   
+!  uvic_icedic  
+   branch => cfg%get_child('uvic_icedic', 'University of Victoria Ice DIC model')   
+   call branch%get(dic_0, 'dic_0', 'initial DIC in water column','mmol/m3', default=2190.0_rk)
+   call branch%get(alk_0 , 'alk_0 ', 'initial TA in water column','mmol/m3', default=2100.0_rk)
+   call branch%get(dic_sw, 'dic_sw', 'dic_sw','', default=2100.0_rk)
+   call branch%get(alk_sw, 'alk_sw', 'alk_sw','', default=2200.0_rk)
+   call branch%get(dic_ice, 'dic_ice', '[DIC] for ice','mmol/m3', default=400.0_rk)  !400.0_rk 
+   call branch%get(alk_ice, 'alk_ice', '[TA] for growing ice','mmol/m3', default=500.0_rk)   !500.0_rk
+   call branch%get(ik_diff, 'ik_diff', 'difference (in melting ice) in [DIC] and 2*[TA] for ice with ikaite precipitaion','', default=50.0_rk)
+   call branch%get(ik_on, 'ik_on', '(0 or 1), turns off ikaite pump','0 or 1', default=1.0_rk)
+   call branch%get(ice_on, 'ice_on', '(0 or 1), turns off ice carbon pump','0 or 1', default=1.0_rk)
+   call branch%get(IA_on, 'IA_on', '(0 or 1), turns off ice algae carbon pump (not used now, bc can do same by ia_0=ia_b=0 in uvic_icealgae)','0 or 1', default=1.0_rk)
+   call branch%get(tplv, 'tplv', 'top of brine-associated DIC depth','m', default=40.0_rk)
+   call branch%get(btlv, 'btlv', 'bottom of brine-associated DIC depth','m', default=50.0_rk)
+   call branch%get(prop2sw, 'prop2sw', 'proportion of DIC rejected that is released into the ocean (the remainder presumably into the atmosphere)','', default=0.99_rk)
+   !prop2sw_melt=0.975_rk !SA: 1.0, def., 0.95, 0.9, 0.5  !right now, default = 0.975
+
+!gotm_fabm.nml  ---  jpnote some of these might be already existing settings --------- !but there isnt actually anywhere in the yaml where they are coming in from 
+   branch => cfg%get_child('gotm_fabm_nml', 'setting from the mortenson gotm_fabm.nml')  !jpnote change label/how it's organized in the yaml? 
+   call branch%get(fabm_calc, 'fabm_calc', 'fabm_calc', default=.false.)
+   call branch%get(cnpar, 'cnpar', 'cnpar','', default=1.0_rk)  
+   !call branch%get(w_adv_discr, 'w_adv_discr', 'w_adv_discr','', default=6.0_rk)
+   !call branch%get(ode_method, 'ode_method', 'ode_method','', default=1.0_rk)
+   !call branch%get(split_factor, 'split_factor', 'split_factor','', default=1.0_rk)
+   call branch%get(bioshade_feedback, 'bioshade_feedback', 'bioshade_feedback', default=.true.)
+   call branch%get(bioalbedo_feedback, 'bioalbedo_feedback', 'bioalbedo_feedback', default=.true.)
+   call branch%get(biodrag_feedback, 'biodrag_feedback', 'biodrag_feedback', default=.true.)
+   call branch%get(repair_state, 'repair_state', 'repair_state', default=.false.)
+   call branch%get(salinity_relaxation_to_freshwater_flux, 'salinity_relaxation_to_freshwater_flux', 'salinity_relaxation_to_freshwater_flux', default=.false.)
+   call branch%get(no_precipitation_dilution, 'no_precipitation_dilution', 'no_precipitation_dilution', default=.false.)
+   call branch%get(save_inputs, 'save_inputs', 'save_inputs', default=.false.)
+
+
+   !call branch%get(, '', '','', default=_rk)
+
+
+
+
+
+
+!fabm_input.nml 
+   
+
+
+
+
+   
    LEVEL2 'done.'
 
    end subroutine configure_gotm_fabm
diff --git a/src/gotm/gotm.F90 b/src/gotm/gotm.F90
index e7b5c168..e103bce6 100644
--- a/src/gotm/gotm.F90
+++ b/src/gotm/gotm.F90
@@ -61,9 +61,9 @@ module gotm
 #ifdef _ICE_
    use ice,         only: init_ice, post_init_ice, do_ice, clean_ice, ice_cover
    use stim_variables,  only: Tice_surface,albedo_ice,transmissivity,nilay,ice_uvic_Tice
+#endif
    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
-#endif
 
    use turbulence,  only: turb_method
    use turbulence,  only: init_turbulence,post_init_turbulence,do_turbulence
@@ -568,6 +568,28 @@ 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 
+#if 0
+      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_bulk_data(standard_variables_fabm%zonal_current,u(1:nlev))
+      call model_fabm%link_bulk_data(standard_variables_fabm%meridional_current,v(1:nlev))
+      !------------------
+#endif
+
       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)
@@ -595,6 +617,10 @@ 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)
+
+
+   !jpnote call init_gotm_fabm_output? 
+
 #endif
 
    if (restart) then
@@ -796,7 +822,13 @@ subroutine integrate_gotm()
 
 !     update temperature and salinity
       if (sprof%method .ne. 0) then
+#ifdef _ICE_
+
          call salinity(nlev,dt,cnpar,nus,gams,ice_uvic_Fs,ice_uvic_Ff) 
+#else
+         call salinity(nlev,dt,cnpar,nus,gams)
+#endif
+
       endif
 
       if (tprof%method .ne. 0) then
diff --git a/src/meanflow/salinity.F90 b/src/meanflow/salinity.F90
index fd1fd4c6..eb30aa42 100644
--- a/src/meanflow/salinity.F90
+++ b/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  !jpnote added  !optional - not needed when not running ice? 
 
 !
 ! !REVISION HISTORY:
@@ -116,8 +116,11 @@ subroutine salinity(nlev,dt,cnpar,nus,gams,Fs,Ff) !jpnote passing in Fs and Ff f
 !  set boundary conditions
    DiffBcup       = Neumann
    DiffBcdw       = Neumann
+#ifdef _ICE_
    DiffSup        = -S(nlev)*(precip%value+evap+Ff)-Fs !NSnote, check signs
-!  DiffSup        = -S(nlev)*(precip%value+evap) 
+#else
+  DiffSup        = -S(nlev)*(precip%value+evap) 
+#endif
    DiffSdw        = _ZERO_
 
    AdvBcup       = oneSided

From 19f0a352c4a55dd9332fb03bbf48fa3af7ba727b Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Fri, 6 Aug 2021 09:44:00 -0700
Subject: [PATCH 02/20] removed _ICE_ condition from salinity

---
 src/gotm/gotm.F90         | 4 ----
 src/meanflow/salinity.F90 | 7 +++----
 2 files changed, 3 insertions(+), 8 deletions(-)

diff --git a/src/gotm/gotm.F90 b/src/gotm/gotm.F90
index e103bce6..4095ea41 100644
--- a/src/gotm/gotm.F90
+++ b/src/gotm/gotm.F90
@@ -822,12 +822,8 @@ subroutine integrate_gotm()
 
 !     update temperature and salinity
       if (sprof%method .ne. 0) then
-#ifdef _ICE_
 
          call salinity(nlev,dt,cnpar,nus,gams,ice_uvic_Fs,ice_uvic_Ff) 
-#else
-         call salinity(nlev,dt,cnpar,nus,gams)
-#endif
 
       endif
 
diff --git a/src/meanflow/salinity.F90 b/src/meanflow/salinity.F90
index eb30aa42..267f6450 100644
--- a/src/meanflow/salinity.F90
+++ b/src/meanflow/salinity.F90
@@ -116,11 +116,10 @@ subroutine salinity(nlev,dt,cnpar,nus,gams,Fs,Ff) !jpnote passing in Fs and Ff f
 !  set boundary conditions
    DiffBcup       = Neumann
    DiffBcdw       = Neumann
-#ifdef _ICE_
    DiffSup        = -S(nlev)*(precip%value+evap+Ff)-Fs !NSnote, check signs
-#else
-  DiffSup        = -S(nlev)*(precip%value+evap) 
-#endif
+
+  !DiffSup        = -S(nlev)*(precip%value+evap) 
+
    DiffSdw        = _ZERO_
 
    AdvBcup       = oneSided

From 223c40628898a7c678b1b585d4a83dda35c9a5a1 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Mon, 9 Aug 2021 14:54:18 -0700
Subject: [PATCH 03/20] adding model uvic to fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index ab8a6ffc..bbfd4720 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit ab8a6ffc5ceecd6e467878fe7e3142adcc49ab3a
+Subproject commit bbfd4720db23295f15ff964b0b54459b4a51903c

From bc861d5642723cf7695353b53abc59d2cc76356c Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 10 Aug 2021 11:21:21 -0700
Subject: [PATCH 04/20] changes for new model

---
 src/fabm/gotm_fabm.F90 | 96 +-----------------------------------------
 src/gotm/gotm.F90      | 10 +++--
 2 files changed, 8 insertions(+), 98 deletions(-)

diff --git a/src/fabm/gotm_fabm.F90 b/src/fabm/gotm_fabm.F90
index b85afa3d..9fe70aef 100644
--- a/src/fabm/gotm_fabm.F90
+++ b/src/fabm/gotm_fabm.F90
@@ -330,102 +330,10 @@ subroutine configure_gotm_fabm(cfg)
 !mortenson FABM VARS from gotm.yaml
    
 !fabm.nml
-!uvic_eco
-   branch => cfg%get_child('uvic_eco', 'University of Victoria Ecosystem model')
-   call branch%get(ac, 'ac', 'light attenuation coefficient','m-1', default=0.03_rk)
-   call branch%get(f_seed, 'f_seed', 'fraction of ice algal fux as ph2 seeding','-', default=0.0_rk)
-   call branch%get(ph1_0, 'ph1_0', 'ph1 initial value','umol/L', default=1.0_rk )
-   call branch%get(ph2_0 , 'ph2_0 ', 'ph2 initial value','umol/L', default=0.5_rk)
-   call branch%get(zo1_0, 'zo1_0', 'zo1 initial value','umol/L', default=0.2_rk)
-   call branch%get(zo2_0, 'zo2_0', 'zo2 initial value','umol/L', default=0.1_rk)
-   call branch%get(nh4_0, 'nh4_0', 'nh4 initial value','umol/L', default=10.0_rk)
-   call branch%get(no3_0, 'no3_0', 'no3 initial value','umol/L', default=10.0_rk)
-   call branch%get(de1_0, 'de1_0', 'de1 initial value ','umol/L', default=1.0_rk)
-   call branch%get(de2_0, 'de2_0', 'de2 initial value','umol/L', default=1.0_rk)
-   call branch%get(bsi_0, 'bsi_0', 'bsi initial value','umol/L', default=1.0_rk)
-   call branch%get(sil_0 , 'sil_0 ', 'sil initial value','umol/L', default=5.0_rk)
-   call branch%get(w1 , 'w1 ', 'de1 sinking rate','m/d', default=6.0_rk)
-   call branch%get(w2 , 'w2 ', 'de2 sinking rate','m/d', default=6.0_rk)
-   call branch%get(mu1, 'mu1', 'ph1 maximum growth rate','1/d', default=2.0_rk)
-   call branch%get(mu2 , 'mu2 ', 'ph2 maximum growth rate','1/d', default=2.0_rk)
-   call branch%get(kn, 'kn', ' no3 & nh4 half saturation constant','umol/L', default=0.1_rk)
-   !call branch%get(al1, 'al1', 'initial slope of P-I curve ([time] is same unit as pmax)','m2/W/time', default=0)
-   !call branch%get(al2 , 'al2 ', 'initial slope of P-I curve ([time] is same unit as pmax)','m2/W/time', default=0)
-   call branch%get(rpp1, 'rpp1', 'maximum photosynthetic rate ([time] is same unit as alpha)','1/time', default=0.05_rk)
-   call branch%get(rpp2 , 'rpp2 ', 'maximum photosynthetic rate ([time] is same unit as alpha)','1/time', default=0.05_rk)
-   call branch%get(mp1, 'mp1', 'ph1 excretion rate ','1/d', default=0.05_rk)
-   call branch%get(mp2 , 'mp2 ', 'ph2 excretion rate','1/d', default=0.05_rk)
-   call branch%get(mz1, 'mz1', 'zo1 excretion rate ','1/d', default=0.1_rk)
-   call branch%get(mz2, 'mz2', 'zo2 excretion rate ','1/d', default=0.3_rk)
-   call branch%get(gz1, 'gz1', 'zo1 maximum grazing rate ','1/d', default=1.3_rk)
-   call branch%get(gz2, 'gz2', 'zo2 maximum grazing rate ','1/d', default=0.8_rk)
-   call branch%get(kz1, 'kz1', 'zo1 grazing half saturation constant ','umol/L', default=0.6_rk)
-   call branch%get(kz2, 'kz2', 'zo2 grazing half saturation constant ','umol/L', default=0.75_rk)
-   call branch%get(az1, 'az1', 'zo1 assimilation fraction ','-', default=0.7_rk)
-   call branch%get(az2, 'az2', 'zo2 assimilation fraction ','-', default=0.1_rk)
-   call branch%get(rc, 'rc', 'closure mortality rate','1/d', default=0.003_rk)
-  ! call branch%get(htlnh4, 'htlnh4', 'closure loss fraction to nh4','-', default=0)
-  ! call branch%get(htldet, 'htldet', 'closure loss fraction to det','-', default=0)
-   call branch%get(pp1, 'pp1', 'ph1 fraction as food for zo1','-', default=1.0_rk)
-   call branch%get(pp2, 'pp2', 'ph2 fraction as food for zo2','-', default=1.0_rk)
-   call branch%get(pd1, 'pd1', 'de1 fraction as food for zo1','-', default=0.5_rk)
-   call branch%get(pd2, 'pd2', 'de2 fraction as food for zo2','-', default=0.5_rk)
-   call branch%get(pz1, 'pz1', 'zo1 fraction as food for zo2','-', default=1.0_rk)
-   call branch%get(rd1, 'rd1', ' de1 remineralization rate','1/d', default=0.1_rk)
-   call branch%get(rd2, 'rd2', 'de2 remineralization rate','1/d', default=0.1_rk)
-   call branch%get(rd3, 'rd3', 'bsi remineralization rate','1/d', default=0.1_rk)
-   call branch%get(rpf, 'rpf', 'scale factor for nitrogen preference function','-', default=0.2_rk )
-   call branch%get(rn0, 'rn0', 'nitrification rate at 0 deg.C','1/d', default=0.03_rk)
-   call branch%get(knt, 'knt', 'nitrification temperature coefficient','1/deg.C', default=0.0693_rk)
-   call branch%get(qp, 'qp', 'ph1 & ph2 Q10 factor','-', default=2.0_rk)
-   call branch%get(qz, 'qz', 'zo1 & zo2 Q10 factor','-', default=3.0_rk )
-   call branch%get(qb, 'qb', 'bacteria Q10 factor','-', default=3.0_rk)
-   call branch%get(agg, 'agg', 'ph2 aggregation rate','1/d', default=0.07_rk)
-   call branch%get(rsin , 'rsin ', 'ph2 Si:N ratio','mol-Si/mol-N', default=2.0_rk)
-   call branch%get(ks, 'ks', 'si half saturation constant','umol/L', default=2.0_rk)
-   call branch%get(pmin, 'pmin', 'background plankton concentration','umol-N/L', default=0.01_rk)
+
 
 !ice_algea
-   branch => cfg%get_child('uvic_icealgae', 'University of Victoria Ice Algae model')
-   call branch%get(r_pond, 'r_pond', 'melt pond drainage rate','', default=0.0175_rk)
-  ! call branch%get(fmethod, 'fmethod', 'method for ice-ocean flux ', '', &
-   !            default=0._rk, pchild=twig)
-   !call twig%get(fmethod_method, 'type', 'method for ice-ocean flux ', &
-    !             options=(/option(0, 'for diffusion', 'diffusion'), option(1, 'for ice growth/melting', 'ice growth/melting'), &
-     !            option(2, 'for both diffusion and growth/melting', 'both diffusion and growth/melting')/), default=0)
-                  !(0 for diffusion; 1 for ice growth/melting; 2 for both diffusion and growth/melting)
-   call branch%get(fmethod, 'fmethod ', 'method for ice-ocean flux','', default=0.0_rk)
-  ! call branch%get(fflush, 'fflush', 'method for flushing', '', &
-   !              default=0._rk, pchild=twig)
-   !call twig%get(fflush , 'fflush ', 'method for flushing', &
-   !              options=(/option(0, 'for no flushing', 'no flushing'), option(1, 'for surface flushing', 'surface flushing'), &
-   !             option(2, 'for surface+internal flushing', 'surface+internal flushing')/), default=0)
-   call branch%get(fflush , 'fflush ', 'method for flushing','', default=0.0_rk)
-   call branch%get(drag , 'drag ', 'drag coefficient at the ice-water interface ','-', default=0.005_rk)
-   call branch%get(f_graze, 'f_graze', 'fraction of ice algal growth lost due to grazing ','-', default=0.1_rk)
-   call branch%get(zia, 'zia', 'ice algal layer thickness ','m', default=0.03_rk)
-   call branch%get(ac_ia, 'ac_ia', 'specific light attenuation coefficient for ice algae','', default=0.03_rk)
-   call branch%get(rnit , 'rnit ', 'nitrification rate ','per day', default=0.1_rk)
-   call branch%get(ia_0 , 'ia_0 ', 'ia initial value ','mmol-N/m3', default=0.16_rk)
-   call branch%get(ia_b , 'ia_b ', 'ia background value ','mmol-N/m3', default=0.01_rk)
-   call branch%get(skno3_0, 'skno3_0', 'no3 initial value ','mmol/m3', default=2.0_rk)
-   call branch%get(sknh4_0, 'sknh4_0', 'nh4 initial value ','mmol/m3', default=0.01_rk)
-   call branch%get(sksil_0, 'sksil_0', 'sil initial value ','mmol/m3', default=5.0_rk)
-   call branch%get(ks_no3, 'ks_no3', 'no3 half-saturation value ','mmol/m3', default=1.0_rk)
-   call branch%get(ks_sil, 'ks_sil', 'sil half-saturation value ','mmol/m3', default=4.0_rk)
-   call branch%get(maxg, 'maxg', 'maximum specific growth rate ','d-1', default=0.8511_rk)
-   call branch%get(mort , 'mort ', 'linear mortality rate','d-1', default=0.05_rk)
-   call branch%get(mort2, 'mort2', 'quadratic mortality rate ','d-1', default=0.05_rk)
-   call branch%get(crit_melt, 'crit_melt', 'critical melt rate [m d-1]','m d-1', default=0.015_rk)
-   call branch%get(lcompp, 'lcompp', '# compensation intensity ','umol m-2 s-1', default=0.4_rk)
-   call branch%get(rpp , 'rpp ', 'ratio of photosynthetic parameters (alpha and pbm) [W m-2]-1','[W m-2]-1', default=0.1_rk)
-   !call branch%get(rpi , 'rpi ', 'ratio of photoinhibition parameters (beta and pbm)', default=0)
-   call branch%get(t_sens , 't_sens ', 'temperature sensitivity ','deg.C-1', default=0.0633_rk)
-   call branch%get(nu , 'nu ', 'kinematic viscosity?','', default=1.86e-6_rk)
-   call branch%get(md_no3, 'md_no3', 'molecular diffusion coefficient for nitrate','', default=0.47e-9_rk)
-   call branch%get(md_sil , 'md_sil ', 'molecular diffusion coefficient for dissolved silica','', default=0.47e-9_rk)
-   call branch%get(chl2n , 'chl2n ', 'chl to nitrogen ratio','', default=2.8_rk)
-   call branch%get(sil2n , 'sil2n ', 'silicon to nitrogen ratio','', default=1.7_rk)
+  
    
 !  uvic_icedic  
    branch => cfg%get_child('uvic_icedic', 'University of Victoria Ice DIC model')   
diff --git a/src/gotm/gotm.F90 b/src/gotm/gotm.F90
index 4095ea41..4400525b 100644
--- a/src/gotm/gotm.F90
+++ b/src/gotm/gotm.F90
@@ -401,6 +401,8 @@ subroutine initialize_gotm()
 #ifdef _FABM_
    if (read_nml) call configure_gotm_fabm_from_nml(namlst, 'gotm_fabm.nml')
 
+   !jpnote ?? look at gotm_fabm nml and fabm_input.nml .. 
+
    ! Allow FABM to create its model tree. After this we know all biogeochemical variables
    ! This must be done before gotm_fabm_input configuration.
    call gotm_fabm_create_model(namlst)
@@ -571,7 +573,7 @@ subroutine initialize_gotm()
 
 
       !ice vars--------- jpnote 
-#if 0
+
       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)
@@ -585,10 +587,10 @@ subroutine initialize_gotm()
       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_bulk_data(standard_variables_fabm%zonal_current,u(1:nlev))
-      call model_fabm%link_bulk_data(standard_variables_fabm%meridional_current,v(1:nlev))
+      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
       !------------------
-#endif
+
 
       if (fluxes_method /= 0) then
          call model_fabm%link_horizontal_data(standard_variables_fabm%surface_specific_humidity,qa)

From 8138ac0ae1258cffd8223069bc21b75b7c4b4756 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Mon, 16 Aug 2021 14:51:48 -0700
Subject: [PATCH 05/20] synch

---
 src/fabm/gotm_fabm.F90 | 48 ++++++++++++++++--------------------------
 1 file changed, 18 insertions(+), 30 deletions(-)

diff --git a/src/fabm/gotm_fabm.F90 b/src/fabm/gotm_fabm.F90
index 003fc421..36429570 100644
--- a/src/fabm/gotm_fabm.F90
+++ b/src/fabm/gotm_fabm.F90
@@ -111,12 +111,13 @@ module gotm_fabm
    !declare fabm variables here? jpnote 
       !Yaml variables 
       !
-   logical                    :: no_precipitation_dilution
-   REALTYPE                   :: r_pond,fmethod,fflush,drag,f_graze,zia,ac_ia,rnit,skno3_0,sknh4_0, &
-                                 sksil_0,ia_0,ia_b,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp, &
-                                 t_sens,nu,md_no3,md_sil,chl2n,sil2n
+!logical :: use_icealgae
+   !logical                    :: no_precipitation_dilution
+   !REALTYPE                   :: r_pond,fmethod,fflush,drag,f_graze,zia,ac_ia,rnit,skno3_0,sknh4_0, &
+    !                             sksil_0,ia_0,ia_b,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp, &
+     !                            t_sens,nu,md_no3,md_sil,chl2n,sil2n
 
-   REALTYPE :: ac,f_seed,ph1_0,ph2_0,zo1_0,zo2_0,no3_0,nh4_0,de1_0,de2_0,bsi_0,sil_0,w1,w2,mu1,mu2,kn,rpp1,rpp2,mp1,mp2,gz1,kz1,az1,az2,mz1,rc,pp1,pp2,pd1,pd2,pz1,gz2,kz2,mz2,rd1,rd2,rd3,rpf,rn0,knt,qp,qz,qb,agg,rsin,ks,pmin
+   !REALTYPE :: ac,f_seed,ph1_0,ph2_0,zo1_0,zo2_0,no3_0,nh4_0,de1_0,de2_0,bsi_0,sil_0,w1,w2,mu1,mu2,kn,rpp1,rpp2,mp1,mp2,gz1,kz1,az1,az2,mz1,rc,pp1,pp2,pd1,pd2,pz1,gz2,kz2,mz2,rd1,rd2,rd3,rpf,rn0,knt,qp,qz,qb,agg,rsin,ks,pmin
    !real(rk) :: drag,f_graze,zia,ac_ia,ia_0,ia_b,rnit,skno3_0,sknh4_0,sksil_0,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp,rpi,t_sens,nu,md_no3,md_sil,chl2n,sil2n
    REALTYPE:: dic_0, alk_0, dic_sw, alk_sw, dic_ice, alk_ice, ik_diff, ik_on, ice_on, IA_on, tplv, btlv, prop2sw, prop2sw_melt
    
@@ -336,36 +337,22 @@ subroutine configure_gotm_fabm(cfg)
   
    
 !  uvic_icedic  
-   branch => cfg%get_child('uvic_icedic', 'University of Victoria Ice DIC model')   
-   call branch%get(dic_0, 'dic_0', 'initial DIC in water column','mmol/m3', default=2190.0_rk)
-   call branch%get(alk_0 , 'alk_0 ', 'initial TA in water column','mmol/m3', default=2100.0_rk)
-   call branch%get(dic_sw, 'dic_sw', 'dic_sw','', default=2100.0_rk)
-   call branch%get(alk_sw, 'alk_sw', 'alk_sw','', default=2200.0_rk)
-   call branch%get(dic_ice, 'dic_ice', '[DIC] for ice','mmol/m3', default=400.0_rk)  !400.0_rk 
-   call branch%get(alk_ice, 'alk_ice', '[TA] for growing ice','mmol/m3', default=500.0_rk)   !500.0_rk
-   call branch%get(ik_diff, 'ik_diff', 'difference (in melting ice) in [DIC] and 2*[TA] for ice with ikaite precipitaion','', default=50.0_rk)
-   call branch%get(ik_on, 'ik_on', '(0 or 1), turns off ikaite pump','0 or 1', default=1.0_rk)
-   call branch%get(ice_on, 'ice_on', '(0 or 1), turns off ice carbon pump','0 or 1', default=1.0_rk)
-   call branch%get(IA_on, 'IA_on', '(0 or 1), turns off ice algae carbon pump (not used now, bc can do same by ia_0=ia_b=0 in uvic_icealgae)','0 or 1', default=1.0_rk)
-   call branch%get(tplv, 'tplv', 'top of brine-associated DIC depth','m', default=40.0_rk)
-   call branch%get(btlv, 'btlv', 'bottom of brine-associated DIC depth','m', default=50.0_rk)
-   call branch%get(prop2sw, 'prop2sw', 'proportion of DIC rejected that is released into the ocean (the remainder presumably into the atmosphere)','', default=0.99_rk)
-   !prop2sw_melt=0.975_rk !SA: 1.0, def., 0.95, 0.9, 0.5  !right now, default = 0.975
+
 
 !gotm_fabm.nml  ---  jpnote some of these might be already existing settings --------- !but there isnt actually anywhere in the yaml where they are coming in from 
-   branch => cfg%get_child('gotm_fabm_nml', 'setting from the mortenson gotm_fabm.nml')  !jpnote change label/how it's organized in the yaml? 
-   call branch%get(fabm_calc, 'fabm_calc', 'fabm_calc', default=.false.)
-   call branch%get(cnpar, 'cnpar', 'cnpar','', default=1.0_rk)  
+   !branch => cfg%get_child('gotm_fabm_nml', 'setting from the mortenson gotm_fabm.nml')  !jpnote change label/how it's organized in the yaml? 
+   !call branch%get(fabm_calc, 'fabm_calc', 'fabm_calc', default=.false.)
+   !call branch%get(cnpar, 'cnpar', 'cnpar','', default=1.0_rk)  
    !call branch%get(w_adv_discr, 'w_adv_discr', 'w_adv_discr','', default=6.0_rk)
    !call branch%get(ode_method, 'ode_method', 'ode_method','', default=1.0_rk)
    !call branch%get(split_factor, 'split_factor', 'split_factor','', default=1.0_rk)
-   call branch%get(bioshade_feedback, 'bioshade_feedback', 'bioshade_feedback', default=.true.)
-   call branch%get(bioalbedo_feedback, 'bioalbedo_feedback', 'bioalbedo_feedback', default=.true.)
-   call branch%get(biodrag_feedback, 'biodrag_feedback', 'biodrag_feedback', default=.true.)
-   call branch%get(repair_state, 'repair_state', 'repair_state', default=.false.)
-   call branch%get(salinity_relaxation_to_freshwater_flux, 'salinity_relaxation_to_freshwater_flux', 'salinity_relaxation_to_freshwater_flux', default=.false.)
-   call branch%get(no_precipitation_dilution, 'no_precipitation_dilution', 'no_precipitation_dilution', default=.false.)
-   call branch%get(save_inputs, 'save_inputs', 'save_inputs', default=.false.)
+   !call branch%get(bioshade_feedback, 'bioshade_feedback', 'bioshade_feedback', default=.true.)
+  ! call branch%get(bioalbedo_feedback, 'bioalbedo_feedback', 'bioalbedo_feedback', default=.true.)
+   !call branch%get(biodrag_feedback, 'biodrag_feedback', 'biodrag_feedback', default=.true.)
+   !call branch%get(repair_state, 'repair_state', 'repair_state', default=.false.)
+   !call branch%get(salinity_relaxation_to_freshwater_flux, 'salinity_relaxation_to_freshwater_flux', 'salinity_relaxation_to_freshwater_flux', default=.false.)
+   !call branch%get(no_precipitation_dilution, 'no_precipitation_dilution', 'no_precipitation_dilution', default=.false.)
+   !call branch%get(save_inputs, 'save_inputs', 'save_inputs', default=.false.)
 
 
    !call branch%get(, '', '','', default=_rk)
@@ -775,6 +762,7 @@ 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)'

From 968421137e894f5337fcf2d58e1475effc9952ec Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 24 Aug 2021 11:16:13 -0700
Subject: [PATCH 06/20] clean

---
 extern/fabm               |  2 +-
 src/fabm/gotm_fabm.F90    | 37 ++++++++-----------------------------
 src/gotm/gotm.F90         |  4 ----
 src/meanflow/salinity.F90 |  2 +-
 4 files changed, 10 insertions(+), 35 deletions(-)

diff --git a/extern/fabm b/extern/fabm
index bbfd4720..378ae847 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit bbfd4720db23295f15ff964b0b54459b4a51903c
+Subproject commit 378ae8476bc1586cd7b56b4b31025ad7bb934fef
diff --git a/src/fabm/gotm_fabm.F90 b/src/fabm/gotm_fabm.F90
index 36429570..80088b73 100644
--- a/src/fabm/gotm_fabm.F90
+++ b/src/fabm/gotm_fabm.F90
@@ -119,7 +119,7 @@ module gotm_fabm
 
    !REALTYPE :: ac,f_seed,ph1_0,ph2_0,zo1_0,zo2_0,no3_0,nh4_0,de1_0,de2_0,bsi_0,sil_0,w1,w2,mu1,mu2,kn,rpp1,rpp2,mp1,mp2,gz1,kz1,az1,az2,mz1,rc,pp1,pp2,pd1,pd2,pz1,gz2,kz2,mz2,rd1,rd2,rd3,rpf,rn0,knt,qp,qz,qb,agg,rsin,ks,pmin
    !real(rk) :: drag,f_graze,zia,ac_ia,ia_0,ia_b,rnit,skno3_0,sknh4_0,sksil_0,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp,rpi,t_sens,nu,md_no3,md_sil,chl2n,sil2n
-   REALTYPE:: dic_0, alk_0, dic_sw, alk_sw, dic_ice, alk_ice, ik_diff, ik_on, ice_on, IA_on, tplv, btlv, prop2sw, prop2sw_melt
+  ! REALTYPE:: dic_0, alk_0, dic_sw, alk_sw, dic_ice, alk_ice, ik_diff, ik_on, ice_on, IA_on, tplv, btlv, prop2sw, prop2sw_melt
    
 
 
@@ -285,7 +285,6 @@ subroutine configure_gotm_fabm(cfg)
    ! Initialize all namelist variables to reasonable default values.
    call cfg%get(fabm_calc, 'use', 'enable FABM', &
                 default=.false.)
-   !fabm_calc ??? same as fabm calc 
    
    call cfg%get(freshwater_impact, 'freshwater_impact', 'enable dilution/concentration by precipitation/evaporation', &
                 default=.true.) ! disable to check mass conservation
@@ -298,6 +297,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'), &
@@ -329,17 +329,8 @@ subroutine configure_gotm_fabm(cfg)
 
    
 !mortenson FABM VARS from gotm.yaml
-   
-!fabm.nml
-
-
-!ice_algea
-  
-   
-!  uvic_icedic  
-
 
-!gotm_fabm.nml  ---  jpnote some of these might be already existing settings --------- !but there isnt actually anywhere in the yaml where they are coming in from 
+!gotm_fabm.nml  jpnote some of these might be already existing settings --------- !but there isnt actually anywhere in the yaml where they are coming in from 
    !branch => cfg%get_child('gotm_fabm_nml', 'setting from the mortenson gotm_fabm.nml')  !jpnote change label/how it's organized in the yaml? 
    !call branch%get(fabm_calc, 'fabm_calc', 'fabm_calc', default=.false.)
    !call branch%get(cnpar, 'cnpar', 'cnpar','', default=1.0_rk)  
@@ -353,22 +344,7 @@ subroutine configure_gotm_fabm(cfg)
    !call branch%get(salinity_relaxation_to_freshwater_flux, 'salinity_relaxation_to_freshwater_flux', 'salinity_relaxation_to_freshwater_flux', default=.false.)
    !call branch%get(no_precipitation_dilution, 'no_precipitation_dilution', 'no_precipitation_dilution', default=.false.)
    !call branch%get(save_inputs, 'save_inputs', 'save_inputs', default=.false.)
-
-
-   !call branch%get(, '', '','', default=_rk)
-
-
-
-
-
-
-!fabm_input.nml 
-   
-
-
-
-
-   
+    
    LEVEL2 'done.'
 
    end subroutine configure_gotm_fabm
@@ -762,7 +738,10 @@ 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 
+   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)'
diff --git a/src/gotm/gotm.F90 b/src/gotm/gotm.F90
index 4400525b..9ae42abc 100644
--- a/src/gotm/gotm.F90
+++ b/src/gotm/gotm.F90
@@ -401,8 +401,6 @@ subroutine initialize_gotm()
 #ifdef _FABM_
    if (read_nml) call configure_gotm_fabm_from_nml(namlst, 'gotm_fabm.nml')
 
-   !jpnote ?? look at gotm_fabm nml and fabm_input.nml .. 
-
    ! Allow FABM to create its model tree. After this we know all biogeochemical variables
    ! This must be done before gotm_fabm_input configuration.
    call gotm_fabm_create_model(namlst)
@@ -571,7 +569,6 @@ subroutine initialize_gotm()
       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)
@@ -621,7 +618,6 @@ subroutine initialize_gotm()
    if (fabm_calc) call init_gotm_fabm_state(nlev)
 
 
-   !jpnote call init_gotm_fabm_output? 
 
 #endif
 
diff --git a/src/meanflow/salinity.F90 b/src/meanflow/salinity.F90
index 267f6450..49ce3d87 100644
--- a/src/meanflow/salinity.F90
+++ b/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, optional, intent(in)                :: Fs,Ff  !jpnote added  !optional - not needed when not running ice? 
+   REALTYPE, optional, intent(in)                :: Fs,Ff  
 
 !
 ! !REVISION HISTORY:

From 3269b2bad4637b0df923aa18737c7ce86b889158 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 31 Aug 2021 17:24:05 -0700
Subject: [PATCH 07/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index 378ae847..a6e3d10f 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 378ae8476bc1586cd7b56b4b31025ad7bb934fef
+Subproject commit a6e3d10f7eb1e3018a35042295c8806b569005a7

From 9958834fb6018ceeaf19124c02e5e0cd7fc48c7a Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 31 Aug 2021 17:24:05 -0700
Subject: [PATCH 08/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index 378ae847..d938a3e2 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 378ae8476bc1586cd7b56b4b31025ad7bb934fef
+Subproject commit d938a3e25ade729fe93d3fba496165badcea6b84

From ac8735f0c870e20eadc9f8dcdd23c24d20c5de5f Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 31 Aug 2021 17:33:45 -0700
Subject: [PATCH 09/20] clean

---
 src/fabm/gotm_fabm.F90 | 37 +------------------------------------
 1 file changed, 1 insertion(+), 36 deletions(-)

diff --git a/src/fabm/gotm_fabm.F90 b/src/fabm/gotm_fabm.F90
index 80088b73..8945afa3 100644
--- a/src/fabm/gotm_fabm.F90
+++ b/src/fabm/gotm_fabm.F90
@@ -108,23 +108,6 @@ module gotm_fabm
    
 
 
-   !declare fabm variables here? jpnote 
-      !Yaml variables 
-      !
-!logical :: use_icealgae
-   !logical                    :: no_precipitation_dilution
-   !REALTYPE                   :: r_pond,fmethod,fflush,drag,f_graze,zia,ac_ia,rnit,skno3_0,sknh4_0, &
-    !                             sksil_0,ia_0,ia_b,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp, &
-     !                            t_sens,nu,md_no3,md_sil,chl2n,sil2n
-
-   !REALTYPE :: ac,f_seed,ph1_0,ph2_0,zo1_0,zo2_0,no3_0,nh4_0,de1_0,de2_0,bsi_0,sil_0,w1,w2,mu1,mu2,kn,rpp1,rpp2,mp1,mp2,gz1,kz1,az1,az2,mz1,rc,pp1,pp2,pd1,pd2,pz1,gz2,kz2,mz2,rd1,rd2,rd3,rpf,rn0,knt,qp,qz,qb,agg,rsin,ks,pmin
-   !real(rk) :: drag,f_graze,zia,ac_ia,ia_0,ia_b,rnit,skno3_0,sknh4_0,sksil_0,ks_no3,ks_sil,maxg,mort,mort2,crit_melt,lcompp,rpp,rpi,t_sens,nu,md_no3,md_sil,chl2n,sil2n
-  ! REALTYPE:: dic_0, alk_0, dic_sw, alk_sw, dic_ice, alk_ice, ik_diff, ik_on, ice_on, IA_on, tplv, btlv, prop2sw, prop2sw_melt
-   
-
-
-
-
    ! Arrays for work, vertical movement, and cross-boundary fluxes
    REALTYPE,allocatable,dimension(:,:) :: ws
    REALTYPE,allocatable,dimension(:)   :: sfl,bfl,Qsour,Lsour,DefaultRelaxTau,cc_old,curh,curnuh,iweights
@@ -328,22 +311,6 @@ subroutine configure_gotm_fabm(cfg)
                 default=-1, display=display_advanced)
 
    
-!mortenson FABM VARS from gotm.yaml
-
-!gotm_fabm.nml  jpnote some of these might be already existing settings --------- !but there isnt actually anywhere in the yaml where they are coming in from 
-   !branch => cfg%get_child('gotm_fabm_nml', 'setting from the mortenson gotm_fabm.nml')  !jpnote change label/how it's organized in the yaml? 
-   !call branch%get(fabm_calc, 'fabm_calc', 'fabm_calc', default=.false.)
-   !call branch%get(cnpar, 'cnpar', 'cnpar','', default=1.0_rk)  
-   !call branch%get(w_adv_discr, 'w_adv_discr', 'w_adv_discr','', default=6.0_rk)
-   !call branch%get(ode_method, 'ode_method', 'ode_method','', default=1.0_rk)
-   !call branch%get(split_factor, 'split_factor', 'split_factor','', default=1.0_rk)
-   !call branch%get(bioshade_feedback, 'bioshade_feedback', 'bioshade_feedback', default=.true.)
-  ! call branch%get(bioalbedo_feedback, 'bioalbedo_feedback', 'bioalbedo_feedback', default=.true.)
-   !call branch%get(biodrag_feedback, 'biodrag_feedback', 'biodrag_feedback', default=.true.)
-   !call branch%get(repair_state, 'repair_state', 'repair_state', default=.false.)
-   !call branch%get(salinity_relaxation_to_freshwater_flux, 'salinity_relaxation_to_freshwater_flux', 'salinity_relaxation_to_freshwater_flux', default=.false.)
-   !call branch%get(no_precipitation_dilution, 'no_precipitation_dilution', 'no_precipitation_dilution', default=.false.)
-   !call branch%get(save_inputs, 'save_inputs', 'save_inputs', default=.false.)
     
    LEVEL2 'done.'
 
@@ -738,10 +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) 
-   
+   call model%link_scalar(standard_variables%timestep,dt)  !jpnote added 
    
-   !jpnote added 
 
    allocate(Qsour(0:nlev),stat=rc)
    if (rc /= 0) stop 'init_var_gotm_fabm(): Error allocating (Qsour)'

From add99fd372ba4ef974d1b7f75e46a6c75ae7c9cc Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 31 Aug 2021 17:35:10 -0700
Subject: [PATCH 10/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index 54420326..c069e29d 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 544203261f1a9c45a203d7c4b77beb2cac870fba
+Subproject commit c069e29db82da530d9e8abf1bd2d442572e67da4

From cfa77d9d8d6b61224e67ec7ae08604293ac245e5 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Thu, 7 Oct 2021 11:23:12 -0700
Subject: [PATCH 11/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index c069e29d..eda64831 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit c069e29db82da530d9e8abf1bd2d442572e67da4
+Subproject commit eda6483104f76297d38a583373252f63045f711e

From 55f1185dc2312a8782321c7cf648a3dff7a84a4a Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Sun, 10 Oct 2021 14:42:58 -0700
Subject: [PATCH 12/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index eda64831..ca96b301 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit eda6483104f76297d38a583373252f63045f711e
+Subproject commit ca96b3019e7622d52c7fd867f42ec39bf5d64972

From b55b4a2b8daab6aef2734e57fd29ca881540009f Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Sun, 10 Oct 2021 14:50:00 -0700
Subject: [PATCH 13/20] adding fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index ca96b301..17d5bc97 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit ca96b3019e7622d52c7fd867f42ec39bf5d64972
+Subproject commit 17d5bc97fb096f8122b74cd7a6b8158194a6b014

From 9815c713cc58a1f556eb7bb61efea2a52867c912 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Mon, 20 Dec 2021 17:20:10 -0800
Subject: [PATCH 14/20] adding extern/fabm,flexout/stim

---
 extern/fabm    | 2 +-
 extern/flexout | 2 +-
 extern/stim    | 2 +-
 3 files changed, 3 insertions(+), 3 deletions(-)

diff --git a/extern/fabm b/extern/fabm
index 17d5bc97..9e7b2f94 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 17d5bc97fb096f8122b74cd7a6b8158194a6b014
+Subproject commit 9e7b2f94f47b68e1ab5d5b991545e4946eb0b010
diff --git a/extern/flexout b/extern/flexout
index f129ad21..3da32b10 160000
--- a/extern/flexout
+++ b/extern/flexout
@@ -1 +1 @@
-Subproject commit f129ad2181576624ae929331388d0a69bc574f7f
+Subproject commit 3da32b1063e1f5ddfce56edaedbc1b29674f4e84
diff --git a/extern/stim b/extern/stim
index a8d2346f..bd929672 160000
--- a/extern/stim
+++ b/extern/stim
@@ -1 +1 @@
-Subproject commit a8d2346f08b61012d4ba015eb5a9d3d1a86305df
+Subproject commit bd9296724cd2911eb30fc2ef6c9d17dd70b76a5c

From 09169bc9d3aca387575347493577c88416ad9612 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 21 Dec 2021 13:23:42 -0800
Subject: [PATCH 15/20] adding fabm - dms changes

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index 9e7b2f94..350c694d 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 9e7b2f94f47b68e1ab5d5b991545e4946eb0b010
+Subproject commit 350c694dac0a3e7ebf33a7e30ceebabe5e8e4e37

From 5d166ae82cba5014b3274b8c22c538f25a27710e Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 21 Dec 2021 14:09:21 -0800
Subject: [PATCH 16/20] flexout

---
 extern/flexout | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/flexout b/extern/flexout
index 3da32b10..4c141f25 160000
--- a/extern/flexout
+++ b/extern/flexout
@@ -1 +1 @@
-Subproject commit 3da32b1063e1f5ddfce56edaedbc1b29674f4e84
+Subproject commit 4c141f253dbc1e37657fbcf00be0c2da611f9ec0

From f25f6444fc0c73dfab4c3c597f009369db5daf1c Mon Sep 17 00:00:00 2001
From: Karsten Bolding <karsten@bolding-bruggeman.com>
Date: Mon, 7 Dec 2020 10:15:44 +0100
Subject: [PATCH 17/20] added changed from commit1
 2015a027d621296bee2174bd0c8dd847f298e7fb

---
 src/airsea/airsea.F90             | 30 +++++++++++++++++-------------
 src/airsea/airsea_variables.F90   |  3 ++-
 src/airsea/longwave_radiation.F90 | 28 +++++++++++++++++++++++-----
 3 files changed, 42 insertions(+), 19 deletions(-)

diff --git a/src/airsea/airsea.F90 b/src/airsea/airsea.F90
index 34b0e832..d51b0227 100644
--- a/src/airsea/airsea.F90
+++ b/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
@@ -490,12 +491,18 @@ subroutine init_airsea_yaml()
                 default=.false.)
    call branch%get(calc_evaporation, 'calc_evaporation', 'calculate evaporation from meteorological conditions', &
                 default=.false.)
-
-   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', &
+!jpnote: before cherry pick  -- keeping commented because unsure if removing I/0 does anything -- also should check yaml. 
+   !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', &
+               ! 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 
+   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', &
@@ -687,7 +694,7 @@ subroutine post_init_airsea(lat,lon)
             LEVEL4 'using Fairall et. all formulation'
          case default
       end select
-      LEVEL3 'net longwave radiation:'
+      LEVEL3 'longwave radiation:'
       select case (ql%method)
          case(0) ! Read from file instead of calculating
             call register_input(ql)
@@ -1017,9 +1024,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)
       end if
 #if 0
       call airsea_fluxes(fluxes_method,rain_impact,calc_evaporation, &
@@ -1069,10 +1075,8 @@ 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)
    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)
diff --git a/src/airsea/airsea_variables.F90 b/src/airsea/airsea_variables.F90
index 813273a7..d53da035 100644
--- a/src/airsea/airsea_variables.F90
+++ b/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]
diff --git a/src/airsea/longwave_radiation.F90 b/src/airsea/longwave_radiation.F90
index cbf32712..cc33b6bc 100644
--- a/src/airsea/longwave_radiation.F90
+++ b/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,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},
@@ -27,16 +28,25 @@ subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
    IMPLICIT NONE
 !
 ! !INPUT PARAMETERS:
-   integer, intent(in)                 :: method
+   integer, intent(in)                 :: method,type
    REALTYPE, intent(in)                :: dlat,tw,ta,cloud
 !
 ! !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,7 +70,6 @@ subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
 
    REALTYPE                  :: ccf
    REALTYPE                  :: x1,x2,x3
-!
 !EOP
 !-----------------------------------------------------------------------
 !BOC
@@ -68,7 +77,16 @@ subroutine longwave_radiation(method,dlat,tw,ta,cloud,ql)
    ccf= cloud_correction_factor(nint(abs(dlat))+1)
 
    select case(method)
-      case(CLARK)
+      !before
+      !case(CLARK)
+      !after
+      case(from_file)
+         select case(type)
+            case(1)
+            case(2)
+               ql = ql-bolz*emiss*(tw**4)
+         end select
+      case(clark)
 !        Clark et al. (1974) formula.
 !        unit of ea is Pascal, must hPa
 !        Black body defect term, clouds, water vapor correction

From 4db20075e2184e33a08dca6062d08675b1e9a028 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 21 Dec 2021 14:31:42 -0800
Subject: [PATCH 18/20] fabm

---
 extern/fabm | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/extern/fabm b/extern/fabm
index 350c694d..717e63e3 160000
--- a/extern/fabm
+++ b/extern/fabm
@@ -1 +1 @@
-Subproject commit 350c694dac0a3e7ebf33a7e30ceebabe5e8e4e37
+Subproject commit 717e63e3f4922e0975d4eebb20af8998ee194d43

From 2590787c06f1da2ec9ac39b197efcac33a96568a Mon Sep 17 00:00:00 2001
From: Karsten Bolding <karsten@bolding-bruggeman.com>
Date: Wed, 27 Jan 2021 08:35:45 +0100
Subject: [PATCH 19/20] changes from commit 2 22dd70983da5d9e98ce2a1675fdfa8a7b

---
 src/airsea/airsea.F90               | 46 ++++++++++++++++-------------
 src/airsea/airsea_variables.F90     |  1 +
 src/airsea/longwave_radiation.F90   | 21 ++++++-------
 src/gotm/register_all_variables.F90 |  3 +-
 4 files changed, 40 insertions(+), 31 deletions(-)

diff --git a/src/airsea/airsea.F90 b/src/airsea/airsea.F90
index d51b0227..0bb9c2c0 100644
--- a/src/airsea/airsea.F90
+++ b/src/airsea/airsea.F90
@@ -77,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
@@ -342,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,19 +491,25 @@ subroutine init_airsea_yaml()
                 default=.false.)
    call branch%get(calc_evaporation, 'calc_evaporation', 'calculate evaporation from meteorological conditions', &
                 default=.false.)
-!jpnote: before cherry pick  -- keeping commented because unsure if removing I/0 does anything -- also should check yaml. 
-   !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')/))
+
+ !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', &
                ! 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 
-   call branch%get(ql, 'longwave_radiation', 'longwave radiation', 'W/m^2', &
+!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(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)
@@ -695,10 +701,10 @@ subroutine post_init_airsea(lat,lon)
          case default
       end select
       LEVEL3 'longwave radiation:'
-      select case (ql%method)
+      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'
@@ -799,10 +805,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
@@ -1024,8 +1030,8 @@ subroutine flux_from_meteo(jul,secs)
       cloud1 = cloud2
 
       call humidity(hum_method,hum,airp,tw,ta)
-      call longwave_radiation(ql%method,longwave_type, &
-                              dlat,tw_k,ta_k,cloud,ql%value)
+      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, &
@@ -1034,7 +1040,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
@@ -1075,11 +1081,11 @@ subroutine flux_from_meteo(jul,secs)
    end if
 
    call humidity(hum_method,hum%value,airp%value,tw,ta)
-   call longwave_radiation(ql%method,longwave_type, &
-                           dlat,tw_k,ta_k,cloud%value,ql%value)
+   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))
diff --git a/src/airsea/airsea_variables.F90 b/src/airsea/airsea_variables.F90
index d53da035..7773fa70 100644
--- a/src/airsea/airsea_variables.F90
+++ b/src/airsea/airsea_variables.F90
@@ -47,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
diff --git a/src/airsea/longwave_radiation.F90 b/src/airsea/longwave_radiation.F90
index cc33b6bc..579a2430 100644
--- a/src/airsea/longwave_radiation.F90
+++ b/src/airsea/longwave_radiation.F90
@@ -5,7 +5,7 @@
 ! !ROUTINE: Calculate the net longwave radiation \label{sec:back-rad}
 !
 ! !INTERFACE:
-   subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,ql)
+   subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,qlobs,ql)
 !
 ! !DESCRIPTION:
 !
@@ -29,7 +29,7 @@ subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,ql)
 !
 ! !INPUT PARAMETERS:
    integer, intent(in)                 :: method,type
-   REALTYPE, intent(in)                :: dlat,tw,ta,cloud
+   REALTYPE, intent(in)                :: dlat,tw,ta,cloud,qlobs
 !
 ! !OUTPUT PARAMETERS:
    REALTYPE, intent(inout)               :: ql
@@ -39,13 +39,13 @@ subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,ql)
 !
 ! !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)
+   !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, &
@@ -83,8 +83,9 @@ subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,ql)
       case(from_file)
          select case(type)
             case(1)
+               ql=qlobs
             case(2)
-               ql = ql-bolz*emiss*(tw**4)
+               ql = qlobs-bolz*emiss*(tw**4)
          end select
       case(clark)
 !        Clark et al. (1974) formula.
diff --git a/src/gotm/register_all_variables.F90 b/src/gotm/register_all_variables.F90
index 62711e54..44185e96 100644
--- a/src/gotm/register_all_variables.F90
+++ b/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')

From 95dc0967aaf6d96a104bfa41b6b4d6da51ea2811 Mon Sep 17 00:00:00 2001
From: Julia Putko <jputko@aeolus.seos.uvic.ca>
Date: Tue, 21 Dec 2021 15:19:48 -0800
Subject: [PATCH 20/20] additional changes regarding longwave: FILE, ql_

---
 src/airsea/airsea.F90             | 8 +++++---
 src/airsea/longwave_radiation.F90 | 6 +++---
 2 files changed, 8 insertions(+), 6 deletions(-)

diff --git a/src/airsea/airsea.F90 b/src/airsea/airsea.F90
index 0bb9c2c0..af4f8d75 100644
--- a/src/airsea/airsea.F90
+++ b/src/airsea/airsea.F90
@@ -775,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, &
diff --git a/src/airsea/longwave_radiation.F90 b/src/airsea/longwave_radiation.F90
index 579a2430..c98e086f 100644
--- a/src/airsea/longwave_radiation.F90
+++ b/src/airsea/longwave_radiation.F90
@@ -24,7 +24,7 @@ subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,qlobs,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:
@@ -80,14 +80,14 @@ subroutine longwave_radiation(method,type,dlat,tw,ta,cloud,qlobs,ql)
       !before
       !case(CLARK)
       !after
-      case(from_file)
+      case(FILE)
          select case(type)
             case(1)
                ql=qlobs
             case(2)
                ql = qlobs-bolz*emiss*(tw**4)
          end select
-      case(clark)
+      case(CLARK)
 !        Clark et al. (1974) formula.
 !        unit of ea is Pascal, must hPa
 !        Black body defect term, clouds, water vapor correction