Updates to source while working on the point version.

Issue wwieder#12 (casa_inout.f90): Column header line added to CASACNP pool initialization/restart files.
Subroutine casa_init(), that reads the initialization/restart file, assumes there is a column header
line and was updated to skip past the first line in the file. Subroutine write_cnpepool_header()
was created and added to casa_inout.f90 to write a column header line to .csv restart files.
Subroutine write_cnpflux_header() was created and added to casa_inout.f90 to write a column header
line .csv output flux files (for example cnpflux_end_*.csv). The .csv flux output files are indeed
output only.

Issue wwieder#13 (casaoffline_driver_clm.f90): Size allocated to CORPSE output structure depends on size of met.nc file.
Read the number of years in the met.nc file to determine the number of timesteps needed in the
CORPSE output file data structures.  Previously the allocation to these output data structures
assumed a max number of 7 years. This was an OK assumption when running a 2-degree grid because
the size of met.nc files could not exceed 2GB, which was 5-7 years, depending on the number of
variables strored in the met.nc files. Seven is not a sufficient number of years for the point
weather input files that may contain many more than 7 years wihout exceeding the maximum netcdf file size.�

Issue wwieder#14 (corpse_cycle.f90, corpse_soil_carbon.f90): Updates to CORPSE fW calculation
Add parameter fWmin, the minimum soil moisture effect on decomposition, to the CORPSE parameter file.
Save CORPSE fW function (soil moisture effect on soil decomposition) to the output file.Issue wwieder#13 (casaoffline_driver_clm.f90): Size allocated to CORPSE output structure depends on size of met.nc file.

Issue wwieder#15 (mimics_cycle.f90, mimics_input.f90): Write MIMICS point output to .csv file
The subroutine WritePointMIMICS was created and added to mimics_inout.f90 to write daily output
from MIMICS point simulations to a .csv file. Some output code formerly in mimics_cycle.f90, was moved
to this subroutine.

SOURCE_CODE_05.13.2017/casa_inout.f90

@@ -830,12 +830,13 @@ SUBROUTINE casa_init(filename_cnpipool,mp,ms,mst)
   integer  np,npt,npz,nl,ns,nland,nlandz,mp,ms,mst
   real(r_2) nyearz,ivtz,istz,latz,lonz,areacellz,glaiz,slaz,isoz
 
-  print *, 'initcasa ', initcasa
+  print *, 'initcasa = ', initcasa
   !phen%phase = 2
 IF (initcasa>=1) then
       write(*,*) 
       write(*,*) "Reading initial CASACNP pool file: ", filename_cnpipool, "..."
       open(99,file=filename_cnpipool)
+      read(99,*)   ! Skip past file header. -mdh 1/30/2018
       do npt =1, mp
 !! Commented out this section (-MDH 6/14/2014)
 !!       Select Case(icycle)
@@ -972,10 +973,99 @@ SUBROUTINE casa_init(filename_cnpipool,mp,ms,mst)
         casabal%psoilocclast  = casapool%psoilocc
         casabal%sumpbal       = 0.0
 
-    Endif
+    endif
 
 end SUBROUTINE casa_init
 
