got GRU 1 to produce the exact same output over 1 year

build/cmake/CMakeLists.txt

@@ -25,6 +25,7 @@ if (SUNDIALS)
     set(GRU_ACTOR_DIR ${ACTORS_DIR}/gru_actor)
     set(SUMMA_DSHARE_DIR ${PARENT_DIR}/build/summa-sundials/build/source/dshare)
     set(SUMMA_ENGINE_DIR ${PARENT_DIR}/build/summa-sundials/build/source/engine)
+    set(SUMMA_NOAHMP_DIR ${PARENT_DIR}/build/summa-sundials/build/source/noah-mp)
 
     set(DIR_SUNDIALS /usr/local/sundials)
 
@@ -34,7 +35,7 @@ if (SUNDIALS)
         ${SUMMA_ENGINE_DIR}/nr_utility.f90)
     
     set(NRPROC
-        ${ENGINE_DIR}/expIntegral.f90
+        ${SUMMA_ENGINE_DIR}/expIntegral.f90
         ${SUMMA_ENGINE_DIR}/spline_int.f90)
     
     SET(HOOKUP
@@ -42,7 +43,7 @@ if (SUNDIALS)
         ${HOOKUP_DIR}/summaActors_FileManager.f90)
     
     SET(DATAMS
-        ${DSHARE_DIR}/multiconst.f90
+        ${SUMMA_DSHARE_DIR}/multiconst.f90
         ${SUMMA_DSHARE_DIR}/var_lookup.f90
         ${DSHARE_DIR}/data_types.f90
         ${DSHARE_DIR}/globalData.f90
@@ -141,10 +142,10 @@ if (SUNDIALS)
         ${ENGINE_DIR}/check_icondActors.f90)
     
     set(NOAHMP
-        ${NOAHMP_DIR}/module_model_constants.F
-        ${NOAHMP_DIR}/module_sf_noahutl.F
-        ${NOAHMP_DIR}/module_sf_noahlsm.F
-        ${NOAHMP_DIR}/module_sf_noahmplsm.F)
+        ${SUMMA_NOAHMP_DIR}/module_model_constants.F
+        ${SUMMA_NOAHMP_DIR}/module_sf_noahutl.F
+        ${SUMMA_NOAHMP_DIR}/module_sf_noahlsm.F
+        ${SUMMA_NOAHMP_DIR}/module_sf_noahmplsm.F)
     
     set(MODRUN
         ${SUMMA_ENGINE_DIR}/indexState.f90

build/source/engine/sundials/coupled_em.f90

@@ -866,6 +866,7 @@ subroutine coupled_em(&
       print*, message
       return
     end if
+    ! print*, 'prog_data%var(iLookPROG%mLayerVolFracLiq)%dat()', prog_data%var(iLookPROG%mLayerVolFracLiq)%dat(1:nSoil)
   
 
     ! process the flag for too much melt
@@ -1037,6 +1038,7 @@ subroutine coupled_em(&
 
   end do  substeps ! (sub-step loop)
 
+
   ! *** add snowfall to the snowpack...
   ! -----------------------------------
 

build/unneeded/systemSolvSundials.f90

@@ -343,7 +343,7 @@ subroutine systemSolvSundials(&
     ! try to accelerate solution for energy
     if(ixCasNrg/=integerMissing) stateVecTrial(ixCasNrg) = stateVecInit(ixCasNrg) + (airtemp - stateVecInit(ixCasNrg))*tempAccelerate
     if(ixVegNrg/=integerMissing) stateVecTrial(ixVegNrg) = stateVecInit(ixVegNrg) + (airtemp - stateVecInit(ixVegNrg))*tempAccelerate
-
+    
     ! compute H_T at the beginning of the data step
     call t2enthalpy_T(&
                     ! input: data structures

build/unneeded/updateVarsSundials.f90

@@ -546,6 +546,7 @@ subroutine updateVarsSundials(&
 
         ! case of energy state OR coupled solution (or adjusting the temperature)
         elseif(do_adjustTemp .or. ( (isNrgState .or. isCoupled) ) )then
+        print*, 'Before mLayerVolFracLiqTrial', mLayerVolFracLiqTrial
 
           ! identify domain type
           select case(ixDomainType)
@@ -638,6 +639,7 @@ subroutine updateVarsSundials(&
           endif
 
         endif  ! if energy state or solution is coupled
+        print*, 'After mLayerVolFracLiqTrial', mLayerVolFracLiqTrial
 
         ! -----
         ! ------------------------