Commit 1b8a59c6 authored by Pillai, Himanshu's avatar Pillai, Himanshu
Browse files

Cleaning code

parent 89f99d7f
......@@ -15,12 +15,11 @@ EXEC_TESTS = CanopyHydrology_kern1_single \
CanopyHydrology_kern1_multiple \
CanopyHydrology_module
.PHONY: links library test
default: all
all: links library $(TESTS)
all: logo links library $(TESTS)
test: $(EXEC_TESTS)
python ../compare_to_gold.py $(TESTS)
......@@ -79,3 +78,5 @@ links:
@echo "making in links"
$(MAKE) -C ../links links
include $(ELM_BASE)/config/Makefile.logo
......@@ -3,6 +3,10 @@
#ifndef ELM_KERNEL_TEST_UTILS_HH_
#define ELM_KERNEL_TEST_UTILS_HH_
#include <iostream>
#include <memory>
#include <type_traits>
namespace ELM {
namespace Utils {
......@@ -10,56 +14,83 @@ namespace Utils {
template<size_t N, typename T=double>
class VectorStatic {
public:
VectorStatic() {}
VectorStatic(T t) { *this = t; }
VectorStatic() :
do_(std::shared_ptr<std::array<T,N>>(new std::array<T,N>())) {
d_ = do_.get();
}
VectorStatic(T t) :
VectorStatic()
{ *this = t; }
VectorStatic(const std::shared_ptr<std::array<T,N> >& d) :
do_(d) {
d_ = do_.get();
}
VectorStatic(std::array<T,N>* d) :
do_(nullptr),
d_(d)
{}
VectorStatic(const VectorStatic<N,T>& other) = default;
T& operator()(size_t i) { assert(0 <= i && i < N); return d_[i]; }
const T& operator()(size_t i) const { assert(0 <= i && i < N); return d_[i]; }
T& operator()(size_t i) { assert(0 <= i && i < N); return (*d_)[i]; }
const T& operator()(size_t i) const { assert(0 <= i && i < N); return (*d_)[i]; }
T& operator[](size_t i) { assert(0 <= i && i < N); return d_[i]; }
const T& operator[](size_t i) const { assert(0 <= i && i < N); return d_[i]; }
T& operator[](size_t i) { assert(0 <= i && i < N); return (*d_)[i]; }
const T& operator[](size_t i) const { assert(0 <= i && i < N); return (*d_)[i]; }
void operator=(T t) {
for (size_t i=0; i!=N; ++i) {
d_[i] = t;
(*d_)[i] = t;
}
}
double const * begin() const { return d_.begin(); }
double const * end() const { return d_.end(); }
double const * begin() const { return &((*d_)[0]); }
double const * end() const { return &((*d_)[N]); }
private:
std::array<T,N> d_;
std::array<T,N>* d_;
std::shared_ptr< std::array<T,N> > do_;
};
template<size_t ROW, size_t COL, typename T=double>
class MatrixStatic {
using Data_t = std::array<std::array<T,COL>, ROW>;
public:
MatrixStatic() {}
MatrixStatic(T t) { *this = t; }
MatrixStatic() :
d_(std::shared_ptr<Data_t>(new Data_t())) {}
MatrixStatic(T t) :
MatrixStatic()
{ *this = t; }
MatrixStatic(const MatrixStatic<ROW,COL,T>& other) = default;
T& operator()(size_t i, size_t j) { assert(0 <= i && i < ROW && 0 <= j && j < COL); return d_[i][j]; }
const T& operator()(size_t i, size_t j) const { assert(0 <= i && i < ROW && 0 <= j && j < COL); return d_[i][j]; }
T& operator()(size_t i, size_t j) { assert(0 <= i && i < ROW && 0 <= j && j < COL); return (*d_)[i][j]; }
const T& operator()(size_t i, size_t j) const { assert(0 <= i && i < ROW && 0 <= j && j < COL); return (*d_)[i][j]; }
VectorStatic<COL,T>& operator[](size_t i) { assert(0 <= i && i < ROW); return d_[i]; }
const VectorStatic<COL,T>& operator[](size_t i) const { assert(0 <= i && i < ROW); return d_[i]; }
VectorStatic<COL,T> operator[](size_t i) { assert(0 <= i && i < ROW); return &((*d_)[i]); }
const VectorStatic<COL,T> operator[](size_t i) const { assert(0 <= i && i < ROW); return &((*d_)[i]); }
void operator=(T t) {
for (size_t i=0; i!=ROW; ++i) {
for (size_t j=0; j!=COL; ++j) {
d_[i][j] = t;
(*d_)[i][j] = t;
}
}
}
double const * begin() const { return &d_[0][0]; }
double const * end() const { return &d_[ROW-1][COL-1] +1; }
double const * begin() const { return &(*d_)[0][0]; }
double const * end() const { return &(*d_)[ROW-1][COL-1] +1; }
private:
std::array<VectorStatic<COL,T>,ROW> d_;
std::shared_ptr<Data_t> d_;
};
......
......@@ -30,9 +30,6 @@ static const int n_max_times = 31 * 24 * 2; // max days per month times hours pe
// day * half hour timestep
static const int n_grid_cells = 24;
// using MatrixState = MatrixStatic<n_grid_cells, n_pfts>;
// using MatrixForc = MatrixStatic<n_max_times,n_grid_cells>;
} // namespace
} // namespace
......@@ -81,15 +78,10 @@ int main(int argc, char ** argv)
ViewMatrixType::HostMirror h_esai = Kokkos::create_mirror_view( esai );
// phenology state
// ELM::Utils::MatrixState elai;
// ELM::Utils::MatrixState esai;
ELM::Utils::read_phenology("../links/surfacedataWBW.nc", n_months, n_pfts, 0, h_elai, h_esai);
ELM::Utils::read_phenology("../links/surfacedataBRW.nc", n_months, n_pfts, n_months, h_elai, h_esai);
// forcing state
// ELM::Utils::MatrixForc forc_rain;
// ELM::Utils::MatrixForc forc_snow;
// ELM::Utils::MatrixForc forc_air_temp;
ViewMatrixType forc_rain( "forc_rain", n_max_times,n_grid_cells );
ViewMatrixType forc_snow( "forc_snow", n_max_times,n_grid_cells );
......@@ -100,23 +92,8 @@ int main(int argc, char ** argv)
const int n_times = ELM::Utils::read_forcing("../links/forcing", n_max_times, 0, n_grid_cells, h_forc_rain, h_forc_snow, h_forc_air_temp);
ViewMatrixType forc_irrig( "forc_irrig", n_max_times,n_grid_cells );
ViewMatrixType::HostMirror h_forc_irrig = Kokkos::create_mirror_view( forc_irrig );
//ELM::Utils::MatrixForc forc_irrig; forc_irrig = 0.;
// output state by the grid cell
// auto qflx_prec_intr = std::array<double,n_grid_cells>();
// auto qflx_irrig = std::array<double,n_grid_cells>();
// auto qflx_prec_grnd = std::array<double,n_grid_cells>();
// auto qflx_snwcp_liq = std::array<double,n_grid_cells>();
// auto qflx_snwcp_ice = std::array<double,n_grid_cells>();
// auto qflx_snow_grnd_patch = std::array<double,n_grid_cells>();
// auto qflx_rain_grnd = std::array<double,n_grid_cells>();
// auto qflx_prec_intr = ELM::Utils::MatrixState();
// auto qflx_irrig = ELM::Utils::MatrixState();
// auto qflx_prec_grnd = ELM::Utils::MatrixState();
// auto qflx_snwcp_liq = ELM::Utils::MatrixState();
// auto qflx_snwcp_ice = ELM::Utils::MatrixState();
// auto qflx_snow_grnd_patch = ELM::Utils::MatrixState();
// auto qflx_rain_grnd = ELM::Utils::MatrixState();
ViewMatrixType qflx_prec_intr( "qflx_prec_intr", n_grid_cells, n_pfts );
ViewMatrixType qflx_irrig( "qflx_irrig", n_grid_cells, n_pfts );
ViewMatrixType qflx_prec_grnd( "qflx_prec_grnd", n_grid_cells, n_pfts );
......@@ -134,44 +111,28 @@ int main(int argc, char ** argv)
// output state by the pft
ViewMatrixType h2o_can( "h2o_can", n_grid_cells, n_pfts );
//auto h2o_can = ELM::Utils::MatrixState();
ViewMatrixType::HostMirror h_h2o_can = Kokkos::create_mirror_view( h2o_can );
//h_h2o_can = 0.;
//auto h2o_can1 = ELM::Utils::MatrixState();
// Array<int64_t, 2> a = h_h2o_can;
//const size_t n0 = h2o_can.extent_0 ();
// const int64_t begin0 = h_h2o_can.begin();
// const int64_t end0= h_h2o_can.end();
Kokkos::deep_copy( elai, h_elai);
Kokkos::deep_copy( esai, h_esai);
Kokkos::deep_copy( forc_rain, h_forc_rain);
Kokkos::deep_copy( forc_snow, h_forc_snow);
//Kokkos::deep_copy( forc_irrig, h_forc_irrig);
Kokkos::deep_copy( forc_air_temp, h_forc_air_temp);
//Kokkos::deep_copy( qflx_prec_intr, h_qflx_prec_intr);
// Kokkos::deep_copy( qflx_irrig, h_qflx_irrig);
// Kokkos::deep_copy( qflx_prec_grnd, h_qflx_prec_grnd);
// Kokkos::deep_copy( qflx_snwcp_liq, h_qflx_snwcp_liq);
// Kokkos::deep_copy( qflx_snwcp_ice, h_qflx_snwcp_ice);
// Kokkos::deep_copy( qflx_snow_grnd_patch, h_qflx_snow_grnd_patch);
// Kokkos::deep_copy( qflx_rain_grnd, h_qflx_rain_grnd);
// Kokkos::deep_copy( h2o_can, h_h2o_can);
double* end = &h_h2o_can(n_grid_cells-1, n_pfts-1) ;
std::ofstream soln_file;
soln_file.open("test_CanopyHydrology_kern1_multiple.soln");
soln_file << "Time\t Total Canopy Water\t Min Water\t Max Water" << std::endl;
std::cout << "Time\t Total Canopy Water\t Min Water\t Max Water" << std::endl;
auto min_max = std::minmax_element(&h_h2o_can(0,0), end+1);//h2o_can1.begin(), h2o_can1.end());
auto min_max = std::minmax_element(&h_h2o_can(0,0), end+1);
soln_file << std::setprecision(16)
<< 0 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.) //h2o_can1.begin(), h2o_can1.end(), 0.)
<< 0 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)
<< "\t" << *min_max.first
<< "\t" << *min_max.second << std::endl;
std::cout << std::setprecision(16)
<< 0 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.) //h2o_can1.begin(), h2o_can1.end(), 0.)
<< 0 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)
<< "\t" << *min_max.first
<< "\t" << *min_max.second << std::endl;
......@@ -189,16 +150,7 @@ int main(int argc, char ** argv)
// Construct 2D MDRangePolicy: lower and upper bounds provided, tile dims defaulted
MDPolicyType_2D mdpolicy_2d( {{0,0}}, {{n_grid_cells,n_pfts}} );
// // Execute parallel_for with rank 2 MDRangePolicy
// Kokkos::parallel_for( "md2d", mdpolicy_2d, ELM::CanopyHydrology_Interception(dtime,
// forc_rain(t,g), forc_snow(t,g), forc_irrig(t,g),
// ltype, ctype, urbpoi, do_capsnow,
// elai(g,p), esai(g,p), dewmx, frac_veg_nosno,
// h2o_can(g,p), n_irrig_steps_left,
// qflx_prec_intr(g,p), qflx_irrig(g,p), qflx_prec_grnd(g,p),
// qflx_snwcp_liq(g,p), qflx_snwcp_ice(g,p),
// qflx_snow_grnd_patch(g,p), qflx_rain_grnd(g,p)); );
Kokkos::parallel_for("md2d",mdpolicy_2d,KOKKOS_LAMBDA (const size_t& g, const size_t& p) {
Kokkos::parallel_for("md2d",mdpolicy_2d,KOKKOS_LAMBDA (const size_t& g, const size_t& p) {
ELM::CanopyHydrology_Interception(dtime,
forc_rain(t,g), forc_snow(t,g), forc_irrig(t,g),
ltype, ctype, urbpoi, do_capsnow,
......@@ -208,21 +160,6 @@ int main(int argc, char ** argv)
qflx_snwcp_liq(g,p), qflx_snwcp_ice(g,p),
qflx_snow_grnd_patch(g,p), qflx_rain_grnd(g,p)); });
// parallel_for (TeamPolicy<> (0,n_grid_cells), KOKKOS_LAMBDA (const size_t& g)
// {
// parallel_for (TeamThreadRange (g, n_pfts),
// [=] (size_t p) {
// ELM::CanopyHydrology_Interception(dtime,
// forc_rain(t,g), forc_snow(t,g), forc_irrig(t,g),
// ltype, ctype, urbpoi, do_capsnow,
// elai(g,p), esai(g,p), dewmx, frac_veg_nosno,
// h2o_can(g,p), n_irrig_steps_left,
// qflx_prec_intr(g,p), qflx_irrig(g,p), qflx_prec_grnd(g,p),
// qflx_snwcp_liq(g,p), qflx_snwcp_ice(g,p),
// qflx_snow_grnd_patch(g,p), qflx_rain_grnd(g,p));
// });
// });
Kokkos::deep_copy( h_qflx_irrig, qflx_irrig);
Kokkos::deep_copy( h_qflx_prec_intr, qflx_prec_intr);
Kokkos::deep_copy( h_qflx_prec_grnd, qflx_prec_grnd);
......@@ -232,13 +169,13 @@ int main(int argc, char ** argv)
Kokkos::deep_copy( h_qflx_rain_grnd, qflx_rain_grnd);
Kokkos::deep_copy( h_h2o_can, h2o_can);
auto min_max = std::minmax_element(&h_h2o_can(0,0), end+1);//h2o_can1.begin(), h2o_can1.end());
auto min_max = std::minmax_element(&h_h2o_can(0,0), end+1);
std::cout << std::setprecision(16)
<< t+1 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)//h2o_can1.begin(), h2o_can1.end(), 0.)
<< t+1 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)
<< "\t" << *min_max.first
<< "\t" << *min_max.second << std::endl;
soln_file << std::setprecision(16)
<< t+1 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)//h2o_can1.begin(), h2o_can1.end(), 0.)
<< t+1 << "\t" << std::accumulate(&h_h2o_can(0,0), end+1, 0.)
<< "\t" << *min_max.first
<< "\t" << *min_max.second << std::endl;
......
......@@ -67,8 +67,7 @@ int main(int argc, char ** argv)
// phenology state
typedef Kokkos::View<double**> ViewMatrixType;
// ELM::Utils::MatrixState elai;
// ELM::Utils::MatrixState esai;
ViewMatrixType elai( "elai", n_months, n_pfts );
ViewMatrixType esai( "esai", n_months, n_pfts );
ViewMatrixType::HostMirror h_elai = Kokkos::create_mirror_view( elai );
......@@ -77,9 +76,7 @@ int main(int argc, char ** argv)
// forcing state
// ELM::Utils::MatrixForc forc_rain;
// ELM::Utils::MatrixForc forc_snow;
// ELM::Utils::MatrixForc forc_air_temp;
ViewMatrixType forc_rain( "forc_rain", n_max_times,1 );
ViewMatrixType forc_snow( "forc_snow", n_max_times,1 );
ViewMatrixType forc_air_temp( "forc_air_temp", n_max_times,1 );
......
......@@ -29,12 +29,6 @@ static const int n_max_times = 31 * 24 * 2; // max days per month times hours pe
static const int n_grid_cells = 24;
static const int n_levels_snow = 5;
// using MatrixStatePFT = MatrixStatic<n_grid_cells, n_pfts>;
// using MatrixStateSoilColumn = MatrixStatic<n_grid_cells, n_levels_snow>;
// using MatrixForc = MatrixStatic<n_max_times,n_grid_cells>;
// using VectorColumn = VectorStatic<n_grid_cells>;
// using VectorColumnInt = VectorStatic<n_grid_cells,int>;
} // namespace
} // namespace
......@@ -84,8 +78,7 @@ int main(int argc, char ** argv)
typedef Kokkos::View<double**> ViewMatrixType;
typedef Kokkos::View<int**> ViewMatrixType1;
typedef Kokkos::View<int*> ViewVectorType1;
// ELM::Utils::MatrixState elai;
// ELM::Utils::MatrixState esai;
ViewMatrixType elai( "elai", n_grid_cells, n_pfts );
ViewMatrixType esai( "esai", n_grid_cells, n_pfts );
ViewMatrixType::HostMirror h_elai = Kokkos::create_mirror_view( elai );
......@@ -101,7 +94,7 @@ int main(int argc, char ** argv)
ViewMatrixType::HostMirror h_forc_snow = Kokkos::create_mirror_view( forc_snow );
ViewMatrixType::HostMirror h_forc_air_temp = Kokkos::create_mirror_view( forc_air_temp );
const int n_times = ELM::Utils::read_forcing("../links/forcing", n_max_times, 0, n_grid_cells, h_forc_rain, h_forc_snow, h_forc_air_temp);
//ELM::Utils::MatrixForc forc_irrig; forc_irrig = 0.;
ViewMatrixType forc_irrig( "forc_irrig", n_max_times,n_grid_cells );
ViewMatrixType::HostMirror h_forc_irrig = Kokkos::create_mirror_view( forc_irrig );
double qflx_floodg = 0.0;
......@@ -111,9 +104,7 @@ int main(int argc, char ** argv)
//
// NOTE: in a real case, these would be populated, but we don't actually
// // need them to be for these kernels. --etc
// auto z = ELM::Utils::MatrixStateSoilColumn(0.);
// auto zi = ELM::Utils::MatrixStateSoilColumn(0.);
// auto dz = ELM::Utils::MatrixStateSoilColumn(0.);
ViewMatrixType z( "z", n_grid_cells, n_levels_snow );
ViewMatrixType zi( "zi", n_grid_cells, n_levels_snow );
ViewMatrixType dz( "dz", n_grid_cells, n_levels_snow );
......@@ -122,12 +113,7 @@ int main(int argc, char ** argv)
ViewMatrixType::HostMirror h_dz = Kokkos::create_mirror_view( dz );
// state variables that require ICs and evolve (in/out)
// auto h2ocan = ELM::Utils::MatrixStatePFT(0.);
// auto swe_old = ELM::Utils::MatrixStateSoilColumn(0.);
// auto h2osoi_liq = ELM::Utils::MatrixStateSoilColumn(0.);
// auto h2osoi_ice = ELM::Utils::MatrixStateSoilColumn(0.);
// auto t_soisno = ELM::Utils::MatrixStateSoilColumn(0.);
// auto frac_iceold = ELM::Utils::MatrixStateSoilColumn(0.);
ViewMatrixType h2ocan( "h2ocan", n_grid_cells, n_pfts );
ViewMatrixType swe_old( "swe_old", n_grid_cells, n_levels_snow );
ViewMatrixType h2osoi_liq( "h2osoi_liq", n_grid_cells, n_levels_snow );
......@@ -141,10 +127,7 @@ int main(int argc, char ** argv)
ViewMatrixType::HostMirror h_t_soisno = Kokkos::create_mirror_view( t_soisno );
ViewMatrixType::HostMirror h_frac_iceold = Kokkos::create_mirror_view( frac_iceold );
// auto t_grnd = ELM::Utils::VectorColumn(0.);
// auto h2osno = ELM::Utils::VectorColumn(0.);
// auto snow_depth = ELM::Utils::VectorColumn(0.);
// auto snl = ELM::Utils::VectorColumnInt(0.); // note this tracks the snow_depth
ViewVectorType t_grnd( "t_grnd", n_grid_cells );
ViewVectorType h2osno( "h2osno", n_grid_cells );
ViewVectorType snow_depth( "snow_depth", n_grid_cells );
......@@ -163,13 +146,7 @@ int main(int argc, char ** argv)
// output fluxes by pft
// auto qflx_prec_intr = ELM::Utils::MatrixStatePFT();
// auto qflx_irrig = ELM::Utils::MatrixStatePFT();
// auto qflx_prec_grnd = ELM::Utils::MatrixStatePFT();
// auto qflx_snwcp_liq = ELM::Utils::MatrixStatePFT();
// auto qflx_snwcp_ice = ELM::Utils::MatrixStatePFT();
// auto qflx_snow_grnd_patch = ELM::Utils::MatrixStatePFT();
// auto qflx_rain_grnd = ELM::Utils::MatrixStatePFT();
ViewMatrixType qflx_prec_intr( "qflx_prec_intr", n_grid_cells, n_pfts );
ViewMatrixType qflx_irrig( "qflx_irrig", n_grid_cells, n_pfts );
ViewMatrixType qflx_prec_grnd( "qflx_prec_grnd", n_grid_cells, n_pfts );
......@@ -187,15 +164,12 @@ int main(int argc, char ** argv)
// FIXME: I have no clue what this is... it is inout on WaterSnow. For now I
// am guessing the data structure. Ask Scott. --etc
//auto integrated_snow = ELM::Utils::VectorColumn(0.);
ViewVectorType integrated_snow( "integrated_snow", n_grid_cells );
ViewVectorType::HostMirror h_integrated_snow = Kokkos::create_mirror_view( integrated_snow);
// output fluxes, state by the column
// auto qflx_snow_grnd_col = ELM::Utils::VectorColumn();
// auto qflx_snow_h2osfc = ELM::Utils::VectorColumn();
// auto qflx_h2osfc2topsoi = ELM::Utils::VectorColumn();
// auto qflx_floodc = ELM::Utils::VectorColumn();
ViewVectorType qflx_snow_grnd_col( "qflx_snow_grnd_col", n_grid_cells );
ViewVectorType qflx_snow_h2osfc( "qflx_snow_h2osfc", n_grid_cells );
ViewVectorType qflx_h2osfc2topsoi( "qflx_h2osfc2topsoi", n_grid_cells );
......@@ -205,21 +179,13 @@ int main(int argc, char ** argv)
ViewVectorType::HostMirror h_qflx_h2osfc2topsoi = Kokkos::create_mirror_view( qflx_h2osfc2topsoi);
ViewVectorType::HostMirror h_qflx_floodc = Kokkos::create_mirror_view( qflx_floodc);
// auto frac_sno_eff = ELM::Utils::VectorColumn();
// auto frac_sno = ELM::Utils::VectorColumn();
ViewVectorType frac_sno_eff( "frac_sno_eff", n_grid_cells );
ViewVectorType frac_sno( "frac_sno", n_grid_cells );
ViewVectorType::HostMirror h_frac_sno_eff = Kokkos::create_mirror_view( frac_sno_eff);
ViewVectorType::HostMirror h_frac_sno = Kokkos::create_mirror_view( frac_sno);
// std::cout << "Time\t Total Canopy Water\t Min Water\t Max Water" << std::endl;
// auto min_max = std::minmax_element(h2ocan.begin(), h2ocan.end());
// std::cout << std::setprecision(16)
// << 0 << "\t" << std::accumulate(h2ocan.begin(), h2ocan.end(), 0.)
// << "\t" << *min_max.first
// << "\t" << *min_max.second << std::endl;
Kokkos::deep_copy( elai, h_elai);
Kokkos::deep_copy( esai, h_esai);
Kokkos::deep_copy( forc_rain, h_forc_rain);
......@@ -264,21 +230,7 @@ int main(int argc, char ** argv)
// Construct 2D MDRangePolicy: lower and upper bounds provided, tile dims defaulted
MDPolicyType_2D mdpolicy_2d( {{0,0}}, {{n_grid_cells,n_pfts}} );
// Kokkos::parallel_for("md2d",mdpolicy_2d,KOKKOS_LAMBDA (const size_t& g, const size_t& p) {
// ELM::CanopyHydrology_Interception(dtime,
// forc_rain(t,g), forc_snow(t,g), forc_irrig(t,g),
// ltype, ctype, urbpoi, do_capsnow,
// elai(g,p), esai(g,p), dewmx, frac_veg_nosno,
// h2ocan(g,p), n_irrig_steps_left,
// qflx_prec_intr(g,p), qflx_irrig(g,p), qflx_prec_grnd(g,p),
// qflx_snwcp_liq(g,p), qflx_snwcp_ice(g,p),
// qflx_snow_grnd_patch(g,p), qflx_rain_grnd(g,p));
// double fwet = 0., fdry = 0.;
// ELM::CanopyHydrology_FracWet(frac_veg_nosno, h2ocan(g,p), elai(g,p), esai(g,p), dewmx, fwet, fdry);
// });
// Column level operations
// NOTE: this is effectively an accumulation kernel/task! --etc
typedef Kokkos::TeamPolicy<> team_policy ;
......@@ -325,47 +277,12 @@ int main(int argc, char ** argv)
h2osno(g), h2osfc(g), h2osoi_liq(g,0), frac_sno(g), frac_sno_eff(g),
qflx_h2osfc2topsoi(g), frac_h2osfc(g));
});
// Kokkos::parallel_reduce( Kokkos::RangePolicy<execution_space>(0,n_grid_cells), KOKKOS_LAMBDA (const size_t& g, double& upd) {
// upd += qflx_snow_grnd_patch(g,p);
// }, sum);
// qflx_snow_grnd_col(g) = sum ;
// Calculate ?water balance? on the snow column, adding throughfall,
// removing melt, etc.
//
// local outputs
// Kokkos::parallel_for (n_grid_cells,
// KOKKOS_LAMBDA (const size_t& g) {
// int newnode;
// ELM::CanopyHydrology_SnowWater(dtime, qflx_floodg,
// ltype, ctype, urbpoi, do_capsnow, oldfflag,
// forc_air_temp(t,g), t_grnd(g),
// qflx_snow_grnd_col(g), qflx_snow_melt, n_melt, frac_h2osfc(g),
// snow_depth(g), h2osno(g), integrated_snow(g), Kokkos::subview(swe_old, g , Kokkos::ALL),
// Kokkos::subview(h2osoi_liq, g , Kokkos::ALL), Kokkos::subview(h2osoi_ice, g , Kokkos::ALL), Kokkos::subview(t_soisno, g , Kokkos::ALL), Kokkos::subview(frac_iceold, g , Kokkos::ALL),
// snow_level(g), Kokkos::subview(dz, g , Kokkos::ALL), Kokkos::subview(z, g , Kokkos::ALL), Kokkos::subview(zi, g , Kokkos::ALL), newnode,
// qflx_floodc(g), qflx_snow_h2osfc(g), frac_sno_eff(g), frac_sno(g));
// // Calculate Fraction of Water to the Surface?
// //
// // FIXME: Fortran black magic... h2osoi_liq is a vector, but the
// // interface specifies a single double. For now passing the 0th
// // entry. --etc
// ELM::CanopyHydrology_FracH2OSfc(dtime, min_h2osfc, ltype, micro_sigma,
// h2osno(g), h2osfc(g), h2osoi_liq(g,0), frac_sno(g), frac_sno_eff(g),
// qflx_h2osfc2topsoi(g), frac_h2osfc(g));
// }); // end grid cell loop
// auto min_max = std::minmax_element(h2ocan.begin(), h2ocan.end());
// std::cout << std::setprecision(16)
// << t+1 << "\t" << std::accumulate(h2ocan.begin(), h2ocan.end(), 0.)
// << "\t" << *min_max.first
// << "\t" << *min_max.second << std::endl;
Kokkos::deep_copy( h_qflx_irrig, qflx_irrig);
Kokkos::deep_copy( h_qflx_prec_intr, qflx_prec_intr);
Kokkos::deep_copy( h_qflx_prec_grnd, qflx_prec_grnd);
......
/* ---------------------------------------------
Makefile constructed configuration:
Thu Aug 1 11:41:20 EDT 2019
Thu Aug 15 11:40:05 EDT 2019
----------------------------------------------*/
#if !defined(KOKKOS_MACROS_HPP) || defined(KOKKOS_CORE_CONFIG_H)
#error "Do not include KokkosCore_config.h directly; include Kokkos_Macros.hpp instead."
......
/* ---------------------------------------------
Makefile constructed configuration:
Thu Aug 1 12:19:53 EDT 2019
Thu Aug 15 13:01:42 EDT 2019
----------------------------------------------*/
#if !defined(KOKKOS_MACROS_HPP) || defined(KOKKOS_CORE_CONFIG_H)
#error "Do not include KokkosCore_config.h directly; include Kokkos_Macros.hpp instead."
......
......@@ -3,6 +3,10 @@
#ifndef ELM_KERNEL_TEST_UTILS_HH_
#define ELM_KERNEL_TEST_UTILS_HH_
#include <iostream>
#include <memory>
#include <type_traits>
namespace ELM {
namespace Utils {
......@@ -10,53 +14,83 @@ namespace Utils {
template<size_t N, typename T=double>
class VectorStatic {
public:
VectorStatic() {}
VectorStatic(T t) { *this = t; }
VectorStatic() :
do_(std::shared_ptr<std::array<T,N>>(new std::array<T,N>())) {
d_ = do_.get();
}
VectorStatic(T t) :
VectorStatic()
{ *this = t; }
VectorStatic(const std::shared_ptr<std::array<T,N> >& d) :
do_(d) {
d_ = do_.get();
}
VectorStatic(std::array<T,N>* d) :
do_(nullptr),
d_(d)
{}
VectorStatic(const VectorStatic<N,T>& other) = default;
T& operator()(size_t i) { return d_[i]; }
const T& operator()(size_t i) const { return d_[i]; }
T& operator()(size_t i) { assert(0 <= i && i < N); return (*d_)[i]; }
const T& operator()(size_t i) const { assert(0 <= i && i < N); return (*d_)[i]; }
T& operator[](size_t i) { return d_[i]; }
const T& operator[](size_t i) const { return d_[i]; }
T& operator[](size_t i) { assert(0 <= i && i < N); return (*d_)[i]; }
const T& operator[](size_t i) const { assert(0 <= i && i < N); return (*d_)[i]; }
void operator=(T t) {
for (size_t i=0; i!=N; ++i) {
d_[i] = t;
(*d_)[i] = t;
}
}
double const * begin() const { return &((*d_)[0]); }
double const * end() const { return &((*d_)[N]); }
private:
std::array<T,N> d_;
std::array<T,N>* d_;
std::shared_ptr< std::array<T,N> > do_;
};
template<size_t ROW, size_t COL, typename T=double>
class MatrixStatic {
using Data_t = std::array<std::array<T,COL>, ROW>;