Merge pull request #22 from ACME-Climate/bishtgautam/fixes-for-gcc5

Fixes for gcc5

Author: jinyun1tang

src/betr/betr_core/TracerParamsMod.F90

@@ -570,7 +570,6 @@ subroutine calc_bunsen_coeff(bounds, lbj, ubj, jtops, numf, filter, &
            c = filter(fc)
            if(n>=jtops(c))then
              bunsencef_col(c,n, k)= henrycef_col(c,n,k)*t_soisno(c,n)/12.2_r8
-!             if(n==1)print*,k,betrtracer_vars%tracernames(trcid),henrycef_col(c,n,k)
              !add the pH effect for tracers that can exist in multiple aqueous phases
              if(is_h2o(trcid))then
                !for water isotopes
@@ -1376,8 +1375,6 @@ subroutine calc_tracer_infiltration(bounds, jtops, numf, filter, bunsencef_topso
            tracer_flx_infl(c,j) = bunsencef_topsoi(c,betrtracer_vars%volatilegroupid(j)) * &
                 tracerboundarycond_vars%tracer_gwdif_concflux_top_col(c,1,j) * qflx_adv(c,0)
 
-!         print*,'trcflx',j,betrtracer_vars%tracernames(j),tracer_flx_infl(c,j),&
-!            bunsencef_topsoi(c,betrtracer_vars%volatilegroupid(j))
          else
            tracer_flx_infl(c,j) = 0._r8
          endif
@@ -1691,7 +1688,6 @@ subroutine calc_aerecond(bounds, col, pft, num_soilp, filter_soilp, jwt, &
     do fp = 1, num_soilp
       p = filter_soilp (fp)
       c = pft%column(p)
-
       ! Calculate aerenchyma diffusion
       if (j > jwt(c) .and. t_soisno(c,j) > tfrz .and. pftvarcon%is_grass_patch(pft%itype(p))) then
         ! Attn EK: This calculation of aerenchyma properties is very uncertain. Let's check in once all

src/betr/betr_core/TracerStateType.F90

@@ -319,7 +319,6 @@ function int_mass_mobile_col(this, lbj, ubj, c, j, dz, bstatus)result(int_mass)
   call bstatus%reset()
 
   int_mass = dot_sum(this%tracer_conc_mobile_col(c,lbj:ubj,j), dz, bstatus)
-
   end function int_mass_mobile_col
 
   !-----------------------------------------------------------------------
@@ -431,7 +430,6 @@ subroutine retrieve_hist(this, bounds, lbj, ubj, state_2d, state_1d, betrtracer_
     state_1d(begc:endc, addone(idtemp1d))=this%tracer_soi_molarmass_col(begc:endc, jj)
     state_1d(begc:endc, addone(idtemp1d))=this%end_tracer_molarmass_col(begc:endc, jj)
   enddo
-
   end associate
   end subroutine retrieve_hist
 

src/betr/betr_dtype/BeTR_biogeophysInputType.F90

@@ -15,8 +15,6 @@ module BeTR_biogeophysInputType
     !waterstate
     real(r8), pointer :: h2osoi_liq_col(:,:)    => null()    !liquid water (kg/m2) (new) (-nlevsno+1:nlevgrnd)
     real(r8), pointer :: h2osoi_ice_col(:,:)    => null()    !ice lens (kg/m2) (new) (-nlevsno+1:nlevgrnd)
-    real(r8), pointer :: h2osoi_liq_old(:,:)    => null()    !liquid water (kg/m2) (old) (-nlevsno+1:nlevgrnd)
-    real(r8), pointer :: h2osoi_ice_old(:,:)    => null()    !ice lens (kg/m2) (old) (-nlevsno+1:nlevgrnd)
     real(r8), pointer :: h2osoi_liqvol_col(:,:) => null()   !volumetric liquid water content
     real(r8), pointer :: h2osoi_icevol_col(:,:) => null()   !volumetric ice water content
     real(r8), pointer :: h2osoi_vol_col(:,:)    => null()    !volumetric water content, total
@@ -182,8 +180,6 @@ subroutine InitAllocate(this, bounds)
   allocate (this%finundated_col(    begc:endc         ) ) ! fraction of column that is inundated, this is for bgc caclulation in betr
   allocate (this%h2osoi_liq_col(    begc:endc,lbj:ubj ) ) !liquid water (kg/m2) (new) (-nlevsno+1:nlevgrnd)
   allocate (this%h2osoi_ice_col(    begc:endc,lbj:ubj ) ) !ice lens (kg/m2) (new) (-nlevsno+1:nlevgrnd)
-  allocate (this%h2osoi_liq_old(    begc:endc,lbj:ubj ) ) !liquid water (kg/m2) (old) (-nlevsno+1:nlevgrnd)
-  allocate (this%h2osoi_ice_old(    begc:endc,lbj:ubj ) ) !ice lens (kg/m2) (old) (-nlevsno+1:nlevgrnd)
   allocate (this%h2osoi_liqvol_col( begc:endc,lbj:ubj ) ) !volumetric liquid water content
   allocate (this%h2osoi_icevol_col( begc:endc,lbj:ubj ) ) !volumetric ice water content
   allocate (this%h2osoi_vol_col(    begc:endc,lbj:ubj ) ) !volumetric water content, total

src/betr/betr_main/BetrBGCMod.F90

@@ -617,7 +617,6 @@ subroutine tracer_gw_transport(betr_time, bounds, lbj, ubj, col, pft, jtops, num
          betrtracer_vars,                                       &
          tracercoeff_vars, tracerstate_vars,bstatus)
     if(bstatus%check_status())return
-!x    print*,'do_transport'
     !do diffusive and advective transport, assuming aqueous and gaseous phase are in equilbrium
     do kk = 1 , 2
        if (transp_pathway(kk) == diffusion_scheme .and. diffusion_on) then
@@ -734,7 +733,6 @@ subroutine do_tracer_advection(betr_time, bounds, lbj, ubj, col, pft, &
     character(len=255)   :: subname = 'do_tracer_advection'
     character(len=betr_errmsg_len) :: msg
 
-!x    print*,'do_advection'
     call bstatus%reset()
     SHR_ASSERT_ALL((ubound(jtops) == (/bounds%endc/))      , errMsg(mod_filename,__LINE__),bstatus)
     if(bstatus%check_status())return
@@ -884,7 +882,7 @@ subroutine do_tracer_advection(betr_time, bounds, lbj, ubj, col, pft, &
                  tracer_conc_mobile_col(bounds%begc:bounds%endc, lbj:ubj,adv_trc_group(1:ntrcs)), &
                  trc_conc_out(:,:,1:ntrcs),                                                       &
                  leaching_mass(bounds%begc:bounds%endc,1:ntrcs), seep_mass(bounds%begc:bounds%endc, 1:ntrcs))
-!x            print*,'semilag',bstatus%check_status(),trim(bstatus%print_msg())
+
             if(bstatus%check_status())return
             !do soil-root tracer exchange
             do k = 1, ntrcs
@@ -917,7 +915,6 @@ subroutine do_tracer_advection(betr_time, bounds, lbj, ubj, col, pft, &
                endif
             enddo
 
-            !print*,'!do error budget and tracer flux update'
             do k = 1, ntrcs
                trcid = adv_trc_group(k)
                do fc = 1, num_soilc
@@ -1250,8 +1247,9 @@ subroutine do_tracer_gw_diffusion(bounds, lbj, ubj, jtops, num_soilc, filter_soi
                                 diff_surf(c,k) * dtime_loc(c)
                         endif
                      else
-                        write(msg,*) 'mass bal error dif '//trim(tracernames(trcid))//new_line('A'), 'mass 0/1',mass0, mass1,'col=',c, &
-                               new_line('A')//'err=', err_tracer(c,k), 'dmass=',dmass(c,k), ' dif=', diff_surf(c,k)*dtime_loc(c), &
+                        write(msg,*) 'mass bal error dif '//trim(tracernames(trcid))//new_line('A'), 'mass 0/1',mass0, &
+                             mass1,'col=',c, &
+                             new_line('A')//'err=', err_tracer(c,k), 'dmass=',dmass(c,k), ' dif=', diff_surf(c,k)*dtime_loc(c), &
                              ' prod=',dot_sum(x=local_source(c,jtops(c):ubj,k),y=dz(c,jtops(c):ubj),bstatus=bstatus)*dtime_loc(c)
 
                         if(bstatus%check_status())return

src/betr/betr_main/TracerBalanceMod.F90

@@ -142,7 +142,6 @@ subroutine betr_tracer_massbalance_check(betr_time, bounds, col,  numf, filter,
               else
                  err_rel = errtracer(c,kk)/max(abs(beg_tracer_molarmass(c,kk)),abs(end_tracer_molarmass(c,kk)))
               endif
-
               if(abs(err_rel)>err_min_rel .and. do_betr_output)then
                  write(msg,*)'error exceeds the tolerance for tracer '//tracernames(kk), &
                       new_line('A'),'err=',errtracer(c,kk), ' col=',c, &
@@ -154,6 +153,8 @@ subroutine betr_tracer_massbalance_check(betr_time, bounds, col,  numf, filter,
                       new_line('A'),errMsg(mod_filename, __LINE__)
                  call tracerflux_vars%flux_display(c,kk,betrtracer_vars, msg1)
                  msg = trim(msg)//new_line('A')//trim(msg1)
+                 print*,trim(msg)
+                 pause
                  call betr_status%set_msg(msg=msg, err=-1)
                  return
               endif
@@ -229,6 +230,7 @@ subroutine betr_tracer_mass_summary(bounds, col, lbj, ubj, numf, filter, betrtra
                       frozenid(jj),dz(c,1:nlevtrc_soil),betr_status)
                  if(betr_status%check_status())return
               endif
+
            enddo
         enddo
       end associate

src/betr/betr_math/MathfuncMod.F90

@@ -31,7 +31,7 @@ module MathfuncMod
   public :: pd_decomp
   public :: num2str
   public :: fpmax
-
+  public :: bisnan
   interface cumsum
      module procedure cumsum_v, cumsum_m
   end interface cumsum
@@ -586,4 +586,17 @@ function fpmax(inval)result(ans)
   ans = max(inval, 0._r8)
   return
   end function fpmax  
+
+  !-------------------------------------------------------------------------------
+  function bisnan(inval)result(ans)
+  
+  !DESCRIPTION
+  !determine if the variable is nan
+  implicit none
+  real(r8), intent(in) :: inval
+
+  logical :: ans
+
+  ans = (inval/=inval)
+  end function bisnan 
 end module MathfuncMod

src/betr/betr_rxns/MockBGCReactionsType.F90

@@ -250,7 +250,6 @@ subroutine set_boundary_conditions(this, bounds, num_soilc, filter_soilc, dz_top
       do fc = 1, num_soilc
          c = filter_soilc(fc)
          irt = 1.e3_r8/(forc_tbot(c)*rgas)
-
          !eventually, the following code will be implemented using polymorphism
          tracerboundarycond_vars%tracer_gwdif_concflux_top_col(c,1:2,betrtracer_vars%id_trc_n2)   = forc_pbot(c)*0.78084_r8*irt  !mol m-3, contant boundary condition, as concentration
          tracerboundarycond_vars%tracer_gwdif_concflux_top_col(c,1:2,betrtracer_vars%id_trc_o2)   = forc_pbot(c)*0.20946_r8*irt  !mol m-3, contant boundary condition, as concentration

src/driver/alm/ALMbetrNLMod.F90

@@ -40,7 +40,7 @@ subroutine betr_readNL(NLFilename)
     character(len=1), parameter  :: quote = ''''
     namelist / betr_inparm / reaction_method, &
       advection_on, diffusion_on, reaction_on, ebullition_on
-
+    logical :: esm_on
     ! ----------------------------------------------------------------------
     ! Read namelist from standard input.
     ! ----------------------------------------------------------------------
@@ -49,8 +49,8 @@ subroutine betr_readNL(NLFilename)
     advection_on    = .true.
     diffusion_on    = .true.
     reaction_on     = .true.
-    ebullition_on   =.true.
-
+    ebullition_on   = .true.
+    esm_on          = .true.
     if ( masterproc )then
 
        unitn = getavu()
@@ -75,6 +75,7 @@ subroutine betr_readNL(NLFilename)
 
     write(betr_namelist_buffer,*) '&betr_parameters'//new_line('A'), &
       ' reaction_method='//quote//trim(reaction_method)//quote//new_line('A'), &
+      ' esm_on=',trim(log2str(esm_on)),new_line('A'),&
       ' advection_on=',trim(log2str(advection_on)),new_line('A'), &
       ' diffusion_on=',trim(log2str(diffusion_on)),new_line('A'), &
       ' reaction_on=',trim(log2str(reaction_on)),new_line('A'), &

src/driver/alm/BeTRSimulationALM.F90

@@ -34,12 +34,12 @@ module BeTRSimulationALM
      procedure, public :: StepWithDrainage          => ALMStepWithDrainage
      procedure, public :: SetBiophysForcing         => ALMSetBiophysForcing
      !unique subroutines
-     procedure, public :: DiagnoseDtracerFreezeThaw => ALMDiagnoseDtracerFreezeThaw
      procedure, public :: CalcDewSubFlux            => ALMCalcDewSubFlux
      procedure, public :: SoilFluxStateRecv         => ALMBetrSoilFluxStateRecv
      procedure, public :: CalcSmpL                  => ALMCalcSmpL
      procedure, public :: PlantSoilBGCSend          => ALMBetrPlantSoilBGCSend
      procedure, public :: PlantSoilBGCRecv          => ALMBetrPlantSoilBGCRecv
+     procedure, public :: set_active                => ALMset_active
   end type betr_simulation_alm_type
 
   public :: create_betr_simulation_alm
@@ -207,6 +207,24 @@ subroutine ALMStepWithoutDrainage(this, bounds,  col, pft)
       call endrun(msg=this%bsimstatus%print_msg())
   end subroutine ALMStepWithoutDrainage
 
+  !---------------------------------------------------------------------------------
+  subroutine ALMset_active(this,bounds,col)
+  
+  !
+  !DESCRIPTION
+  !activate columuns that are active in alm
+  use ColumnType     , only : column_type
+  implicit none
+  ! !ARGUMENTS:
+  class(betr_simulation_alm_type) , intent(inout) :: this
+  type(bounds_type)               , intent(in)    :: bounds
+  type(column_type)               , intent(in)    :: col ! column type
+  
+  integer :: c
+  do c = bounds%begc, bounds%endc
+    this%active_col(c) = (this%active_col(c) .and. col%active(c))
+  enddo
+  end subroutine ALMset_active  
   !---------------------------------------------------------------------------------
   subroutine ALMStepWithDrainage(this, bounds,  col)
    !DESCRIPTION
@@ -523,43 +541,6 @@ subroutine ALMBetrPlantSoilBGCRecv(this, bounds, num_soilc,  filter_soilc,&
   end subroutine ALMBetrPlantSoilBGCRecv
   !------------------------------------------------------------------------
 
-  subroutine ALMDiagnoseDtracerFreezeThaw(this, bounds, num_nolakec, filter_nolakec, col, lun)
-  !
-  ! DESCRIPTION
-  ! aqueous tracer partition based on freeze-thaw
-  !
-  ! USES
-  use ColumnType            , only : column_type
-  use LandunitType          , only : landunit_type
-  use WaterStateType        , only : waterstate_type
-  implicit none
-  !
-  ! Arguments
-  class(betr_simulation_alm_type), intent(inout)   :: this
-  type(bounds_type)     , intent(in) :: bounds
-  integer               , intent(in) :: num_nolakec                        ! number of column non-lake points in column filter
-  integer               , intent(in) :: filter_nolakec(:)                  ! column filter for non-lake points
-  type(landunit_type)   , intent(in) :: lun
-!  type(waterstate_type), intent(in) :: waterstate_vars
-  type(column_type)     , intent(in) :: col                                ! column type
-
-  !temporary variables
-  type(betr_bounds_type)     :: betr_bounds
-  integer :: fc, c
-
-  call this%BeTRSetBounds(betr_bounds)
-
-  call this%BeTRSetcps(bounds, col)
-
-  do fc = 1, num_nolakec
-    c = filter_nolakec(fc)
-    if(.not. this%active_col(c))cycle
-    call this%betr(c)%diagnose_dtracer_freeze_thaw(betr_bounds, this%num_soilc, this%filter_soilc,  &
-      this%biophys_forc(c))
-  enddo
-  end subroutine ALMDiagnoseDtracerFreezeThaw
-
-  !------------------------------------------------------------------------
   subroutine ALMCalcDewSubFlux(this,  &
        bounds, col, num_hydrologyc, filter_soilc_hydrologyc)
    !DESCRIPTION
@@ -716,6 +697,7 @@ subroutine ALMSetBiophysForcing(this, bounds, col, pft, carbonflux_vars, waterst
   type(soilstate_type)        , optional, intent(in) :: soilstate_vars
 
   integer :: p, pi, c
+  integer :: npft_loc
 
   call this%BeTRSetBiophysForcing(bounds, col, pft, 1, nlevsoi, carbonflux_vars, waterstate_vars, &
       waterflux_vars, temperature_vars, soilhydrology_vars, atm2lnd_vars, canopystate_vars, &
@@ -724,15 +706,23 @@ subroutine ALMSetBiophysForcing(this, bounds, col, pft, carbonflux_vars, waterst
 
   !the following will be ALM specific
   !big leaf model
+  !set profiles autotrohpic respiration
   do c = bounds%begc, bounds%endc
-    do pi = 1, betr_maxpatch_pft
-      if (pi <= col%npfts(c)) then
-        p = col%pfti(c) + pi - 1
-        if (pft%active(p)) then
-          this%biophys_forc(c)%rr_patch(pi,1:nlevsoi) = carbonflux_vars%rr_patch(p) !* root_prof(p,1:nlevsoi)
+    npft_loc = ubound(carbonflux_vars%rr_patch,1)-lbound(carbonflux_vars%rr_patch,1)+1
+    if(npft_loc /= col%npfts(c) .and. col%pfti(c) /= lbound(carbonflux_vars%rr_patch,1)) then
+      do pi = 1, betr_maxpatch_pft
+        this%biophys_forc(c)%rr_patch(pi,1:nlevsoi) = 0._r8
+      enddo
+    else
+      do pi = 1, betr_maxpatch_pft
+        if (pi <= col%npfts(c)) then
+          p = col%pfti(c) + pi - 1
+          if (pft%active(p)) then
+            this%biophys_forc(c)%rr_patch(pi,1:nlevsoi) = carbonflux_vars%rr_patch(p) !* root_prof(p,1:nlevsoi)
+          endif
         endif
-      endif
-    enddo
+      enddo
+    endif
   enddo
   !dvgm
   !