Option to not cycle met drivers for CO2, plus multiple minor bug fixes

BRAMS/src/ed2/edcp_driver.F90

@@ -66,7 +66,8 @@ subroutine ed_coup_driver()
    integer                     :: jd2
    integer                     :: ierr
    integer                     :: igr
-   integer                     :: ping 
+   integer                     :: ping
+   logical                     :: new_day
    real                        :: wtime1
    real                        :: wtime2
    real                        :: wtime_start     ! wall time

ED/src/driver/ed_driver.F90

@@ -182,7 +182,8 @@ subroutine ed_driver()
 
 
    !---------------------------------------------------------------------------------------!
-   if (trim(runtype) == 'HISTORY' ) then
+   select case (trim(runtype))
+   case ('HISTORY')
       !------------------------------------------------------------------------------------!
       !      Initialize the model state as a replicate image of a previous  state.         !
       !------------------------------------------------------------------------------------!
@@ -207,15 +208,15 @@ subroutine ed_driver()
       if (nnodetot /= 1 ) call MPI_Barrier(MPI_COMM_WORLD,ierr)
 #endif
       !------------------------------------------------------------------------------------!
-   else
+   case default
 
       !------------------------------------------------------------------------------------!
       !      Initialize state properties of polygons/sites/patches/cohorts.                !
       !------------------------------------------------------------------------------------!
       if (mynum == nnodetot) write (unit=*,fmt='(a)') ' [+] Load_Ecosystem_State...'
       call load_ecosystem_state()
       !------------------------------------------------------------------------------------!
-   end if
+   end select
 
    !---------------------------------------------------------------------------------------!
    !      In case the runs is going to produce detailed output, we eliminate all patches   !
@@ -291,11 +292,14 @@ subroutine ed_driver()
    !      Initialise some derived variables.  Skip this in case the simulation is resuming !
    ! from HISTORY.                                                                         !
    !---------------------------------------------------------------------------------------!
-   if (trim(runtype) /= 'HISTORY' ) then
+   select case (trim(runtype))
+   case ('HISTORY')
+      continue
+   case default
       do ifm=1,ngrids
          call update_derived_props(edgrid_g(ifm))
       end do
-   end if
+   end select
    !---------------------------------------------------------------------------------------!
 
 
@@ -304,11 +308,14 @@ subroutine ed_driver()
    !      Initialise drought phenology.  This should be done after the soil moisture has   !
    ! been set up.                                                                          !
    !---------------------------------------------------------------------------------------!
-   if (runtype /= 'HISTORY') then
+   select case (trim(runtype))
+   case ('HISTORY')
+      continue
+   case default
       do ifm=1,ngrids
          call first_phenology(edgrid_g(ifm))
       end do
-   end if
+   end select
    !---------------------------------------------------------------------------------------!
 
 

ED/src/dynamics/disturbance.f90

@@ -1563,7 +1563,7 @@ subroutine site_disturbance_rates(year, cgrid)
                !---------------------------------------------------------------------------!
             case default
                !------ Read anthropogenic disturbance from external data set. -------------!
-               if (clutime%landuse(12) < 0 .or. clutime%landuse(14) < 0) then
+               if (clutime%landuse(12) < 0. .or. clutime%landuse(14) < 0.) then
                   find_target                         = .true.
                   cpoly%primary_harvest_target  (isi) = 0.
                   cpoly%secondary_harvest_target(isi) = 0.
