Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALInference 4.1.9.1-8a6b96f
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
test_cubic_interp.c
Go to the documentation of this file.
1/*
2 * Copyright (C) 2015-2017 Leo Singer
3 *
4 * This program is free software; you can redistribute it and/or modify
5 * it under the terms of the GNU General Public License as published by
6 * the Free Software Foundation; either version 2 of the License, or
7 * (at your option) any later version.
8 *
9 * This program is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
13 *
14 * You should have received a copy of the GNU General Public License
15 * along with with program; see the file COPYING. If not, write to the
16 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17 * MA 02110-1301 USA
18 */
19
20
21#include <lal/cubic_interp.h>
22#include <gsl/gsl_test.h>
23#include <gsl/gsl_math.h>
24
25#include <lal/XLALError.h>
26
27#ifdef __GNUC__
28#define UNUSED __attribute__ ((unused))
29#else
30#define UNUSED
31#endif
32
33
34int main(int UNUSED argc, char UNUSED **argv)
35{
36 {
37 static const double data[] = {0, 0, 0, 0};
38 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
39 XLAL_CHECK_EXIT(interp);
40 for (double t = 0; t <= 1; t += 0.01)
41 {
42 const double result = cubic_interp_eval(interp, t);
43 const double expected = 0;
44 gsl_test_abs(result, expected, 0,
45 "testing cubic interpolant for zero input");
46 }
47 cubic_interp_free(interp);
48 }
49
50 {
51 static const double data[] = {1, 1, 1, 1};
52 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
53 XLAL_CHECK_EXIT(interp);
54 for (double t = 0; t <= 1; t += 0.01)
55 {
56 const double result = cubic_interp_eval(interp, t);
57 const double expected = 1;
58 gsl_test_abs(result, expected, 0,
59 "testing cubic interpolant for unit input");
60 }
61 cubic_interp_free(interp);
62 }
63
64 {
65 static const double data[] = {1, 0, 1, 4};
66 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
67 XLAL_CHECK_EXIT(interp);
68 for (double t = 0; t <= 1; t += 0.01)
69 {
70 const double result = cubic_interp_eval(interp, t);
71 const double expected = gsl_pow_2(t);
72 gsl_test_abs(result, expected, 10 * GSL_DBL_EPSILON,
73 "testing cubic interpolant for quadratic input");
74 }
75 cubic_interp_free(interp);
76 }
77
78 {
79 static const double data[] = {
80 GSL_POSINF, GSL_POSINF, GSL_POSINF, GSL_POSINF};
81 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
82 XLAL_CHECK_EXIT(interp);
83 for (double t = 0; t <= 1; t += 0.01)
84 {
85 const double result = cubic_interp_eval(interp, t);
86 const double expected = GSL_POSINF;
87 gsl_test_abs(result, expected, 0,
88 "testing cubic interpolant for +inf input");
89 }
90 cubic_interp_free(interp);
91 }
92
93 {
94 static const double data[] = {
95 0, GSL_POSINF, GSL_POSINF, GSL_POSINF};
96 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
97 XLAL_CHECK_EXIT(interp);
98 for (double t = 0; t <= 1; t += 0.01)
99 {
100 const double result = cubic_interp_eval(interp, t);
101 const double expected = GSL_POSINF;
102 gsl_test_abs(result, expected, 0,
103 "testing cubic interpolant for +inf input");
104 }
105 cubic_interp_free(interp);
106 }
107
108 {
109 static const double data[] = {
110 GSL_POSINF, GSL_POSINF, GSL_POSINF, 0};
111 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
112 XLAL_CHECK_EXIT(interp);
113 for (double t = 0; t <= 1; t += 0.01)
114 {
115 const double result = cubic_interp_eval(interp, t);
116 const double expected = GSL_POSINF;
117 gsl_test_abs(result, expected, 0,
118 "testing cubic interpolant for +inf input");
119 }
120 cubic_interp_free(interp);
121 }
122
123 {
124 static const double data[] = {
125 0, GSL_POSINF, GSL_POSINF, 0};
126 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
127 XLAL_CHECK_EXIT(interp);
128 for (double t = 0; t <= 1; t += 0.01)
129 {
130 const double result = cubic_interp_eval(interp, t);
131 const double expected = GSL_POSINF;
132 gsl_test_abs(result, expected, 0,
133 "testing cubic interpolant for +inf input");
134 }
135 cubic_interp_free(interp);
136 }
137
138 {
139 static const double data[] = {
140 0, 0, GSL_POSINF, 0};
141 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
142 XLAL_CHECK_EXIT(interp);
143 for (double t = 0.01; t <= 1; t += 0.01)
144 {
145 const double result = cubic_interp_eval(interp, t);
146 const double expected = 0;
147 gsl_test_abs(result, expected, 0,
148 "testing cubic interpolant for +inf input");
149 }
150 cubic_interp_free(interp);
151 }
152
153 {
154 static const double data[] = {
155 0, GSL_NEGINF, GSL_POSINF, 0};
156 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
157 XLAL_CHECK_EXIT(interp);
158 const double result = cubic_interp_eval(interp, 1);
159 cubic_interp_free(interp);
160 const double expected = GSL_POSINF;
161 gsl_test_abs(result, expected, 0,
162 "testing cubic interpolant for +inf input");
163 }
164
165 {
166 static const double data[] = {
167 0, GSL_POSINF, GSL_NEGINF, 0};
168 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
169 XLAL_CHECK_EXIT(interp);
170 const double result = cubic_interp_eval(interp, 0);
171 cubic_interp_free(interp);
172 const double expected = GSL_POSINF;
173 gsl_test_abs(result, expected, 0,
174 "testing cubic interpolant for +inf input");
175 }
176
177 {
178 static const double data[] = {
179 0, GSL_NEGINF, GSL_NEGINF, GSL_NEGINF};
180 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
181 XLAL_CHECK_EXIT(interp);
182 for (double t = 0; t <= 1; t += 0.01)
183 {
184 const double result = cubic_interp_eval(interp, t);
185 const double expected = GSL_NEGINF;
186 gsl_test_abs(result, expected, 0,
187 "testing cubic interpolant for -inf input");
188 }
189 cubic_interp_free(interp);
190 }
191
192 {
193 static const double data[] = {
194 GSL_NEGINF, GSL_NEGINF, GSL_NEGINF, 0};
195 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
196 XLAL_CHECK_EXIT(interp);
197 for (double t = 0; t <= 1; t += 0.01)
198 {
199 const double result = cubic_interp_eval(interp, t);
200 const double expected = GSL_NEGINF;
201 gsl_test_abs(result, expected, 0,
202 "testing cubic interpolant for -inf input");
203 }
204 cubic_interp_free(interp);
205 }
206
207 {
208 static const double data[] = {
209 0, GSL_NEGINF, GSL_NEGINF, 0};
210 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
211 XLAL_CHECK_EXIT(interp);
212 for (double t = 0; t <= 1; t += 0.01)
213 {
214 const double result = cubic_interp_eval(interp, t);
215 const double expected = GSL_NEGINF;
216 gsl_test_abs(result, expected, 0,
217 "testing cubic interpolant for -inf input");
218 }
219 cubic_interp_free(interp);
220 }
221
222 {
223 static const double data[] = {
224 0, 0, GSL_NEGINF, 0};
225 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
226 XLAL_CHECK_EXIT(interp);
227 for (double t = 0.01; t <= 1; t += 0.01)
228 {
229 const double result = cubic_interp_eval(interp, t);
230 const double expected = 0;
231 gsl_test_abs(result, expected, 0,
232 "testing cubic interpolant for -inf input");
233 }
234 cubic_interp_free(interp);
235 }
236
237 {
238 static const double data[] = {
239 0, GSL_NEGINF, GSL_POSINF, 0};
240 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
241 XLAL_CHECK_EXIT(interp);
242 const double result = cubic_interp_eval(interp, 0);
243 cubic_interp_free(interp);
244 const double expected = GSL_NEGINF;
245 gsl_test_abs(result, expected, 0,
246 "testing cubic interpolant for -inf input");
247 }
248
249 {
250 static const double data[] = {
251 0, GSL_POSINF, GSL_NEGINF, 0};
252 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
253 XLAL_CHECK_EXIT(interp);
254 const double result = cubic_interp_eval(interp, 1);
255 cubic_interp_free(interp);
256 const double expected = GSL_NEGINF;
257 gsl_test_abs(result, expected, 0,
258 "testing cubic interpolant for -inf input");
259 }
260
261 {
262 static const double data[] = {
263 0, GSL_NEGINF, GSL_POSINF, 0};
264 cubic_interp *interp = cubic_interp_init(data, 4, -1, 1);
265 XLAL_CHECK_EXIT(interp);
266 for (double t = 0.01; t < 1; t += 0.01)
267 {
268 const double result = cubic_interp_eval(interp, t);
269 const double expected = GSL_NEGINF;
270 gsl_test_abs(result, expected, 0,
271 "testing cubic interpolant for indeterminate input");
272 }
273 cubic_interp_free(interp);
274 }
275
276
277 {
278 static const double constants[] = {
279 0, 1, GSL_POSINF, GSL_NEGINF, GSL_NAN};
280 for (unsigned k = 0; k < sizeof(constants) / sizeof(*constants); k ++)
281 {
282 double data[4][4];
283 for (int i = 0; i < 4; i ++)
284 for (int j = 0; j < 4; j ++)
285 data[i][j] = constants[k];
287 *data, 4, 4, -1, -1, 1, 1);
288 for (double s = -5; s <= 2; s += 0.1)
289 {
290 for (double t = -5; t <= 1; t += 0.1)
291 {
292 const double result = bicubic_interp_eval(interp, s, t);
293 const double expected = constants[k];
294 gsl_test_abs(result, expected, 0,
295 "testing bicubic interpolant for constant %g input",
296 constants[k]);
297 }
298 }
299 XLAL_CHECK_EXIT(interp);
300 bicubic_interp_free(interp);
301 }
302 }
303
304 for (int k = 1; k < 3; k ++)
305 {
306 {
307 double data[4][4];
308 for (int i = 0; i < 4; i ++)
309 for (int j = 0; j < 4; j ++)
310 data[i][j] = gsl_pow_int(i - 1, k);
312 *data, 4, 4, -1, -1, 1, 1);
313 for (double s = 0; s <= 1; s += 0.1)
314 {
315 for (double t = 0; t <= 1; t += 0.1)
316 {
317 const double result = bicubic_interp_eval(interp, s, t);
318 const double expected = gsl_pow_int(s, k);
319 gsl_test_abs(result, expected, 10 * GSL_DBL_EPSILON,
320 "testing bicubic interpolant for s^%d input", k);
321 }
322 }
323 XLAL_CHECK_EXIT(interp);
324 bicubic_interp_free(interp);
325 }
326
327 {
328 double data[4][4];
329 for (int i = 0; i < 4; i ++)
330 for (int j = 0; j < 4; j ++)
331 data[i][j] = gsl_pow_int(j - 1, k);
333 *data, 4, 4, -1, -1, 1, 1);
334 for (double s = 0; s <= 1; s += 0.1)
335 {
336 for (double t = 0; t <= 1; t += 0.1)
337 {
338 const double result = bicubic_interp_eval(interp, s, t);
339 const double expected = gsl_pow_int(t, k);
340 gsl_test_abs(result, expected, 10 * GSL_DBL_EPSILON,
341 "testing bicubic interpolant for t^%d input", k);
342 }
343 }
344 XLAL_CHECK_EXIT(interp);
345 bicubic_interp_free(interp);
346 }
347
348 {
349 double data[4][4];
350 for (int i = 0; i < 4; i ++)
351 for (int j = 0; j < 4; j ++)
352 data[i][j] = gsl_pow_int(i - 1, k) + gsl_pow_int(j - 1, k);
354 *data, 4, 4, -1, -1, 1, 1);
355 for (double s = 0; s <= 1; s += 0.1)
356 {
357 for (double t = 0; t <= 1; t += 0.1)
358 {
359 const double result = bicubic_interp_eval(interp, s, t);
360 const double expected = gsl_pow_int(s, k)
361 + gsl_pow_int(t, k);
362 gsl_test_abs(result, expected, 10 * GSL_DBL_EPSILON,
363 "testing bicubic interpolant for s^%d + t^%d input",
364 k, k);
365 }
366 }
367 XLAL_CHECK_EXIT(interp);
368 bicubic_interp_free(interp);
369 }
370 }
371
372 return gsl_test_summary();
373}
int j
int k
int s
bicubic_interp * bicubic_interp_init(const double *data, int ns, int nt, double smin, double tmin, double ds, double dt)
Definition: cubic_interp.c:136
void cubic_interp_free(cubic_interp *interp)
Definition: cubic_interp.c:120
void bicubic_interp_free(bicubic_interp *interp)
Definition: cubic_interp.c:187
double bicubic_interp_eval(const bicubic_interp *interp, double s, double t)
Definition: cubic_interp.c:193
cubic_interp * cubic_interp_init(const double *data, int n, double tmin, double dt)
Definition: cubic_interp.c:95
double cubic_interp_eval(const cubic_interp *interp, double t)
Definition: cubic_interp.c:126
sigmaKerr data[0]
#define XLAL_CHECK_EXIT(assertion)
int main(int UNUSED argc, char UNUSED **argv)