soil_nutcarb_init.f90 Source File


This file depends on

sourcefile~~soil_nutcarb_init.f90~~EfferentGraph sourcefile~soil_nutcarb_init.f90 soil_nutcarb_init.f90 sourcefile~basin_module.f90 basin_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~basin_module.f90 sourcefile~carbon_module.f90 carbon_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~carbon_module.f90 sourcefile~hru_module.f90 hru_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~hru_module.f90 sourcefile~organic_mineral_mass_module.f90 organic_mineral_mass_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~organic_mineral_mass_module.f90 sourcefile~soil_data_module.f90 soil_data_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~soil_data_module.f90 sourcefile~soil_module.f90 soil_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~soil_module.f90 sourcefile~tillage_data_module.f90 tillage_data_module.f90 sourcefile~soil_nutcarb_init.f90->sourcefile~tillage_data_module.f90 sourcefile~organic_mineral_mass_module.f90->sourcefile~carbon_module.f90 sourcefile~soil_module.f90->sourcefile~carbon_module.f90

Source Code

     subroutine soil_nutcarb_init (isol)

!!    ~ ~ ~ PURPOSE ~ ~ ~
!!    this subroutine initializes soil chemical properties
!!    and intial soil layer bmix efficiency

      use hru_module, only : hru, ihru, sol_plt_ini
      use soil_module
      use soil_data_module
      use basin_module
      use organic_mineral_mass_module
      use carbon_module
      use tillage_data_module

      implicit none 
      
      integer, intent (in)  :: isol     !none       |unused soil (hru) number
      integer :: nly = 0
      integer :: ly = 0
      integer :: isolt = 0              !counter    |soil plant initialization file pointer
      integer :: isol_pl = 0            !counter    |soil nutrient initialization pointer (nutrients.sol)
      real :: wt1 = 0.                  !kg/ha      |weight of the soil layer
      real :: dep_frac = 0.             !0-1        |fraction of surface concentration at depth
      real :: frac_hum_active = 0.      !0-1        |fraction of humus in active pool - old SWAT
      real :: actp = 0.
      real :: solp = 0.
      real :: ssp = 0.
      real :: psp = 0.  
      real :: mathers_frac              !frac       !The fraction of carbon that is slow humas pool
      real :: tot_mass                  !kg/ha      |total mass of the soil layer

      !! suppress unused variable warning
      if (isol < 0) continue

      nly = soil(ihru)%nly

      !! set soil nutrient initialization from nutrients.sol
      isol_pl = hru(ihru)%dbs%soil_plant_init
      isolt = sol_plt_ini(isol_pl)%nut          ! isolt = 0 = default in type
      
      do ly = 1, nly
        if (ly == 1) then
          soil1(ihru)%cbn(ly) = max(0.001, soil(ihru)%phys(ly)%cbn)    !! assume 0.001% carbon if zero
        else
          soil1(ihru)%cbn(ly) = soil(ihru)%phys(ly)%cbn    
        endif
      enddo

      !! calculate initial nutrient contents of layers, profile and
      !! average in soil for the entire watershed

      do ly = 1, nly
        soil(ihru)%phys(ly)%conv_wt = soil(ihru)%phys(ly)%bd * soil(ihru)%phys(ly)%thick / 100.    ! mg/kg => kg/ha
        wt1 = soil(ihru)%phys(ly)%conv_wt
        
        !! set initial mineral pools - no3
        dep_frac = Exp(-solt_db(isolt)%exp_co * soil(ihru)%phys(ly)%d)
        if (solt_db(isolt)%nitrate > 1.e-9) then
          soil1(ihru)%mn(ly)%no3 = solt_db(isolt)%nitrate * dep_frac
        else
          soil1(ihru)%mn(ly)%no3 = 7. * dep_frac
        end if
        soil1(ihru)%mn(ly)%no3 =  soil1(ihru)%mn(ly)%no3 * wt1      !! mg/kg => kg/ha

        !set initial labile P pool
        if (solt_db(isolt)%lab_p > 1.e-9) then
          soil1(ihru)%mp(ly)%lab = solt_db(isolt)%lab_p * dep_frac
        else
          !! assume initial concentration of 5 mg/kg
          soil1(ihru)%mp(ly)%lab = 5. * dep_frac
        end if
        soil1(ihru)%mp(ly)%lab = soil1(ihru)%mp(ly)%lab * wt1   !! mg/kg => kg/ha

        !! set active mineral P pool based on dynamic PSP MJW
        if (bsn_cc%sol_P_model == 1) then 
          !! Allow Dynamic PSP Ratio
          !! convert to concentration
          solp = soil1(ihru)%mp(ly)%lab / wt1
          !! PSP = -0.045*log (% clay) + 0.001*(Solution P, mg kg-1) - 0.035*(% Organic C) + 0.43
          if (soil(ihru)%phys(ly)%clay > 0.) then
            psp = -0.045 * log(soil(ihru)%phys(ly)%clay) + (0.001 * solp) 
            psp = psp - (0.035 * soil1(ihru)%cbn(ly)) + 0.43 
          endif         
          !! Limit PSP range
          if (psp < .10) then
            psp = 0.10
          else if (psp > 0.7) then
            psp = 0.7
          end if
        else
          psp = hru(ihru)%nut%psp
        end if
        soil1(ihru)%mp(ly)%act = soil1(ihru)%mp(ly)%lab * (1. - psp) / psp

        !! Set Stable pool based on dynamic coefficient
        if (bsn_cc%sol_P_model == 1) then  !! From White et al 2009 
          !! convert to concentration for ssp calculation
          actp = soil1(ihru)%mp(ly)%act / wt1
          solp = soil1(ihru)%mp(ly)%lab / wt1
          !! estimate Total Mineral P in this soil based on data from sharpley 2004
          ssp = 25.044 * (actp + solp)** (-0.3833)
          !!limit SSP Range
          if (ssp > 7.) ssp = 7.
          if (ssp < 1.) ssp = 1.            
          soil1(ihru)%mp(ly)%sta = ssp * (soil1(ihru)%mp(ly)%act + soil1(ihru)%mp(ly)%lab)
        else
          !! the original code
          soil1(ihru)%mp(ly)%sta = 4. * soil1(ihru)%mp(ly)%act
        end if
      end do

      !! set initial organic pools - originally by Zhang
      do ly = 1, nly

        !!initialize total soil organic pool - no litter
        !!kg/ha = mm * t/m3 * m/1,000 mm * 1,000 kg/t * 10,000 m2/ha
        tot_mass = 10000. * soil(ihru)%phys(ly)%thick * soil(ihru)%phys(ly)%bd
        !! total mass of soil organic matter - cbn in %, and assume SOM is 58% carbon
        soil1(ihru)%tot(ly)%m = tot_mass * (soil1(ihru)%cbn(ly) / 100.) / 0.58
        soil1(ihru)%tot(ly)%c = tot_mass * (soil1(ihru)%cbn(ly) / 100.)

        ! soil1(ihru)%tot(ly)%c = soil1(ihru)%tot(ly)%m * soil1(ihru)%cbn(ly) / 100.
        soil1(ihru)%tot(ly)%n = soil1(ihru)%tot(ly)%c / 10.     !assume 10:1 C:N ratio
        soil1(ihru)%tot(ly)%p = soil1(ihru)%tot(ly)%c / 100.    !assume 100:1 C:P ratio
           
        if (bsn_cc%cswat == 0) then
          !set active humus fraction for original SWAT model
          if (solt_db(isolt)%fr_hum_act < 1.e-9) solt_db(isolt)%fr_hum_act = .02
          frac_hum_active = solt_db(isolt)%fr_hum_act
        
          !initialize original SWAT active and stable organic pools (from EPIC)
          !initialize active humus pool
          soil1(ihru)%hact(ly)%m = frac_hum_active * soil1(ihru)%tot(ly)%m
          soil1(ihru)%hact(ly)%c = frac_hum_active * soil1(ihru)%tot(ly)%c
          soil1(ihru)%hact(ly)%n = soil1(ihru)%hact(ly)%c / solt_db(isolt)%hum_c_n
          soil1(ihru)%hact(ly)%p = soil1(ihru)%hact(ly)%c / solt_db(isolt)%hum_c_p
            
          !initialize stable humus pool
          soil1(ihru)%hsta(ly)%m = (1. - frac_hum_active) * soil1(ihru)%tot(ly)%m
          soil1(ihru)%hsta(ly)%c = (1. - frac_hum_active) * soil1(ihru)%tot(ly)%c
          soil1(ihru)%hsta(ly)%n = soil1(ihru)%hsta(ly)%c / solt_db(isolt)%hum_c_n
          soil1(ihru)%hsta(ly)%p = soil1(ihru)%hsta(ly)%c / solt_db(isolt)%hum_c_p
        end if
        
        if (bsn_cc%cswat == 2 ) then
          !! initialize CENTURY organic pools - set soil humus fractions for CENTURY from DSSAT
          !! frac_not_seq is the residue/litter share derived from the user-supplied init_seq (carbon.bsn)
          org_frac%frac_not_seq = 1.0 - org_frac%frac_seq

          !!initialize passive humus pool
          soil1(ihru)%hp(ly)%m = org_frac%frac_seq * org_frac%frac_hum_passive * soil1(ihru)%tot(ly)%m
          soil1(ihru)%hp(ly)%c = org_frac%frac_seq * org_frac%frac_hum_passive * soil1(ihru)%tot(ly)%c
          soil1(ihru)%hp(ly)%n = soil1(ihru)%hp(ly)%c / 10.                   !assume 10:1 C:N ratio
          soil1(ihru)%hp(ly)%p = soil1(ihru)%hp(ly)%c / 80.                   !assume 80:1 C:P ratio
              
          !!initialize slow humus pool original method.  This is the default method
          if (org_frac% mathers_method .eqv. .false.) then
            soil1(ihru)%hs(ly)%m = org_frac%frac_seq * org_frac%frac_hum_slow * soil1(ihru)%tot(ly)%m
            soil1(ihru)%hs(ly)%c = org_frac%frac_seq * org_frac%frac_hum_slow * soil1(ihru)%tot(ly)%c
            soil1(ihru)%hs(ly)%n = soil1(ihru)%hs(ly)%c / 10.                   !assume 10:1 C:N ratio
            soil1(ihru)%hs(ly)%p = soil1(ihru)%hs(ly)%c / 80.                   !assume 80:1 C:P ratio
          endif
              
          ! initialize slow humus pool by Mathers approach ref: "Updating carbon pool initialization with DSSAT-CENTURY"
          if (org_frac%mathers_method .eqv. .true.) then
            if ((soil(ihru)%phys(ly)%clay + soil(ihru)%phys(ly)%silt) >= 0.35 ) then
              mathers_frac = (.41 + 0.0053 * 100.0 * (soil(ihru)%phys(ly)%clay + soil(ihru)%phys(ly)%silt)) / 100.0 
              soil1(ihru)%hs(ly)%c = org_frac%frac_seq * mathers_frac * soil1(ihru)%tot(ly)%c
              soil1(ihru)%hs(ly)%m = soil1(ihru)%hs(ly)%c / 0.58
              soil1(ihru)%hs(ly)%n = soil1(ihru)%hs(ly)%c / 10.                   !assume 10:1 C:N ratio
              soil1(ihru)%hs(ly)%p = soil1(ihru)%hs(ly)%c / 80.                   !assume 80:1 C:P ratio
            else
              mathers_frac = (.069 + 0.015 * 100.0 * (soil(ihru)%phys(ly)%clay + soil(ihru)%phys(ly)%silt)) / 100.0
              soil1(ihru)%hs(ly)%c = org_frac%frac_seq * mathers_frac * soil1(ihru)%tot(ly)%c  
              soil1(ihru)%hs(ly)%m = soil1(ihru)%hs(ly)%c / 0.58
              soil1(ihru)%hs(ly)%n = soil1(ihru)%hs(ly)%c / 10.                   !assume 10:1 C:N ratio
              soil1(ihru)%hs(ly)%p = soil1(ihru)%hs(ly)%c / 80.                   !assume 80:1 C:P ratio
            endif
          endif


          !!initialize microbial pool
          soil1(ihru)%microb(ly)%m = org_frac%frac_seq * org_frac%frac_hum_microb * soil1(ihru)%tot(ly)%m
          soil1(ihru)%microb(ly)%c = org_frac%frac_seq * org_frac%frac_hum_microb * soil1(ihru)%tot(ly)%c
          soil1(ihru)%microb(ly)%n = soil1(ihru)%microb(ly)%c / 8.            !assume 8:1 C:N ratio
          soil1(ihru)%microb(ly)%p = soil1(ihru)%microb(ly)%c / 80.           !assume 80:1 C:P ratio
            
          !! metabolic residue
          ! soil1(ihru)%meta(ly) = plt_mass_z
          soil1(ihru)%meta(ly)%m = org_frac%frac_not_seq * 0.85 * soil1(ihru)%tot(ly)%m
          soil1(ihru)%meta(ly)%c = org_frac%frac_not_seq * 0.85 * soil1(ihru)%tot(ly)%c              
          soil1(ihru)%meta(ly)%n = soil1(ihru)%meta(ly)%c / 10.               
          soil1(ihru)%meta(ly)%p = soil1(ihru)%meta(ly)%c / 100.   
            
          ! structural residue
          ! soil1(ihru)%str(ly) = plt_mass_z
          soil1(ihru)%str(ly)%m = org_frac%frac_not_seq * 0.15 * soil1(ihru)%tot(ly)%m
          soil1(ihru)%str(ly)%c = org_frac%frac_not_seq * 0.15 * soil1(ihru)%tot(ly)%c               
          soil1(ihru)%str(ly)%n = soil1(ihru)%str(ly)%c / 150.                !assume 150:1 C:N ratio (EPIC)
          soil1(ihru)%str(ly)%p = soil1(ihru)%str(ly)%c / 1500.
          
          !! lignin residue
          ! soil1(ihru)%lig(ly) = plt_mass_z
          soil1(ihru)%lig(ly)%m = 0.8 * soil1(ihru)%str(ly)%m
          soil1(ihru)%lig(ly)%c = 0.8 * soil1(ihru)%str(ly)%c                 !assume 80% structural c is lig
          soil1(ihru)%lig(ly)%n = 0.2 * soil1(ihru)%str(ly)%n
          soil1(ihru)%lig(ly)%p = 0.02 * soil1(ihru)%str(ly)%p
            
          !! non-lignin residue
          ! soil1(ihru)%lig(ly) = plt_mass_z
          soil1(ihru)%nonlig(ly)%m = soil1(ihru)%str(ly)%m - soil1(ihru)%lig(ly)%m 
          soil1(ihru)%nonlig(ly)%c = soil1(ihru)%str(ly)%c - soil1(ihru)%lig(ly)%c 
          soil1(ihru)%nonlig(ly)%n = soil1(ihru)%str(ly)%n - soil1(ihru)%lig(ly)%n 
          soil1(ihru)%nonlig(ly)%p = soil1(ihru)%str(ly)%p - soil1(ihru)%lig(ly)%p 
        end if

        ! set the initial biomix for each soil layer
        if (bmix_depth >= soil(ihru)%phys(ly)%d) then
          soil(ihru)%ly(ly)%init_bmix = bmix_eff
        else if (bmix_depth > soil(ihru)%phys(ly-1)%d .and. bmix_depth <= soil(ihru)%phys(ly)%d ) then
          !interpolate
          soil(ihru)%ly(ly)%init_bmix = bmix_eff * (bmix_depth - soil(ihru)%phys(ly-1)%d )/soil(ihru)%phys(ly)%thick 
        else
          soil(ihru)%ly(ly)%init_bmix = 0.0
        endif

        soil1(ihru)%tot(ly) = soil1(ihru)%str(ly) + soil1(ihru)%meta(ly) + soil1(ihru)%hs(ly) + soil1(ihru)%hp(ly) + soil1(ihru)%microb(ly)
        soil1(ihru)%seq(ly) = soil1(ihru)%hs(ly) + soil1(ihru)%hp(ly) + soil1(ihru)%microb(ly)

      end do   !! end soil layer loop 

      return
      end subroutine soil_nutcarb_init