+!--------------------------------------------------------------------------------
+! Write a header to the casacnp end-of-simulation pool (restart) file.
+! This subroutine assumes that the calling routine has already opened the output
+! file with unit number nout. The filename is for reference only.
+!
+
+SUBROUTINE write_cnpepool_header(icycle, nout, filename_cnpepool)
+    implicit none
+    integer, intent(in) :: icycle, nout
+    character(len=100), intent(in) :: filename_cnpepool
+
+    !                              10        20        30        40        50        60        70        80
+    !                     123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890
+    WRITE(nout,'(a72$)') 'iYrCnt,npt,veg%iveg,soil%isoilm,casamet%isorder,casamet%lat,casamet%lon,'
+    WRITE(nout,'(a82$)') 'casamet%areacell,casamet%glai,casabiome%sla(veg%iveg),phen%phase,casapool%clabile,'
+    WRITE(nout,'(a67$)') 'casapool%cplant(LEAF),casapool%cplant(WOOD),casapool%cplant(FROOT),'
+    WRITE(nout,'(a67$)') 'casapool%clitter(METB),casapool%clitter(STR),casapool%clitter(CWD),'
+    WRITE(nout,'(a62$)') 'casapool%csoil(MIC),casapool%csoil(SLOW),casapool%csoil(PASS),'
+    WRITE(nout,'(a67$)') 'casapool%nplant(LEAF),casapool%nplant(WOOD),casapool%nplant(FROOT),'
+    WRITE(nout,'(a67$)') 'casapool%nlitter(METB),casapool%nlitter(STR),casapool%nlitter(CWD),'
+    WRITE(nout,'(a80$)') 'casapool%nsoil(MIC),casapool%nsoil(SLOW),casapool%nsoil(PASS),casapool%nsoilmin,'
+    WRITE(nout,'(a67$)') 'casapool%pplant(LEAF),casapool%pplant(WOOD),casapool%pplant(FROOT),'
+    WRITE(nout,'(a67$)') 'casapool%plitter(METB),casapool%plitter(STR),casapool%plitter(CWD),'
+    WRITE(nout,'(a62$)') 'casapool%psoil(MIC),casapool%psoil(SLOW),casapool%psoil(PASS),'  
+    WRITE(nout,'(a55$)') 'casapool%psoillab,casapool%psoilsorb,casapool%psoilocc,'
+    WRITE(nout,'(a47)')  'casabal%sumcbal,casabal%sumnbal,casabal%sumpbal'
+
+end SUBROUTINE write_cnpepool_header
+
+!--------------------------------------------------------------------------------
+! This subroutine assumes that the calling routine has already opened the output 
+! file with unit number nout. The filename is for reference only.
+!
+
+SUBROUTINE write_cnpflux_header(icycle, nout, filename_cnpflux)
+    implicit none
+    integer, intent(in) :: icycle, nout
+    character(len=100), intent(in) :: filename_cnpflux
+
+
+! LEAF    = 1
+! WOOD    = 2
+! FROOT   = 3
+! 
+! METB    = 1
+! STR     = 2
+! CWD     = 3
+! 
+! MIC     = 1
+! SLOW    = 2
+! PASS    = 3
+! 
+! PLAB    = 1
+! PSORB   = 2
+! POCC    = 3
+
+    Select Case(icycle)
+    Case(1)
+    !                                  10        20        30        40        50        60        70        80
+    !                         123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890
+        WRITE(nout,'(a47$)') 'myear,npt,veg%iveg,soil%isoilm,casamet%isorder,'
+        WRITE(nout,'(a41$)') 'casamet%lat,casamet%lon,casamet%areacell,'
+        WRITE(nout,'(a70$)') 'casabal%Fcnppyear,casabal%FCrsyear,casabal%FCneeyear,casabal%FCrpyear,'
+        WRITE(nout,'(a58$)') 'clitterinput(LEAF),clitterinput(WOOD),clitterinput(FROOT),'
+        WRITE(nout,'(a50)')  'csoilinput(METB),csoilinput(METB),csoilinput(METB)'
+
+!ATTENTION: Update for N and P
+    ! Case(2)
+        ! write(nout,??) myear,npt,veg%iveg(npt),soil%isoilm(npt),casamet%isorder(npt), &
+        ! casamet%lat(npt),casamet%lon(npt),casamet%areacell(npt)*(1.0e-9), &
+        ! casabal%FCnppyear(npt),casabal%FCrsyear(npt),casabal%FCneeyear(npt),casabal%FCrpyear(npt),&
+        ! clitterinput(npt,:),csoilinput(npt,:), &
+        ! casabal%FNdepyear(npt),casabal%FNfixyear(npt),  casabal%FNsnetyear(npt), &
+        ! casabal%FNupyear(npt), casabal%FNleachyear(npt),casabal%FNlossyear(npt)
+
+    ! Case(3)
+        ! write(nout,??) myear,npt,veg%iveg(npt),soil%isoilm(npt),casamet%isorder(npt), &
+        ! casamet%lat(npt),casamet%lon(npt),casamet%areacell(npt)*(1.0e-9), &
+        ! casabal%FCnppyear(npt),casabal%FCrsyear(npt),casabal%FCneeyear(npt),casabal%FCrpyear(npt),&
+        ! clitterinput(npt,:),csoilinput(npt,:), &
+        ! casabal%FNdepyear(npt),casabal%FNfixyear(npt),  casabal%FNsnetyear(npt), &
+        ! casabal%FNupyear(npt), casabal%FNleachyear(npt),casabal%FNlossyear(npt), &
+        ! casabal%FPweayear(npt),casabal%FPdustyear(npt), casabal%FPsnetyear(npt), &
+        ! casabal%FPupyear(npt), casabal%FPleachyear(npt),casabal%FPlossyear(npt)
+
+    END Select 
+
+end SUBROUTINE write_cnpflux_header
+
 
 !--------------------------------------------------------------------------------
 SUBROUTINE casa_poolout(filename_cnpepool,iYrCnt,myear,writeToRestartCSVfile)
