1 ////////////////////////////////////////////////////////////////////////////////////// 2 // This file is distributed under the University of Illinois/NCSA Open Source License. 3 // See LICENSE file in top directory for details. 4 // 5 // Copyright (c) 2019 QMCPACK developers. 6 // 7 // File developed by: Peter Doak, doakpw@ornl.gov, Oak Ridge National Lab 8 // Mark Dewing, markdewing@gmail.com, University of Illinois at Urbana-Champaign 9 // 10 // File created by: Mark Dewing, markdewing@gmail.com, University of Illinois at Urbana-Champaign 11 ////////////////////////////////////////////////////////////////////////////////////// 12 13 14 #include "catch.hpp" 15 16 #include "einspline/bspline_base.h" 17 #include "einspline/bspline_create.h" 18 #include "einspline/bspline_eval_d.h" 19 #include "einspline/multi_bspline_create.h" 20 #include "einspline/multi_bspline_eval_d.h" 21 #include "einspline/multi_bspline_eval_s.h" 22 23 #include <stdio.h> 24 #include <vector> 25 #include "config/stdlib/Constants.h" 26 27 TEST_CASE("double_1d_natural", "[einspline]") 28 { 29 // two point case 30 Ugrid x_grid; 31 x_grid.start = 1.0; 32 x_grid.end = 10.0; 33 x_grid.num = 2; 34 35 std::vector<double> data(2); 36 data[0] = 2.0; data[1] = 3.0; 37 38 BCtype_d xBC; 39 xBC.lCode = NATURAL; 40 xBC.rCode = NATURAL; 41 UBspline_1d_d* s = create_UBspline_1d_d(x_grid, xBC, data.data()); 42 43 REQUIRE(s); 44 45 double val; 46 47 eval_UBspline_1d_d(s, 1.0, &val); 48 REQUIRE(val == Approx(2.0)); 49 50 eval_UBspline_1d_d(s, 9.9999999, &val); 51 REQUIRE(val == Approx(3.0)); 52 53 // This should assert 54 // eval_UBspline_1d_d(s, 10.0, &val); 55 // REQUIRE(val == Approx(3.0)); 56 57 eval_UBspline_1d_d(s, 5.5, &val); 58 REQUIRE(val == Approx(2.5)); 59 60 destroy_Bspline(s); 61 62 // three point case 63 x_grid.start = 1.0; 64 x_grid.end = 10.0; 65 x_grid.num = 3; 66 67 data.resize(3); 68 data[0] = 2.0; data[1] = 2.7; data[2] = 3.0; 69 70 xBC.lCode = NATURAL; 71 xBC.rCode = NATURAL; 72 s = create_UBspline_1d_d(x_grid, xBC, data.data()); 73 74 REQUIRE(s); 75 76 eval_UBspline_1d_d(s, 1.0, &val); 77 REQUIRE(val == Approx(2.0)); 78 79 eval_UBspline_1d_d(s, 9.9999999, &val); 80 REQUIRE(val == Approx(3.0)); 81 82 eval_UBspline_1d_d(s, 5.5, &val); 83 REQUIRE(val == Approx(2.7)); 84 85 destroy_Bspline(s); 86 87 } 88 89 90 TEST_CASE("double_1d_multi", "[einspline]") 91 { 92 Ugrid x_grid; 93 x_grid.start = 1.0; 94 x_grid.end = 10.0; 95 x_grid.num = 2; 96 97 double data[2]; 98 data[0] = 2.0; 99 data[1] = 3.0; 100 101 BCtype_d xBC; 102 xBC.lCode = NATURAL; 103 xBC.rCode = NATURAL; 104 multi_UBspline_1d_d* s = create_multi_UBspline_1d_d(x_grid, xBC, 1); 105 REQUIRE(s); 106 107 set_multi_UBspline_1d_d(s, 0, data); 108 109 double val; 110 eval_multi_UBspline_1d_d(s, 1.0, &val); 111 REQUIRE(val == Approx(2.0)); 112 } 113 114 TEST_CASE("double_1d_periodic", "[einspline]") 115 { 116 Ugrid x_grid; 117 x_grid.start = 0.0; 118 x_grid.end = 1.0; 119 //Enough grid points are required to do the micro evaluation test. 120 constexpr int N = 12; 121 x_grid.num = N; 122 double delta = (x_grid.end - x_grid.start)/x_grid.num; 123 124 double tpi = 2*M_PI; 125 double data[N]; 126 for (int i = 0; i < N; i++) 127 { 128 double x = delta*i; 129 data[i] = sin(tpi*x); 130 } 131 132 BCtype_d bc; 133 bc.lCode = PERIODIC; 134 bc.rCode = PERIODIC; 135 136 UBspline_1d_d* s = create_UBspline_1d_d(x_grid, bc, data); 137 138 REQUIRE(s); 139 140 double val; 141 eval_UBspline_1d_d(s, 0.0, &val); 142 REQUIRE(val == Approx(0.0)); 143 144 eval_UBspline_1d_d(s, delta, &val); 145 REQUIRE(val == Approx(data[1])); 146 147 double micro_delta = delta / 4.0; 148 int micro_N = N * 4; 149 double micro_data[N*4]; 150 for (int i = 0 ; i < micro_N; i++) 151 { 152 double x = micro_delta * i; 153 micro_data[i] = sin(tpi*x); 154 } 155 eval_UBspline_1d_d(s, micro_delta * 3, &val); 156 REQUIRE(val == Approx(micro_data[3]).epsilon(0.001)); 157 158 eval_UBspline_1d_d(s, micro_delta * 17, &val); 159 REQUIRE(val == Approx(micro_data[17]).epsilon(0.001)); 160 161 eval_UBspline_1d_d(s, micro_delta * 31, &val); 162 REQUIRE(val == Approx(micro_data[31]).epsilon(0.001)); 163 164 destroy_Bspline(s); 165 } 166 167 #ifdef QMC_CUDA 168 void test_multi(multi_UBspline_1d_s *cpuSpline, float *pos, float *vals_cuda); 169 170 TEST_CASE("multi_cuda_wrapper", "[einspline]") 171 { 172 Ugrid x_grid; 173 // GPU versions require the grid to start at zero 174 x_grid.start = 0.0; 175 x_grid.end = 10.0; 176 x_grid.num = 2; 177 178 float data[2]; 179 data[0] = 2.0; 180 data[1] = 3.0; 181 182 BCtype_s xBC; 183 xBC.lCode = NATURAL; 184 xBC.rCode = NATURAL; 185 multi_UBspline_1d_s* s = create_multi_UBspline_1d_s(x_grid, xBC, 1); 186 REQUIRE(s); 187 set_multi_UBspline_1d_s(s, 0, data); 188 189 float pos[1]; 190 pos[0] = 0.0; 191 192 // Check the CPU value 193 float cpu_val[1]; 194 eval_multi_UBspline_1d_s(s, pos[0], cpu_val); 195 REQUIRE(cpu_val[0] == 2.0); 196 197 //this would assert in debug and is an illegal value for pos[0] 198 //pos[0] = 11.0; 199 pos[0] = 9.99999999; 200 // Check the CPU value 201 eval_multi_UBspline_1d_s(s, pos[0], cpu_val); 202 //std::cout << std::setprecision(24) << cpu_val[0] << " == " << 3.0000f << '\n'; 203 //With power 9/ clang 8 3.0000f is 3 but cpu_val[0] = 3.0000002384185791015625 204 REQUIRE(cpu_val[0] == Approx(3.0000f).epsilon(0.00000025)); 205 206 // Check the GPU value 207 pos[0] = 0.0; 208 float vals_output[1]; 209 test_multi(s, pos, vals_output); 210 REQUIRE(vals_output[0] == 2.0); 211 212 pos[0] = 9.99999999; 213 test_multi(s, pos, vals_output); 214 REQUIRE(vals_output[0] == 3.0); 215 } 216 217 #endif 218