@@ -2265,6 +2265,9 @@ subroutine increment_patch_vars(csite,np,cp,area_fac,cb_enthalpy,can_exner,cb_ma
       csite%fmean_sfcw_mass      (np) = csite%fmean_sfcw_mass      (np)                    &
                                       + csite%fmean_sfcw_mass      (cp)                    &
                                       * area_fac
+      csite%fmean_snowfac        (np) = csite%fmean_snowfac        (np)                    &
+                                      + csite%fmean_snowfac        (cp)                    &
+                                      * area_fac
       csite%fmean_rshort_gnd     (np) = csite%fmean_rshort_gnd     (np)                    &
                                       + csite%fmean_rshort_gnd     (cp)                    &
                                       * area_fac
@@ -2487,6 +2490,9 @@ subroutine increment_patch_vars(csite,np,cp,area_fac,cb_enthalpy,can_exner,cb_ma
          csite%dmean_sfcw_fliq      (    np) = csite%dmean_sfcw_fliq      (    np)         &
                                              + csite%dmean_sfcw_fliq      (    cp)         &
                                              * area_fac
+         csite%dmean_snowfac        (    np) = csite%dmean_snowfac        (    np)         &
+                                             + csite%dmean_snowfac        (    cp)         &
+                                             * area_fac
          csite%dmean_rshort_gnd     (    np) = csite%dmean_rshort_gnd     (    np)         &
                                              + csite%dmean_rshort_gnd     (    cp)         &
                                              * area_fac
@@ -2770,6 +2776,9 @@ subroutine increment_patch_vars(csite,np,cp,area_fac,cb_enthalpy,can_exner,cb_ma
          csite%mmean_sfcw_fliq      (    np) = csite%mmean_sfcw_fliq      (    np)         &
                                              + csite%mmean_sfcw_fliq      (    cp)         &
                                              * area_fac
+         csite%mmean_snowfac        (    np) = csite%mmean_snowfac        (    np)         &
+                                             + csite%mmean_snowfac        (    cp)         &
+                                             * area_fac
          csite%mmean_rshort_gnd     (    np) = csite%mmean_rshort_gnd     (    np)         &
                                              + csite%mmean_rshort_gnd     (    cp)         &
                                              * area_fac
@@ -3065,6 +3074,9 @@ subroutine increment_patch_vars(csite,np,cp,area_fac,cb_enthalpy,can_exner,cb_ma
          csite%qmean_sfcw_fliq      (  :,np) = csite%qmean_sfcw_fliq      (  :,np)         &
                                              + csite%qmean_sfcw_fliq      (  :,cp)         &
                                              * area_fac
+         csite%qmean_snowfac        (  :,np) = csite%qmean_snowfac        (  :,np)         &
+                                             + csite%qmean_snowfac        (  :,cp)         &
+                                             * area_fac
          csite%qmean_soil_energy    (:,:,np) = csite%qmean_soil_energy    (:,:,np)         &
                                              + csite%qmean_soil_energy    (:,:,cp)         &
                                              * area_fac

ED/src/dynamics/euler_driver.f90

@@ -97,7 +97,7 @@ subroutine euler_timestep(cgrid)
             !     Update the monthly rainfall.                                             !
             !------------------------------------------------------------------------------!
             imon                             = current_time%month
-            cpoly%avg_monthly_pcpg(imon,isi) = cpoly%avg_monthly_pcpg(imon,isi)            &
+            cpoly%avg_monthly_accp(imon,isi) = cpoly%avg_monthly_accp(imon,isi)            &
                                              + cmet%pcpg * dtlsm
             !------------------------------------------------------------------------------!
 

ED/src/dynamics/farq_katul.f90

@@ -247,9 +247,6 @@ subroutine katul_lphys(ib,can_prss,can_rhos,can_shv,can_co2,ipft,leaf_par,leaf_t
 
       !update photosynthetic parameters with water stress
       select case (h2o_plant_lim)
-      case (0,1,2,3)
-          ! use fsw to account for water stress in photosyn_driv
-          water_stress_factor = 1.
       case (4)
           ! leaf water potential will influence stomata optimization 
           ! at two different scales
@@ -268,6 +265,9 @@ subroutine katul_lphys(ib,can_prss,can_rhos,can_shv,can_co2,ipft,leaf_par,leaf_t
                                         1. / (1. +      &
                                 0.1 * (leaf_psi / leaf_psi_tlp(ipft)) ** 6.0)))
           lambda8 = lambda8 * dble(exp(stoma_beta(ipft) * dmax_leaf_psi))
+      case default
+          ! use fsw to account for water stress in photosyn_driv
+          water_stress_factor = 1.
       end select
 
       !thispft(ib)%vm0  = thispft(ib)%vm0 * dble(water_stress_factor)
@@ -848,8 +848,8 @@ subroutine photosynthesis_stomata_solver8(ib,gsc,limit_case,
         real(kind=8)                :: k1,k2        !! Variable used in photosynthesis equation
         real(kind=8)                :: a,b,c        !! Coefficients of the quadratic equation to solve ci
         real(kind=8)                :: rad          !! sqrt(b2-4ac)
-        real(kind=8)                :: dbdg,dcdg    !! derivatives of b,c wrt. gsc
-
+        real(kind=8)                :: dbdg         !! derivatives of b wrt. gsc
+        real(kind=8)                :: dcdg         !! derivatives of c wrt. gsc
         !------------------------------------------------------------------------------------!
 
 

ED/src/dynamics/fire.f90

@@ -62,7 +62,7 @@ subroutine fire_frequency(cgrid)
       real                           :: fuel
       real                           :: ignition_rate
       real                           :: mean_fire_intensity
-      real                           :: sum_pcpg
+      real                           :: sum_accp
       logical                        :: people_around
       !------------------------------------------------------------------------------------!
 
@@ -97,8 +97,8 @@ subroutine fire_frequency(cgrid)
             !     Find the total rainfall of the past year and reset the counter for this  !
             ! month.                                                                       !
             !------------------------------------------------------------------------------!
-            sum_pcpg                         = sum(cpoly%avg_monthly_pcpg(:,isi))
-            cpoly%avg_monthly_pcpg(imon,isi) = 0.
+            sum_accp                         = sum(cpoly%avg_monthly_accp(:,isi))
+            cpoly%avg_monthly_accp(imon,isi) = 0.
             !------------------------------------------------------------------------------!
 
 

ED/src/dynamics/forestry.f90

@@ -142,6 +142,10 @@ subroutine find_lambda_harvest(cpoly,isi,onsp,lambda_harvest)
          ! croplands and pastures.                                                         !
          !---------------------------------------------------------------------------------!
          select case(ilu)
+         case (1,8)
+            !---- Pasture or cropland. Do nothing. ----------------------------------------!
+            continue
+            !------------------------------------------------------------------------------!
          case (2:7)
             hcoh_loop: do ico=1,cpatch%ncohorts
                ipft = cpatch%pft(ico)

ED/src/dynamics/heun_driver.f90

@@ -98,7 +98,7 @@ subroutine heun_timestep(cgrid)
             !     Update the monthly rainfall.                                             !
             !------------------------------------------------------------------------------!
             imon                             = current_time%month
-            cpoly%avg_monthly_pcpg(imon,isi) = cpoly%avg_monthly_pcpg(imon,isi)            &
+            cpoly%avg_monthly_accp(imon,isi) = cpoly%avg_monthly_accp(imon,isi)            &
                                              + cmet%pcpg * dtlsm
             !------------------------------------------------------------------------------!
 

ED/src/dynamics/hybrid_driver.f90

@@ -107,7 +107,7 @@ subroutine hybrid_timestep(cgrid)
            !     Update the monthly rainfall.                                    !
            !---------------------------------------------------------------------!
            imon                             = current_time%month
-           cpoly%avg_monthly_pcpg(imon,isi) = cpoly%avg_monthly_pcpg(imon,isi)   &
+           cpoly%avg_monthly_accp(imon,isi) = cpoly%avg_monthly_accp(imon,isi)   &
                                             + cmet%pcpg * dtlsm
            !---------------------------------------------------------------------!
 

ED/src/dynamics/reproduction.f90

@@ -179,6 +179,7 @@ subroutine reproduction_driver(cgrid,month,veget_dyn_on)
 
          !----- The big loops start here. -------------------------------------------------!
          polyloop: do ipy = 1,cgrid%npolygons
+            cpoly => cgrid%polygon(ipy)
 
             !------------------------------------------------------------------------------!
             !     Check whether this is late spring/early summer.  This is needed for      !
@@ -187,8 +188,8 @@ subroutine reproduction_driver(cgrid,month,veget_dyn_on)
             !------------------------------------------------------------------------------!
             late_spring = (cgrid%lat(ipy) >= 0.0 .and. month == 6) .or.                    &
                           (cgrid%lat(ipy) < 0.0 .and. month == 12)
+            !------------------------------------------------------------------------------!
 
-            cpoly => cgrid%polygon(ipy)
             siteloop_sort: do isi = 1,cpoly%nsites
                csite => cpoly%site(isi)
 
@@ -1286,8 +1287,8 @@ subroutine seed_dispersal(cpoly,late_spring)
                   ! of them will land again in this patch, and we correct for this further !
                   ! down.                                                                  !
                   !------------------------------------------------------------------------!
-                  csite%cbudget_seedrain(donpa) = csite%cbudget_seedrain(donpa)            &
-                                                - bseed_maygo * frqsumi
+                  donsite%cbudget_seedrain(donpa) = donsite%cbudget_seedrain(donpa)        &
+                                                  - bseed_maygo * frqsumi
                   !------------------------------------------------------------------------!
 
 
@@ -1314,7 +1315,7 @@ subroutine seed_dispersal(cpoly,late_spring)
                   ! (4) RPY = DPA * AD * AR           (1->3)                               !
                   ! (5) RPA = DPA * AD                (4->2, regardless of the patch)      !
                   !------------------------------------------------------------------------!
-                  bseed_xpatch = bseed_maygo * csite%area(donpa) * cpoly%area(donsi)
+                  bseed_xpatch = bseed_maygo * donsite%area(donpa) * cpoly%area(donsi)
                   !------------------------------------------------------------------------!
 
 
@@ -1343,8 +1344,8 @@ subroutine seed_dispersal(cpoly,late_spring)
                         ! subtracted all the non-local dispersal outside the receptor site !
                         ! loop.                                                            !
                         !------------------------------------------------------------------!
-                        csite%cbudget_seedrain(recpa) = csite%cbudget_seedrain(recpa)      &
-                                                      + bseed_xpatch * frqsumi
+                        recsite%cbudget_seedrain(recpa) = recsite%cbudget_seedrain(recpa)  &
+                                                        + bseed_xpatch * frqsumi
                         !------------------------------------------------------------------!
 
 

ED/src/dynamics/rk4_copy_patch.f90

@@ -2092,6 +2092,8 @@ subroutine initp2modelp(hdid,initp,csite,ipa,nighttime,wbudget_loss2atm,ebudget_
                                     + csite%ground_shv        (ipa) * dtlsm_o_frqsum
       csite%fmean_can_ggnd    (ipa) = csite%fmean_can_ggnd    (ipa)                        &
                                     + csite%ggnet             (ipa) * dtlsm_o_frqsum
+      csite%fmean_snowfac     (ipa) = csite%fmean_snowfac     (ipa)                        &
+                                    + csite%snowfac           (ipa) * dtlsm_o_frqsum
       !------------------------------------------------------------------------------------!
       !       Snow/pounding layers.  We keep track of the total, not individual layers.    !
       ! Energy will be integrated as an extensive variable, we will convert it by the      !

ED/src/dynamics/rk4_driver.F90

@@ -112,7 +112,7 @@ subroutine rk4_timestep(cgrid)
             !     Update the monthly rainfall.                                             !
             !------------------------------------------------------------------------------!
             imon                             = current_time%month
-            cpoly%avg_monthly_pcpg(imon,isi) = cpoly%avg_monthly_pcpg(imon,isi)            &
+            cpoly%avg_monthly_accp(imon,isi) = cpoly%avg_monthly_accp(imon,isi)            &
                                              + cmet%pcpg * dtlsm
             !------------------------------------------------------------------------------!
 

ED/src/dynamics/vegetation_dynamics.f90

@@ -34,7 +34,7 @@ subroutine veg_dynamics_driver(new_month,new_year,gr_tfact0,veget_dyn_on)
                                       , yr_day                        ! ! intent(in)
       use mem_polygons         , only : maxpatch                      ! ! intent(in)
       use average_utils        , only : normalize_ed_today_vars       & ! sub-routine
-                                      , normalize_ed_todaynpp_vars    & ! sub-routine
+                                      , copy_today_to_dmean_vars      & ! sub-routine
                                       , zero_ed_today_vars            ! ! sub-routine
       use canopy_radiation_coms, only : ihrzrad                       ! ! intent(in)
       use hrzshade_utils       , only : split_hrzshade                & ! sub-routine
@@ -129,7 +129,7 @@ subroutine veg_dynamics_driver(new_month,new_year,gr_tfact0,veget_dyn_on)
 
 
          !------  update dmean and mmean values for NPP allocation terms ------------------!
-         call normalize_ed_todayNPP_vars(cgrid)
+         call copy_today_to_dmean_vars(cgrid)
          !---------------------------------------------------------------------------------!
 
 

ED/src/init/ed_params.f90

@@ -525,7 +525,7 @@ subroutine init_decomp_params()
    rh0          = 0.700 ! 0.701 ! 0.425
    rh_q10       = 1.500 ! 1.500 ! 1.893
    rh_p_smoist  = 1.600 ! 0.836 ! 0.606
-   rh_p_oxygen  = 0.600 ! 0.404 ! 0.164
+   rh_p_oxygen  = 0.450 ! 0.404 ! 0.164
    !---------------------------------------------------------------------------------------!