@@ -1055,6 +1145,8 @@ SUBROUTINE casa_poolout(filename_cnpepool,iYrCnt,myear,writeToRestartCSVfile)
      nout=103
      open(nout,file=filename_cnpepool)
 
+    ! Write a header to the casacnp end-of-simulation pool (restart) file. -mdh 1/30/2018
+    call write_cnpepool_header(icycle, nout, filename_cnpepool)
 
 !      write(*,91) nyear,cplantsum,clittersum,csoilsum 
       casabal%sumcbal=min(9999.0,max(-9999.0,casabal%sumcbal))
@@ -1209,6 +1301,9 @@ SUBROUTINE casa_fluxout(filename_cnpflux,myear,clitterinput,csoilinput,writeToRe
   if (writeToRestartCSVfile) then
       nout=104
       open(nout,file=filename_cnpflux)
+      ! Write a header to the casacnp end-of-simulation flux file. -mdh 1/30/2018
+      call write_cnpflux_header(icycle, nout, filename_cnpflux)
+
       do npt =1,mp
          Select Case(icycle)
          Case(1)
@@ -4300,7 +4395,7 @@ SUBROUTINE casacnpdriver(filename_cnpmet, filename_cnpepool, filename_cnpflux, f
              call corpse_poolfluxout(filename_corpseepool,mp,writeToRestartCSVfile)
              if (casafile%iptToSaveIndx > 0) then
                  !Write point-specific output to sPtFileNameCORPSE
-                 call WritePointCorpse(sPtFileNameCORPSE,casafile%iptToSaveIndx,mp)
+                 call WritePointCORPSE(sPtFileNameCORPSE,casafile%iptToSaveIndx,mp)
              endif
              ! Output current year's CORPSE results for transient run (-mdh 5/16/2016)
              if (mdaily == 1) then
@@ -4402,7 +4497,7 @@ SUBROUTINE WritePointFileHeaders(dirPtFile,mp)
       write(ptstr, '(i10)') iptToSave_corpse
       ptstr = adjustl(ptstr)
       sPtFileNameCORPSE = trim(dirPtFile) // 'TEST_daily_corpse_' // trim(ptstr) // '.csv' 
-      ! Subroutine WritePointCorpse is called from subroutine casacnpdriver a the end of the
+      ! Subroutine WritePointCORPSE is called from subroutine casacnpdriver a the end of the
       ! simulation to write the header and contents to sPtFileNameCORPSE.
   endif
 

SOURCE_CODE_05.13.2017/casaoffline_driver_clm.f90

@@ -33,7 +33,7 @@ PROGRAM offline_casacnp
   USE corpsevariable
 
   IMPLICIT NONE
-  integer mst, mvt, mloop, nloop, mdaily
+  integer mst, mvt, mloop, nloop, mdaily, myear
   integer tyear1, tyear2, ityear, idx, idx2
   integer mreps, irep, iYrCnt, nctime
   character(len=100) :: filename_cnppoint,filename_phen, &
@@ -245,10 +245,12 @@ PROGRAM offline_casacnp
 
   else if (isomModel == CORPSE) then
 
-      print *, 'calling corpse_init'
+      call GetMetNcFileDim(filename_cnpmet, ms, myear)
+
       !! ATTENTION: The allocation of output variables may need to be moved after met.nc file is read to get the exact # simulation years.
       if (initcasa < 2) then
-          maxSteps = mloop * 365 * 7                          ! 7 is the maximum number of years in a spinup file (daily time step - mdh 3/21/2016)
+          !maxSteps = mloop * 365 * 7                         ! 7 is the maximum number of years in a spinup file (daily time step - mdh 3/21/2016)
+          maxSteps = mloop * 365 * myear                      ! myear is the number of years in a spinup file (daily time step - mdh 1/30/2018)
       else
           maxSteps = mloop * 365 * abs(tyear2 - tyear1 + 1)   ! Assumes each transient year met.nc file has 365 days  (daily time step - mdh 3/21/2016)
       endif

SOURCE_CODE_05.13.2017/corpse_cycle.f90

@@ -256,11 +256,12 @@ SUBROUTINE corpse_soil(mp,idoy,cleaf2met,cleaf2str,croot2met,croot2str,cwd2str,c
                   theta_frzn(npt) = min(1.0, casamet%frznmoistavg(npt)/soil%ssat(npt)) ! fraction of frozen water-filled pore space (0.0 - 1.0)
                   air_filled_porosity = max(0.0, 1.0-theta_liq(npt)-theta_frzn(npt))
 
+                  ! fW and fT are output variables only.  To update actual function, go to function Resp in corpse_soil_carbon.f90. -mdh 12/18/2017
                   ! fW(npt)=(theta_liq(npt)**3+0.001)*max((air_filled_porosity)**gas_diffusion_exp,min_anaerobic_resp_factor)
-                  fW(npt) = (theta_liq(npt)**3+0.001)*max((air_filled_porosity)**2.5,0.003)
-!                 fW(npt) = max(0.05*0.022600567942709, fW(npt))  !WW added 12.14.2017 to put lower limit on CORPSE, similar to MIMICS
-                  fW(npt) = max(0.0001                , fW(npt))  !WW added 12.16.2017 to put lower limit on CORPSE, similar to MIMICS
-                  fT(npt) = 0.0 ! placeholder for future output, if needed
+                  fW(npt) = max(pt(npt)%soil(jj)%fWmin, &
+                           (theta_liq(npt)**3+0.001)*max((air_filled_porosity)**pt(npt)%soil(jj)%gas_diffusion_exp, &
+                                                          pt(npt)%soil(jj)%min_anaerobic_resp_factor))
+                  fT(npt) = 0.0 ! placeholder for future output variable, if needed
 
                   call update_pool(pool=pt(npt)%soil(jj), &
                                T=T, &

SOURCE_CODE_05.13.2017/corpse_inout.f90

@@ -12,7 +12,7 @@
 !     SUBROUTINE save_output_line - records carbon sums from one pool into the output data container 
 !                for that pool
 !     SUBROUTINE corpse_caccum - accumulate daily C fluxes and pool values and annual mean each year
-!     SUBROUTINE WritePointCorpse - write all saved output variables for a sigle point to .csv file
+!     SUBROUTINE WritePointCORPSE - write all saved output variables for a sigle point to .csv file
 !                Recordtime in corpse_params.nml file determines the number of timesteps that are saved.
 !     SUBROUTINE corpse_poolfluxout - write all saved output variables for the entire grid to restart 
 !                .csv output file for the final time saved.
@@ -552,7 +552,7 @@ END SUBROUTINE corpse_caccum
 ! Write ALL SAVED output variables and ALL SAVED timesteps for A SINGLE POINT .csv file filenamePtCORPSE.
 ! Recordtime in corpse_params.nml file determines the number of timesteps that are saved.
 !
-SUBROUTINE WritePointCorpse(filenamePtCORPSE,iptToSaveIndx,mp)
+SUBROUTINE WritePointCORPSE(filenamePtCORPSE,iptToSaveIndx,mp)
     USE casavariable
     USE define_types
     USE corpsevariable
@@ -617,7 +617,7 @@ SUBROUTINE WritePointCorpse(filenamePtCORPSE,iptToSaveIndx,mp)
 
     if (verbose .ge. 0) print *, "Done writing output to file ", trim(filenamePtCORPSE), "..."
 
-end subroutine WritePointCorpse
+end subroutine WritePointCORPSE
 
 !--------------------------------------------------------------------------------------------------------------
 ! Write END-OF-SIMULATION values for all saved output variables for the entire grid to filename_corpseepool

SOURCE_CODE_05.13.2017/corpse_soil_carbon.f90

@@ -87,6 +87,7 @@ module corpse_soil_carbon
     real::Qmax              !Pool DOC sorption capacity (See Mayes et al 2012)
     type(litterCohort),allocatable::litterCohorts(:)
     real::dissolved_carbon(nspecies)
+    real::fWmin, gas_diffusion_exp, min_anaerobic_resp_factor
 end type soil_carbon_pool
 
 !For new functionality: Separate cohorts just for rhizosphere and bulk soil
@@ -124,25 +125,25 @@ module corpse_soil_carbon
 logical                  :: microbe_driven_protection=.FALSE.  !Whether to use microbial biomass in protection rate
 real                     :: protected_carbon_decomp_factor=0.0 !vmaxref for protected carbon is multiplied by this (0.0 for total protection)
 real,dimension(nspecies) :: turnover_factor=(/1.0,1.0,1.0/)    !Factor by which each C species changes microbial turnover time (between 0 and 1)
+real                     :: fWmin=0.0                          !Minimum value of f(W) soil mositure effect on respiration (default=0.0)
 
 
 namelist /soil_carbon_nml/ &
             Ea,vmaxref,kC,Tmic,et,eup,minMicrobeC,combineThreshold,soilMaxCohorts,gas_diffusion_exp,&
             tol,enzfrac,tProtected,protection_rate,protection_species,leaching_solubility,&
             flavor_relative_solubility,DOC_deposition_rate,&
             litterDensity,CLASSIC_DECOMP,protected_relative_solubility,min_anaerobic_resp_factor,&
-            microbe_driven_protection,protected_carbon_decomp_factor,turnover_factor
+            microbe_driven_protection,protected_carbon_decomp_factor,turnover_factor,fWmin
 
 contains ! -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
 
 !-------------------------------------------------------------------------------------------------------------------------------------
 ! Allocate 2 (max(BULK,RHIZ)) litterCohorts and  
 ! set protection_rate, Qmax, and max_cohorts for the pool.
 
-subroutine init_soil_carbon(pool,protectionRate,Qmax,max_cohorts)
+subroutine init_soil_carbon(pool,Qmax,max_cohorts)
     implicit none
     type(soil_carbon_pool),intent(inout)::pool
-    real,optional,intent(in) :: protectionRate
     real,optional,intent(in) :: Qmax
     integer,optional,intent(in) :: max_cohorts
 
@@ -152,13 +153,20 @@ subroutine init_soil_carbon(pool,protectionRate,Qmax,max_cohorts)
 
     ! Default initializations for optional function arguments
     pool%protection_rate=protection_rate
+    pool%gas_diffusion_exp=gas_diffusion_exp
+    pool%min_anaerobic_resp_factor=min_anaerobic_resp_factor
+    pool%fWmin=fWmin
     pool%Qmax=0.0
     pool%max_cohorts=soilMaxCohorts
 
+!   write(*,'(1x,a36,f9.6)') 'init_soil_carbon: protection_rate = ', protection_rate
+!   write(*,'(1x,a38,f9.6)') 'init_soil_carbon: gas_diffusion_exp = ', gas_diffusion_exp
+!   write(*,'(1x,a46,f9.6)') 'init_soil_carbon: min_anaerobic_resp_factor = ', min_anaerobic_resp_factor
+!   write(*,'(1x,a26,f9.6)') 'init_soil_carbon: fWmin = ', fWmin
+
     pool%dissolved_carbon=0.0
 
     IF (present(max_cohorts)) pool%max_cohorts=min(max_cohorts,soilMaxCohorts)
-    IF (present(protectionRate)) pool%protection_rate=protectionRate
     IF (present(Qmax)) pool%Qmax=Qmax
 
     IF (allocated(pool%litterCohorts)) THEN
@@ -512,6 +520,7 @@ function Resp(Ctotal,Chet,T,theta_liq,air_filled_porosity,vmax_factor)
     real,dimension(nspecies)::Resp
     !real,dimension(nspecies)::tempresp
     real::enz,Cavail(nspecies),vmax_multiplier
+    real :: fW
 
     vmax_multiplier=1.0
     if (present(vmax_factor)) vmax_multiplier=vmax_factor
@@ -533,8 +542,12 @@ function Resp(Ctotal,Chet,T,theta_liq,air_filled_porosity,vmax_factor)
     ! Put a lower limit of 0.001 on theta_liq^3 (w.wieder, 11/7/2016).
 !   Resp=Vmax(T)*vmax_multiplier*theta_liq**3*(Cavail)*enz/(sum(Cavail)*kC+enz) &
 !        *max((air_filled_porosity)**gas_diffusion_exp,min_anaerobic_resp_factor)
-    Resp=Vmax(T)*vmax_multiplier*(theta_liq**3+0.001)*(Cavail)*enz/(sum(Cavail)*kC+enz) &
-         *max((air_filled_porosity)**gas_diffusion_exp,min_anaerobic_resp_factor)
+!   Resp=Vmax(T)*vmax_multiplier*(theta_liq**3+0.001)*(Cavail)*enz/(sum(Cavail)*kC+enz) &
+!        *max((air_filled_porosity)**gas_diffusion_exp,min_anaerobic_resp_factor)
+    !Set minimum soil water effect to namelist parameter fWmin. -mdh 12/18/2017
+    fW = max(fWmin,(theta_liq**3+0.001)*max((air_filled_porosity)**gas_diffusion_exp,min_anaerobic_resp_factor))
+    Resp=Vmax(T)*vmax_multiplier*fW*(Cavail)*enz/(sum(Cavail)*kC+enz)
+
 
     !ox_avail=oxygen_concentration(Ox,sum(tempresp)/sum(Cavail)*theta_liq*oxPerC)
     !print *,sum(tempresp)/sum(Cavail)