Commit 583eae0d authored by Pillai, Himanshu's avatar Pillai, Himanshu
Browse files

new Makefile for Kokkos and updated Evaluation table

parent a75f711d
Evaluation Evaluation
================ ================
|:-------------:|:------------:|
| Tables | FORTRAN | C++ | | KERNEL |
| -------------- |:-------------:| ------------: | |----------------|:-------------:|:------------:|
| DRIVER | FORTRAN | C++ |
|----------------|:-------------:|:------------:|
| FORTRAN | Works | Works | | FORTRAN | Works | Works |
| C++ | | Works | | C++ | X | Works |
| LEGION Fortran | N/A | 1 | | LEGION Fortran | X | X |
| LEGION C++ |Partially | Works | | LEGION C++ | Partially | Works |
| KOKKOS Fortran | X | X | | KOKKOS Fortran | X | X |
| KOKKOS C++ | X | Debugging mode| | KOKKOS C++ | X | Works |
| CUDA | X | Debugging mode| | CUDA | X | Debugging |
| OpenMP/OpenACC | Works | On it | | OpenMP/OpenACC | N/A | N/A |
| FleCSI | N/A | On it | | FleCSI | N/A | Writing |
| Tpetra | N/A | N/A | | Tpetra | N/A | N/A |
| HPX | N/A | N/A | | HPX | N/A | N/A |
|----------------|:-------------:|:------------:|
SRCDIR = . SRCDIR = .
OBJECT = . OBJECT = .
.PHONY: cc_fortran_wrappers cc_serial cc_kokkos fortran .PHONY: cc_fortran_wrappers cc_serial cc_kokkos fortran cuda
default: all default: all
all: cc_fortran_wrappers fortran cc_serial cc_kokkos all: cc_fortran_wrappers fortran cc_serial cc_kokkos cuda
cc_fortran_wrappers: fortran cc_fortran_wrappers: fortran
$(MAKE) -C cc_fortran_wrappers all $(MAKE) -C cc_fortran_wrappers all
...@@ -19,15 +19,20 @@ cc_kokkos: ...@@ -19,15 +19,20 @@ cc_kokkos:
fortran: fortran:
$(MAKE) -C fortran all $(MAKE) -C fortran all
cuda:
$(MAKE) -C cuda all
clean: clean:
$(MAKE) -C cc clean $(MAKE) -C cc clean
$(MAKE) -C fortran clean $(MAKE) -C fortran clean
$(MAKE) -C cuda clean
allclean: allclean:
$(MAKE) -C cc_serial clean $(MAKE) -C cc_serial clean
$(MAKE) -C cc_kokkos clean $(MAKE) -C cc_kokkos clean
$(MAKE) -C fortran clean $(MAKE) -C fortran clean
$(MAKE) -C cuda clean
...@@ -11,7 +11,7 @@ using namespace std; ...@@ -11,7 +11,7 @@ using namespace std;
namespace ELM { namespace ELM {
template<typename Array_d> template<typename Array_d>
void CanopyHydrology_SnowWater(const double& dtime, void CanopyHydrology_SnowWater(const double& dtime,
const double& qflx_floodg, const double& qflx_floodg,
const int& ltype, const int& ltype,
const int& ctype, const int& ctype,
...@@ -41,10 +41,10 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -41,10 +41,10 @@ void CanopyHydrology_SnowWater(const double& dtime,
double& qflx_snow_h2osfc, double& qflx_snow_h2osfc,
double& frac_sno_eff, double& frac_sno_eff,
double& frac_sno) double& frac_sno)
{ {
//parameters //parameters
double rpi=4.0e0*atan(1.0e0) ; double rpi=4.0e0*atan(1.0e0) ;
double tfrz=273.15; double tfrz=273.15;
double zlnd = 0.010; double zlnd = 0.010;
...@@ -55,33 +55,38 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -55,33 +55,38 @@ void CanopyHydrology_SnowWater(const double& dtime,
// real(r8), intent(inout), dimension(-nlevsno+1:0) :: t_soisno, frac_iceold // real(r8), intent(inout), dimension(-nlevsno+1:0) :: t_soisno, frac_iceold
// real(r8), intent(inout), dimension(-nlevsno+1:0) :: dz, z, zi // real(r8), intent(inout), dimension(-nlevsno+1:0) :: dz, z, zi
//local variables //local variables
double temp_intsnow, temp_snow_depth, z_avg, fmelt, dz_snowf, snowmelt ; double temp_intsnow, temp_snow_depth, z_avg, fmelt, dz_snowf, snowmelt ;
double newsnow, bifall, accum_factor, fsno_new, smr ; double newsnow, bifall, accum_factor, fsno_new, smr ;
int j ; int j ;
//apply gridcell flood water flux to non-lake columns //apply gridcell flood water flux to non-lake columns
if (ctype != icol_sunwall && ctype != icol_shadewall) { if (ctype != icol_sunwall && ctype != icol_shadewall) {
qflx_floodc = qflx_floodg; } qflx_floodc = qflx_floodg; }
else{ else{
qflx_floodc = 0.0; qflx_floodc = 0.0;
} }
//Determine snow height and snow water //Determine snow height and snow water
//Use Alta relationship, Anderson(1976); LaChapelle(1961), //Use Alta relationship, Anderson(1976); LaChapelle(1961),
//U.S.Department of Agriculture Forest Service, Project F, //U.S.Department of Agriculture Forest Service, Project F,
//Progress Rep. 1, Alta Avalanche Study Center:Snow Layer Densification. //Progress Rep. 1, Alta Avalanche Study Center:Snow Layer Densification.
qflx_snow_h2osfc = 0.0; qflx_snow_h2osfc = 0.0;
//set temporary variables prior to updating //set temporary variables prior to updating
temp_snow_depth=snow_depth; temp_snow_depth=snow_depth;
//save initial snow content //save initial snow content
for(j = -nlevsno+1; j < snow_level; j++) { // for(j = -nlevsno+1; j < snow_level; j++) {
swe_old[j] = 0.00; // swe_old[j] = 0.00;
} // }
for(j = snow_level+1; j < 0; j++) { // for(j = snow_level+1; j < 0; j++) {
swe_old[j]=h2osoi_liq[j]+h2osoi_ice[j]; // swe_old[j]=h2osoi_liq[j]+h2osoi_ice[j];
// }
for(j = 0; j < nlevsno; j++) {
if(j < nlevsno+snow_level ) swe_old[j]=0.00;
else swe_old[j] = h2osoi_liq[j]+h2osoi_ice[j];
} }
if (do_capsnow) { if (do_capsnow) {
...@@ -157,7 +162,7 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -157,7 +162,7 @@ void CanopyHydrology_SnowWater(const double& dtime,
if(h2osno < 1.0) { if(h2osno < 1.0) {
frac_sno=std::min(frac_sno,h2osno); frac_sno=std::min(frac_sno,h2osno);
} }
} }
} else { //h2osno == 0 } else { //h2osno == 0
// initialize frac_sno and snow_depth when no snow present initially // initialize frac_sno and snow_depth when no snow present initially
...@@ -210,38 +215,38 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -210,38 +215,38 @@ void CanopyHydrology_SnowWater(const double& dtime,
} else { } else {
frac_sno_eff = 1.; frac_sno_eff = 1.;
} }
} else { } else {
frac_sno_eff = 1.; frac_sno_eff = 1.;
} }
if (ltype==istwet && t_grnd>tfrz) { if (ltype==istwet && t_grnd>tfrz) {
h2osno=0.; h2osno=0.;
snow_depth=0.; snow_depth=0.;
} }
//When the snow accumulation exceeds 10 mm, initialize snow layer //When the snow accumulation exceeds 10 mm, initialize snow layer
//Currently, the water temperature for the precipitation is simply set //Currently, the water temperature for the precipitation is simply set
//as the surface air temperature //as the surface air temperature
newnode = 0 ; //flag for when snow node will be initialized newnode = 0 ; //flag for when snow node will be initialized
if (snow_level == 0 && qflx_snow_grnd_col > 0.00 && frac_sno*snow_depth >= 0.010) { if (snow_level == 0 && qflx_snow_grnd_col > 0.00 && frac_sno*snow_depth >= 0.010) {
newnode = 1; newnode = 1;
snow_level = -1; snow_level = -1;
dz[0] = snow_depth ; //meter dz[nlevsno-1] = snow_depth ; //meter
z[0] = -0.50*dz[0]; z[nlevsno-1] = -0.50*dz[nlevsno-1];
zi[-1] = -dz[0]; zi[nlevsno-2] = -dz[nlevsno-1];
t_soisno[0] = min(tfrz, forc_air_temp) ; //K t_soisno[nlevsno-1] = min(tfrz, forc_air_temp) ; //K
h2osoi_ice[0] = h2osno ; //kg/m2 h2osoi_ice[nlevsno-1] = h2osno ; //kg/m2
h2osoi_liq[0] = 0.0 ; //kg/m2 h2osoi_liq[nlevsno-1] = 0.0 ; //kg/m2
frac_iceold[0] = 1.0; frac_iceold[nlevsno-1] = 1.0;
} }
//The change of ice partial density of surface node due to precipitation. //The change of ice partial density of surface node due to precipitation.
//Only ice part of snowfall is added here, the liquid part will be added //Only ice part of snowfall is added here, the liquid part will be added
//later. //later.
if (snow_level < 0 && newnode == 0) { if (snow_level < 0 && newnode == 0) {
h2osoi_ice[snow_level+1] = h2osoi_ice[snow_level+1]+newsnow; h2osoi_ice[nlevsno-1+snow_level+1] = h2osoi_ice[nlevsno-1+snow_level+1]+newsnow;
dz[snow_level+1] = dz[snow_level+1]+dz_snowf*dtime; dz[nlevsno-1+snow_level+1] = dz[nlevsno-1+snow_level+1]+dz_snowf*dtime;
} }
} }
} }
...@@ -18,13 +18,14 @@ FC_LIB_ROOT = /usr/lib/gcc/x86_64-linux-gnu/8 #/usr/local/Cellar/gcc/8.2.0/lib/g ...@@ -18,13 +18,14 @@ FC_LIB_ROOT = /usr/lib/gcc/x86_64-linux-gnu/8 #/usr/local/Cellar/gcc/8.2.0/lib/g
# assumes you have a working c++ compiler # assumes you have a working c++ compiler
CC = g++ CC = g++
CPP_FLAGS += -g -Wall -Wshadow -std=c++11 CPP_FLAGS += -g -Wall -Wshadow -std=c++11
CXXFLAGS += -g -Wall -Wshadow -std=c++11
STD_LIB_ROOT = /usr STD_LIB_ROOT = /usr
# assumes you have a working CUDA # assumes you have a working CUDA
NVCC=nvcc NVCC=nvcc
CUDA_FLAGS= -std=c++11 -c -arch=sm_60 CUDA_FLAGS= -std=c++11 -c -arch=sm_60
CUDA_LIBS= -lopenblas -lpthread -lcudart -lcublas CUDA_LIBS= -lopenblas -lpthread -lcudart -lcublas
CUDA_LIBDIRS=-L/usr/local/cuda-10.0/lib64 CUDA_LIBDIRS=/usr/local/cuda-10.0/lib64
CUDA_INCDIRS=-I/usr/local/cuda-10.0/include CUDA_INCDIRS=-I/usr/local/cuda-10.0/include
...@@ -36,7 +37,8 @@ CUDA_INCDIRS=-I/usr/local/cuda-10.0/include ...@@ -36,7 +37,8 @@ CUDA_INCDIRS=-I/usr/local/cuda-10.0/include
# linking flags # linking flags
CC_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(STD_LIB_ROOT)/lib -lstdc++ -L$(FC_LIB_ROOT) -lgfortran CC_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(STD_LIB_ROOT)/lib -lstdc++ -L$(FC_LIB_ROOT) -lgfortran
FC_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdff FC_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdff
CUDA_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(CUDA_LIBDIRS) -lstdc++ -lopenblas -lpthread -lcudart -lcublas CXX_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf
CUDA_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(CUDA_LIBDIRS) -lstdc++ -lpthread -lcudart -lcublas
# #
# rules # rules
...@@ -50,7 +52,7 @@ CUDA_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(CUDA_LIBDIRS) - ...@@ -50,7 +52,7 @@ CUDA_LD_FLAGS += -L$(SRCDIR) -lelm -L$(NETCDF_ROOT) -lnetcdf -L$(CUDA_LIBDIRS) -
%.cc.o: %.cc %.cc.o: %.cc
$(CC) $(CC_FLAGS) $(CPP_FLAGS) $(INC_FLAGS) -c $< -o $@ $(CC) $(CC_FLAGS) $(CPP_FLAGS) $(INC_FLAGS) -c $< -o $@
%.cpp.o: %.cpp %.cpp.o: %.cpp
$(CC) $(CC_FLAGS) $(CPP_FLAGS) $(INC_FLAGS) -c $< -o $@ $(CC) $(CC_FLAGS) $(CXXFLAGS) $(INC_FLAGS) -c $< -o $@
%.cu.o: %.cu %.cu.o: %.cu
$(NVCC) $(CUDA_FLAGS) $(INC_FLAGS) -c $< -o $@ $(NVCC) $(CUDA_FLAGS) $(INC_FLAGS) -c $< -o $@
......
...@@ -33,14 +33,14 @@ __global__ void CanopyHydrology_FracWet(const int& frac_veg_nosno, ...@@ -33,14 +33,14 @@ __global__ void CanopyHydrology_FracWet(const int& frac_veg_nosno,
template<typename Array_d> template<typename Array_d>
__global__ void CanopyHydrology_SnowWater(const double& dtime, __device__ void CanopyHydrology_SnowWater(const double& dtime,
const double& qflx_floodg, const double& qflx_floodg,
const int& ltype, const int& ltype,
const int& ctype, const int& ctype,
const bool& urbpoi, const bool& urbpoi,
const bool& do_capsnow, const bool& do_capsnow,
const int& oldfflag, const int& oldfflag,
const double& forc_t, const double& forc_air_temp,
const double& t_grnd, const double& t_grnd,
const double& qflx_snow_grnd_col, const double& qflx_snow_grnd_col,
const double& qflx_snow_melt, const double& qflx_snow_melt,
...@@ -49,15 +49,15 @@ __global__ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -49,15 +49,15 @@ __global__ void CanopyHydrology_SnowWater(const double& dtime,
double& snow_depth, double& snow_depth,
double& h2osno, double& h2osno,
double& int_snow, double& int_snow,
Array_d& swe_old, Array_d swe_old,
Array_d& h2osoi_liq, Array_d h2osoi_liq,
Array_d& h2osoi_ice, Array_d h2osoi_ice,
Array_d& t_soisno, Array_d t_soisno,
Array_d& frac_iceold, Array_d frac_iceold,
int& snl, int& snow_level,
Array_d& dz, Array_d dz,
Array_d& z, Array_d z,
Array_d& zi, Array_d zi,
int& newnode, int& newnode,
double& qflx_floodc, double& qflx_floodc,
double& qflx_snow_h2osfc, double& qflx_snow_h2osfc,
......
...@@ -4,21 +4,21 @@ ...@@ -4,21 +4,21 @@
#include "column_varcon.h" #include "column_varcon.h"
#include "clm_varpar.h" #include "clm_varpar.h"
#include "clm_varctl.h" #include "clm_varctl.h"
#include "CanopyHydrology_cpp.hh" #include "CanopyHydrology.cuh"
using namespace std; using namespace std;
namespace ELM { namespace ELM {
template<typename Array_d> template<typename Array_d>
void CanopyHydrology_SnowWater(const double& dtime, __device__ void CanopyHydrology_SnowWater(const double& dtime,
const double& qflx_floodg, const double& qflx_floodg,
const int& ltype, const int& ltype,
const int& ctype, const int& ctype,
const bool& urbpoi, const bool& urbpoi,
const bool& do_capsnow, const bool& do_capsnow,
const int& oldfflag, const int& oldfflag,
const double& forc_t, const double& forc_air_temp,
const double& t_grnd, const double& t_grnd,
const double& qflx_snow_grnd_col, const double& qflx_snow_grnd_col,
const double& qflx_snow_melt, const double& qflx_snow_melt,
...@@ -27,15 +27,15 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -27,15 +27,15 @@ void CanopyHydrology_SnowWater(const double& dtime,
double& snow_depth, double& snow_depth,
double& h2osno, double& h2osno,
double& int_snow, double& int_snow,
Array_d& swe_old, Array_d swe_old,
Array_d& h2osoi_liq, Array_d h2osoi_liq,
Array_d& h2osoi_ice, Array_d h2osoi_ice,
Array_d& t_soisno, Array_d t_soisno,
Array_d& frac_iceold, Array_d frac_iceold,
int& snl, int& snow_level,
Array_d& dz, Array_d dz,
Array_d& z, Array_d z,
Array_d& zi, Array_d zi,
int& newnode, int& newnode,
double& qflx_floodc, double& qflx_floodc,
double& qflx_snow_h2osfc, double& qflx_snow_h2osfc,
...@@ -77,11 +77,10 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -77,11 +77,10 @@ void CanopyHydrology_SnowWater(const double& dtime,
//set temporary variables prior to updating //set temporary variables prior to updating
temp_snow_depth=snow_depth; temp_snow_depth=snow_depth;
//save initial snow content //save initial snow content
for(j = -nlevsno+1; j < snl; j++) {
swe_old[j] = 0.00; for(j = 0; j < nlevsno; j++) {
} if(j < nlevsno+snow_level ) swe_old[j]=0.00;
for(j = snl+1; j < 0; j++) { else swe_old[j] = h2osoi_liq[j]+h2osoi_ice[j];
swe_old[j]=h2osoi_liq[j]+h2osoi_ice[j];
} }
if (do_capsnow) { if (do_capsnow) {
...@@ -91,10 +90,10 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -91,10 +90,10 @@ void CanopyHydrology_SnowWater(const double& dtime,
int_snow = 5.e2; int_snow = 5.e2;
} else { } else {
if (forc_t > tfrz + 2.) { if (forc_air_temp > tfrz + 2.) {
bifall=50. + 1.7*std::pow((17.0),1.5); bifall=50. + 1.7*pow((17.0),1.5);
} else if (forc_t > tfrz - 15.) { } else if (forc_air_temp > tfrz - 15.) {
bifall=50. + 1.7*std::pow((forc_t - tfrz + 15.),1.5); bifall=50. + 1.7*pow((forc_air_temp - tfrz + 15.),1.5);
} else { } else {
bifall=50.; bifall=50.;
} }
...@@ -103,7 +102,7 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -103,7 +102,7 @@ void CanopyHydrology_SnowWater(const double& dtime,
newsnow = (1. - frac_h2osfc) * qflx_snow_grnd_col * dtime; newsnow = (1. - frac_h2osfc) * qflx_snow_grnd_col * dtime;
// update int_snow // update int_snow
int_snow = std::max(int_snow,h2osno) ; //h2osno could be larger due to frost int_snow = max(int_snow,h2osno) ; //h2osno could be larger due to frost
// snowmelt from previous time step * dtime // snowmelt from previous time step * dtime
snowmelt = qflx_snow_melt * dtime; snowmelt = qflx_snow_melt * dtime;
...@@ -118,9 +117,9 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -118,9 +117,9 @@ void CanopyHydrology_SnowWater(const double& dtime,
// first compute change from melt during previous time step // first compute change from melt during previous time step
if(snowmelt > 0.) { if(snowmelt > 0.) {
smr=std::min(1.,(h2osno)/(int_snow)); smr=min(1.,(h2osno)/(int_snow));
frac_sno = 1. - std::pow((acos(min(1.,(2.*smr - 1.)))/rpi),(n_melt)) ; frac_sno = 1. - pow((acos(min(1.,(2.*smr - 1.)))/rpi),(n_melt)) ;
} }
...@@ -130,8 +129,8 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -130,8 +129,8 @@ void CanopyHydrology_SnowWater(const double& dtime,
frac_sno = fsno_new; frac_sno = fsno_new;
// reset int_snow after accumulation events // reset int_snow after accumulation events
temp_intsnow= (h2osno + newsnow) / (0.5*(cos(rpi*std::pow((1.0-std::max(frac_sno,1e-6)),(1.0/n_melt))+1.0))) ; temp_intsnow= (h2osno + newsnow) / (0.5*(cos(rpi*pow((1.0-max(frac_sno,1e-6)),(1.0/n_melt))+1.0))) ;
int_snow = std::min(1.e8,temp_intsnow) ; int_snow = min(1.e8,temp_intsnow) ;
} }
//==================================================================== //====================================================================
...@@ -152,10 +151,10 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -152,10 +151,10 @@ void CanopyHydrology_SnowWater(const double& dtime,
if (oldfflag == 1) { if (oldfflag == 1) {
// snow cover fraction in Niu et al. 2007 // snow cover fraction in Niu et al. 2007
if(snow_depth > 0.0) { if(snow_depth > 0.0) {
frac_sno = tanh(snow_depth/(2.5*zlnd*std::pow((std::min(800.0,(h2osno+ newsnow)/snow_depth)/100.0),1.0)) ) ; frac_sno = tanh(snow_depth/(2.5*zlnd*pow((min(800.0,(h2osno+ newsnow)/snow_depth)/100.0),1.0)) ) ;
} }
if(h2osno < 1.0) { if(h2osno < 1.0) {
frac_sno=std::min(frac_sno,h2osno); frac_sno=min(frac_sno,h2osno);
} }
} }
...@@ -168,8 +167,8 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -168,8 +167,8 @@ void CanopyHydrology_SnowWater(const double& dtime,
// make int_snow consistent w/ new fsno, h2osno // make int_snow consistent w/ new fsno, h2osno
int_snow = 0. ;//reset prior to adding newsnow below int_snow = 0. ;//reset prior to adding newsnow below
temp_intsnow= (h2osno + newsnow) / (0.5*(cos(rpi*std::pow((1.0-std::max(frac_sno,1e-6)),(1.0/n_melt)))+1.0)); temp_intsnow= (h2osno + newsnow) / (0.5*(cos(rpi*pow((1.0-max(frac_sno,1e-6)),(1.0/n_melt)))+1.0));
int_snow = std::min(1.e8,temp_intsnow); int_snow = min(1.e8,temp_intsnow);
// update snow_depth and h2osno to be consistent with frac_sno, z_avg // update snow_depth and h2osno to be consistent with frac_sno, z_avg
if (subgridflag ==1 && !urbpoi) { if (subgridflag ==1 && !urbpoi) {
...@@ -181,7 +180,7 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -181,7 +180,7 @@ void CanopyHydrology_SnowWater(const double& dtime,
if (oldfflag == 1) { if (oldfflag == 1) {
// snow cover fraction in Niu et al. 2007 // snow cover fraction in Niu et al. 2007
if(snow_depth > 0.0) { if(snow_depth > 0.0) {
frac_sno = tanh(snow_depth/(2.5*zlnd*std::pow((std::min(800.0,newsnow/snow_depth)/100.0),1.0)) ); frac_sno = tanh(snow_depth/(2.5*zlnd*pow((min(800.0,newsnow/snow_depth)/100.0),1.0)) );
} }
} }
} else { } else {
...@@ -223,24 +222,24 @@ void CanopyHydrology_SnowWater(const double& dtime, ...@@ -223,24 +222,24 @@ void CanopyHydrology_SnowWater(const double& dtime,
//Currently, the water temperature for the precipitation is simply set //Currently, the water temperature for the precipitation is simply set
//as the surface air temperature //as the surface air temperature
newnode = 0 ; //flag for when snow node will be initialized newnode = 0 ; //flag for when snow node will be initialized
if (snl == 0 && qflx_snow_grnd_col > 0.00 && frac_sno*snow_depth >= 0.010) { if (snow_level == 0 && qflx_snow_grnd_col > 0.00 && frac_sno*snow_depth >= 0.010) {
newnode = 1; newnode = 1;
snl = -1; snow_level = -1;
dz[0] = snow_depth ; //meter dz[nlevsno-1] = snow_depth ; //meter
z[0] = -0.50*dz[0]; z[nlevsno-1] = -0.50*dz[nlevsno-1];
zi[-1] = -dz[0]; zi[nlevsno-2] = -dz[nlevsno-1];
t_soisno[0] = min(tfrz, forc_t) ; //K t_soisno[nlevsno-1] = min(tfrz, forc_air_temp) ; //K
h2osoi_ice[0] = h2osno ; //kg/m2 h2osoi_ice[nlevsno-1] = h2osno ; //kg/m2