1 //=============================================================================
2 //
3 // Copyright (c) Kitware, Inc.
4 // All rights reserved.
5 // See LICENSE.txt for details.
6 //
7 // This software is distributed WITHOUT ANY WARRANTY; without even
8 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
9 // PURPOSE. See the above copyright notice for more information.
10 //
11 // Copyright 2015 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
12 // Copyright 2015 UT-Battelle, LLC.
13 // Copyright 2015 Los Alamos National Security.
14 //
15 // Under the terms of Contract DE-NA0003525 with NTESS,
16 // the U.S. Government retains certain rights in this software.
17 // Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
18 // Laboratory (LANL), the U.S. Government retains certain rights in
19 // this software.
20 //
21 //=============================================================================
22 #ifndef vtk_m_testing_TestingMath_h
23 #define vtk_m_testing_TestingMath_h
24
25 #include <vtkm/Math.h>
26
27 #include <vtkm/TypeListTag.h>
28 #include <vtkm/VecTraits.h>
29
30 #include <vtkm/exec/FunctorBase.h>
31
32 #include <vtkm/cont/DeviceAdapterAlgorithm.h>
33
34 #include <vtkm/cont/testing/Testing.h>
35
36 #define VTKM_MATH_ASSERT(condition, message) \
37 if (!(condition)) \
38 { \
39 this->RaiseError(message); \
40 }
41
42 //-----------------------------------------------------------------------------
43 namespace UnitTestMathNamespace
44 {
45
46 class Lists
47 {
48 public:
49 static constexpr vtkm::IdComponent NUM_NUMBERS = 5;
50
NumberList(vtkm::Int32 i)51 VTKM_EXEC_CONT vtkm::Float64 NumberList(vtkm::Int32 i) const
52 {
53 vtkm::Float64 numberList[NUM_NUMBERS] = { 0.25, 0.5, 1.0, 2.0, 3.75 };
54 return numberList[i];
55 }
AngleList(vtkm::Int32 i)56 VTKM_EXEC_CONT vtkm::Float64 AngleList(vtkm::Int32 i) const
57 {
58 vtkm::Float64 angleList[NUM_NUMBERS] = { 0.643501108793284, // angle for 3, 4, 5 triangle.
59 0.78539816339745, // pi/4
60 0.5235987755983, // pi/6
61 1.0471975511966, // pi/3
62 0.0 };
63 return angleList[i];
64 }
OppositeList(vtkm::Int32 i)65 VTKM_EXEC_CONT vtkm::Float64 OppositeList(vtkm::Int32 i) const
66 {
67 vtkm::Float64 oppositeList[NUM_NUMBERS] = { 3.0, 1.0, 1.0, 1.732050807568877 /*sqrt(3)*/, 0.0 };
68 return oppositeList[i];
69 }
AdjacentList(vtkm::Int32 i)70 VTKM_EXEC_CONT vtkm::Float64 AdjacentList(vtkm::Int32 i) const
71 {
72 vtkm::Float64 adjacentList[NUM_NUMBERS] = { 4.0, 1.0, 1.732050807568877 /*sqrt(3)*/, 1.0, 1.0 };
73 return adjacentList[i];
74 }
HypotenuseList(vtkm::Int32 i)75 VTKM_EXEC_CONT vtkm::Float64 HypotenuseList(vtkm::Int32 i) const
76 {
77 vtkm::Float64 hypotenuseList[NUM_NUMBERS] = {
78 5.0, 1.414213562373095 /*sqrt(2)*/, 2.0, 2.0, 1.0
79 };
80 return hypotenuseList[i];
81 }
NumeratorList(vtkm::Int32 i)82 VTKM_EXEC_CONT vtkm::Float64 NumeratorList(vtkm::Int32 i) const
83 {
84 vtkm::Float64 numeratorList[NUM_NUMBERS] = { 6.5, 5.8, 9.3, 77.0, 0.1 };
85 return numeratorList[i];
86 }
DenominatorList(vtkm::Int32 i)87 VTKM_EXEC_CONT vtkm::Float64 DenominatorList(vtkm::Int32 i) const
88 {
89 vtkm::Float64 denominatorList[NUM_NUMBERS] = { 2.3, 1.6, 3.1, 19.0, 0.4 };
90 return denominatorList[i];
91 }
FModRemainderList(vtkm::Int32 i)92 VTKM_EXEC_CONT vtkm::Float64 FModRemainderList(vtkm::Int32 i) const
93 {
94 vtkm::Float64 fModRemainderList[NUM_NUMBERS] = { 1.9, 1.0, 0.0, 1.0, 0.1 };
95 return fModRemainderList[i];
96 }
RemainderList(vtkm::Int32 i)97 VTKM_EXEC_CONT vtkm::Float64 RemainderList(vtkm::Int32 i) const
98 {
99 vtkm::Float64 remainderList[NUM_NUMBERS] = { -0.4, -0.6, 0.0, 1.0, 0.1 };
100 return remainderList[i];
101 }
QuotientList(vtkm::Int32 i)102 VTKM_EXEC_CONT vtkm::Int64 QuotientList(vtkm::Int32 i) const
103 {
104 vtkm::Int64 quotientList[NUM_NUMBERS] = { 3, 4, 3, 4, 0 };
105 return quotientList[i];
106 }
XList(vtkm::Int32 i)107 VTKM_EXEC_CONT vtkm::Float64 XList(vtkm::Int32 i) const
108 {
109 vtkm::Float64 xList[NUM_NUMBERS] = { 4.6, 0.1, 73.4, 55.0, 3.75 };
110 return xList[i];
111 }
FractionalList(vtkm::Int32 i)112 VTKM_EXEC_CONT vtkm::Float64 FractionalList(vtkm::Int32 i) const
113 {
114 vtkm::Float64 fractionalList[NUM_NUMBERS] = { 0.6, 0.1, 0.4, 0.0, 0.75 };
115 return fractionalList[i];
116 }
FloorList(vtkm::Int32 i)117 VTKM_EXEC_CONT vtkm::Float64 FloorList(vtkm::Int32 i) const
118 {
119 vtkm::Float64 floorList[NUM_NUMBERS] = { 4.0, 0.0, 73.0, 55.0, 3.0 };
120 return floorList[i];
121 }
CeilList(vtkm::Int32 i)122 VTKM_EXEC_CONT vtkm::Float64 CeilList(vtkm::Int32 i) const
123 {
124 vtkm::Float64 ceilList[NUM_NUMBERS] = { 5.0, 1.0, 74.0, 55.0, 4.0 };
125 return ceilList[i];
126 }
RoundList(vtkm::Int32 i)127 VTKM_EXEC_CONT vtkm::Float64 RoundList(vtkm::Int32 i) const
128 {
129 vtkm::Float64 roundList[NUM_NUMBERS] = { 5.0, 0.0, 73.0, 55.0, 4.0 };
130 return roundList[i];
131 }
132 };
133
134 //-----------------------------------------------------------------------------
135 template <typename T>
136 struct ScalarFieldTests : public vtkm::exec::FunctorBase
137 {
138 VTKM_EXEC
TestPiScalarFieldTests139 void TestPi() const
140 {
141 // std::cout << "Testing Pi" << std::endl;
142 VTKM_MATH_ASSERT(test_equal(vtkm::Pi(), 3.14159265), "Pi not correct.");
143 VTKM_MATH_ASSERT(test_equal(vtkm::Pif(), 3.14159265f), "Pif not correct.");
144 VTKM_MATH_ASSERT(test_equal(vtkm::Pi<vtkm::Float64>(), 3.14159265),
145 "Pi template function not correct.");
146 }
147
148 VTKM_EXEC
TestArcTan2ScalarFieldTests149 void TestArcTan2() const
150 {
151 // std::cout << "Testing arc tan 2" << std::endl;
152
153 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(0.0), T(1.0)), T(0.0)), "ATan2 x+ axis.");
154 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(1.0), T(0.0)), T(0.5 * vtkm::Pi())),
155 "ATan2 y+ axis.");
156 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(-1.0), T(0.0)), T(-0.5 * vtkm::Pi())),
157 "ATan2 y- axis.");
158
159 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(1.0), T(1.0)), T(0.25 * vtkm::Pi())),
160 "ATan2 Quadrant 1");
161 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(1.0), T(-1.0)), T(0.75 * vtkm::Pi())),
162 "ATan2 Quadrant 2");
163 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(-1.0), T(-1.0)), T(-0.75 * vtkm::Pi())),
164 "ATan2 Quadrant 3");
165 VTKM_MATH_ASSERT(test_equal(vtkm::ATan2(T(-1.0), T(1.0)), T(-0.25 * vtkm::Pi())),
166 "ATan2 Quadrant 4");
167 }
168
169 VTKM_EXEC
TestPowScalarFieldTests170 void TestPow() const
171 {
172 // std::cout << "Running power tests." << std::endl;
173 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS; index++)
174 {
175 T x = static_cast<T>(Lists{}.NumberList(index));
176 T powx = vtkm::Pow(x, static_cast<T>(2.0));
177 T sqrx = x * x;
178 VTKM_MATH_ASSERT(test_equal(powx, sqrx), "Power gave wrong result.");
179 }
180 }
181
182 VTKM_EXEC
TestLog2ScalarFieldTests183 void TestLog2() const
184 {
185 // std::cout << "Testing Log2" << std::endl;
186 VTKM_MATH_ASSERT(test_equal(vtkm::Log2(T(0.25)), T(-2.0)), "Bad value from Log2");
187 VTKM_MATH_ASSERT(test_equal(vtkm::Log2(vtkm::Vec<T, 4>(0.5, 1.0, 2.0, 4.0)),
188 vtkm::Vec<T, 4>(-1.0, 0.0, 1.0, 2.0)),
189 "Bad value from Log2");
190 }
191
192 VTKM_EXEC
TestNonFinitesScalarFieldTests193 void TestNonFinites() const
194 {
195 // std::cout << "Testing non-finites." << std::endl;
196
197 T zero = 0.0;
198 T finite = 1.0;
199 T nan = vtkm::Nan<T>();
200 T inf = vtkm::Infinity<T>();
201 T neginf = vtkm::NegativeInfinity<T>();
202 T epsilon = vtkm::Epsilon<T>();
203
204 // General behavior.
205 VTKM_MATH_ASSERT(nan != vtkm::Nan<T>(), "Nan not equal itself.");
206 VTKM_MATH_ASSERT(!(nan >= zero), "Nan not greater or less.");
207 VTKM_MATH_ASSERT(!(nan <= zero), "Nan not greater or less.");
208 VTKM_MATH_ASSERT(!(nan >= finite), "Nan not greater or less.");
209 VTKM_MATH_ASSERT(!(nan <= finite), "Nan not greater or less.");
210
211 VTKM_MATH_ASSERT(neginf < inf, "Infinity big");
212 VTKM_MATH_ASSERT(zero < inf, "Infinity big");
213 VTKM_MATH_ASSERT(finite < inf, "Infinity big");
214 VTKM_MATH_ASSERT(zero > -inf, "-Infinity small");
215 VTKM_MATH_ASSERT(finite > -inf, "-Infinity small");
216 VTKM_MATH_ASSERT(zero > neginf, "-Infinity small");
217 VTKM_MATH_ASSERT(finite > neginf, "-Infinity small");
218
219 VTKM_MATH_ASSERT(zero < epsilon, "Negative epsilon");
220 VTKM_MATH_ASSERT(finite > epsilon, "Large epsilon");
221
222 // Math check functions.
223 VTKM_MATH_ASSERT(!vtkm::IsNan(zero), "Bad IsNan check.");
224 VTKM_MATH_ASSERT(!vtkm::IsNan(finite), "Bad IsNan check.");
225 VTKM_MATH_ASSERT(vtkm::IsNan(nan), "Bad IsNan check.");
226 VTKM_MATH_ASSERT(!vtkm::IsNan(inf), "Bad IsNan check.");
227 VTKM_MATH_ASSERT(!vtkm::IsNan(neginf), "Bad IsNan check.");
228 VTKM_MATH_ASSERT(!vtkm::IsNan(epsilon), "Bad IsNan check.");
229
230 VTKM_MATH_ASSERT(!vtkm::IsInf(zero), "Bad infinity check.");
231 VTKM_MATH_ASSERT(!vtkm::IsInf(finite), "Bad infinity check.");
232 VTKM_MATH_ASSERT(!vtkm::IsInf(nan), "Bad infinity check.");
233 VTKM_MATH_ASSERT(vtkm::IsInf(inf), "Bad infinity check.");
234 VTKM_MATH_ASSERT(vtkm::IsInf(neginf), "Bad infinity check.");
235 VTKM_MATH_ASSERT(!vtkm::IsInf(epsilon), "Bad infinity check.");
236
237 VTKM_MATH_ASSERT(vtkm::IsFinite(zero), "Bad finite check.");
238 VTKM_MATH_ASSERT(vtkm::IsFinite(finite), "Bad finite check.");
239 VTKM_MATH_ASSERT(!vtkm::IsFinite(nan), "Bad finite check.");
240 VTKM_MATH_ASSERT(!vtkm::IsFinite(inf), "Bad finite check.");
241 VTKM_MATH_ASSERT(!vtkm::IsFinite(neginf), "Bad finite check.");
242 VTKM_MATH_ASSERT(vtkm::IsFinite(epsilon), "Bad finite check.");
243 }
244
245 VTKM_EXEC
TestRemaindersScalarFieldTests246 void TestRemainders() const
247 {
248 // std::cout << "Testing remainders." << std::endl;
249 Lists table;
250 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS; index++)
251 {
252 T numerator = static_cast<T>(table.NumeratorList(index));
253 T denominator = static_cast<T>(table.DenominatorList(index));
254 T fmodremainder = static_cast<T>(table.FModRemainderList(index));
255 T remainder = static_cast<T>(table.RemainderList(index));
256 vtkm::Int64 quotient = table.QuotientList(index);
257
258 VTKM_MATH_ASSERT(test_equal(vtkm::FMod(numerator, denominator), fmodremainder),
259 "Bad FMod remainder.");
260 VTKM_MATH_ASSERT(test_equal(vtkm::Remainder(numerator, denominator), remainder),
261 "Bad remainder.");
262 vtkm::Int64 q;
263 VTKM_MATH_ASSERT(test_equal(vtkm::RemainderQuotient(numerator, denominator, q), remainder),
264 "Bad remainder-quotient remainder.");
265 VTKM_MATH_ASSERT(test_equal(q, quotient), "Bad reminder-quotient quotient.");
266 }
267 }
268
269 VTKM_EXEC
TestRoundScalarFieldTests270 void TestRound() const
271 {
272 // std::cout << "Testing round." << std::endl;
273 Lists table;
274 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS; index++)
275 {
276 T x = static_cast<T>(table.XList(index));
277 T fractional = static_cast<T>(table.FractionalList(index));
278 T floor = static_cast<T>(table.FloorList(index));
279 T ceil = static_cast<T>(table.CeilList(index));
280 T round = static_cast<T>(table.RoundList(index));
281
282 T intPart;
283 VTKM_MATH_ASSERT(test_equal(vtkm::ModF(x, intPart), fractional),
284 "ModF returned wrong fractional part.");
285 VTKM_MATH_ASSERT(test_equal(intPart, floor), "ModF returned wrong integral part.");
286 VTKM_MATH_ASSERT(test_equal(vtkm::Floor(x), floor), "Bad floor.");
287 VTKM_MATH_ASSERT(test_equal(vtkm::Ceil(x), ceil), "Bad ceil.");
288 VTKM_MATH_ASSERT(test_equal(vtkm::Round(x), round), "Bad round.");
289 }
290 }
291
292 VTKM_EXEC
TestIsNegativeScalarFieldTests293 void TestIsNegative() const
294 {
295 // std::cout << "Testing SignBit and IsNegative." << std::endl;
296 T x = 0;
297 VTKM_MATH_ASSERT(vtkm::SignBit(x) == 0, "SignBit wrong for 0.");
298 VTKM_MATH_ASSERT(!vtkm::IsNegative(x), "IsNegative wrong for 0.");
299
300 x = 20;
301 VTKM_MATH_ASSERT(vtkm::SignBit(x) == 0, "SignBit wrong for 20.");
302 VTKM_MATH_ASSERT(!vtkm::IsNegative(x), "IsNegative wrong for 20.");
303
304 x = -20;
305 VTKM_MATH_ASSERT(vtkm::SignBit(x) != 0, "SignBit wrong for -20.");
306 VTKM_MATH_ASSERT(vtkm::IsNegative(x), "IsNegative wrong for -20.");
307
308 x = 0.02f;
309 VTKM_MATH_ASSERT(vtkm::SignBit(x) == 0, "SignBit wrong for 0.02.");
310 VTKM_MATH_ASSERT(!vtkm::IsNegative(x), "IsNegative wrong for 0.02.");
311
312 x = -0.02f;
313 VTKM_MATH_ASSERT(vtkm::SignBit(x) != 0, "SignBit wrong for -0.02.");
314 VTKM_MATH_ASSERT(vtkm::IsNegative(x), "IsNegative wrong for -0.02.");
315 }
316
317 VTKM_EXEC
operatorScalarFieldTests318 void operator()(vtkm::Id) const
319 {
320 this->TestPi();
321 this->TestArcTan2();
322 this->TestPow();
323 this->TestLog2();
324 this->TestNonFinites();
325 this->TestRemainders();
326 this->TestRound();
327 this->TestIsNegative();
328 }
329 };
330
331 template <typename Device>
332 struct TryScalarFieldTests
333 {
334 template <typename T>
operatorTryScalarFieldTests335 void operator()(const T&) const
336 {
337 vtkm::cont::DeviceAdapterAlgorithm<Device>::Schedule(ScalarFieldTests<T>(), 1);
338 }
339 };
340
341 //-----------------------------------------------------------------------------
342 template <typename VectorType>
343 struct ScalarVectorFieldTests : public vtkm::exec::FunctorBase
344 {
345 using Traits = vtkm::VecTraits<VectorType>;
346 using ComponentType = typename Traits::ComponentType;
347 enum
348 {
349 NUM_COMPONENTS = Traits::NUM_COMPONENTS
350 };
351
352 VTKM_EXEC
TestTriangleTrigScalarVectorFieldTests353 void TestTriangleTrig() const
354 {
355 // std::cout << "Testing normal trig functions." << std::endl;
356 Lists table;
357 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS - NUM_COMPONENTS + 1; index++)
358 {
359 VectorType angle;
360 VectorType opposite;
361 VectorType adjacent;
362 VectorType hypotenuse;
363 for (vtkm::IdComponent componentIndex = 0; componentIndex < NUM_COMPONENTS; componentIndex++)
364 {
365 Traits::SetComponent(angle,
366 componentIndex,
367 static_cast<ComponentType>(table.AngleList(componentIndex + index)));
368 Traits::SetComponent(
369 opposite,
370 componentIndex,
371 static_cast<ComponentType>(table.OppositeList(componentIndex + index)));
372 Traits::SetComponent(
373 adjacent,
374 componentIndex,
375 static_cast<ComponentType>(table.AdjacentList(componentIndex + index)));
376 Traits::SetComponent(
377 hypotenuse,
378 componentIndex,
379 static_cast<ComponentType>(table.HypotenuseList(componentIndex + index)));
380 }
381
382 VTKM_MATH_ASSERT(test_equal(vtkm::Sin(angle), opposite / hypotenuse), "Sin failed test.");
383 VTKM_MATH_ASSERT(test_equal(vtkm::Cos(angle), adjacent / hypotenuse), "Cos failed test.");
384 VTKM_MATH_ASSERT(test_equal(vtkm::Tan(angle), opposite / adjacent), "Tan failed test.");
385
386 VTKM_MATH_ASSERT(test_equal(vtkm::ASin(opposite / hypotenuse), angle),
387 "Arc Sin failed test.");
388
389 #if defined(VTKM_ICC)
390 // When the intel compiler has vectorization enabled ( -O2/-O3 ) it converts the
391 // `adjacent/hypotenuse` divide operation into reciprocal (rcpps) and
392 // multiply (mulps) operations. This causes a change in the expected result that
393 // is larger than the default tolerance of test_equal.
394 //
395 VTKM_MATH_ASSERT(test_equal(vtkm::ACos(adjacent / hypotenuse), angle, 0.0004),
396 "Arc Cos failed test.");
397 #else
398 VTKM_MATH_ASSERT(test_equal(vtkm::ACos(adjacent / hypotenuse), angle),
399 "Arc Cos failed test.");
400 #endif
401 VTKM_MATH_ASSERT(test_equal(vtkm::ATan(opposite / adjacent), angle), "Arc Tan failed test.");
402 }
403 }
404
405 VTKM_EXEC
TestHyperbolicTrigScalarVectorFieldTests406 void TestHyperbolicTrig() const
407 {
408 // std::cout << "Testing hyperbolic trig functions." << std::endl;
409
410 const VectorType zero(0);
411 const VectorType half(0.5);
412 Lists table;
413 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS - NUM_COMPONENTS + 1; index++)
414 {
415 VectorType x;
416 for (vtkm::IdComponent componentIndex = 0; componentIndex < NUM_COMPONENTS; componentIndex++)
417 {
418 Traits::SetComponent(
419 x, componentIndex, static_cast<ComponentType>(table.AngleList(componentIndex + index)));
420 }
421
422 const VectorType minusX = zero - x;
423
424 VTKM_MATH_ASSERT(test_equal(vtkm::SinH(x), half * (vtkm::Exp(x) - vtkm::Exp(minusX))),
425 "SinH does not match definition.");
426 VTKM_MATH_ASSERT(test_equal(vtkm::CosH(x), half * (vtkm::Exp(x) + vtkm::Exp(minusX))),
427 "SinH does not match definition.");
428 VTKM_MATH_ASSERT(test_equal(vtkm::TanH(x), vtkm::SinH(x) / vtkm::CosH(x)),
429 "TanH does not match definition");
430
431 VTKM_MATH_ASSERT(test_equal(vtkm::ASinH(vtkm::SinH(x)), x), "SinH not inverting.");
432 VTKM_MATH_ASSERT(test_equal(vtkm::ACosH(vtkm::CosH(x)), x), "CosH not inverting.");
433 VTKM_MATH_ASSERT(test_equal(vtkm::ATanH(vtkm::TanH(x)), x), "TanH not inverting.");
434 }
435 }
436
437 template <typename FunctionType>
RaiseToTestScalarVectorFieldTests438 VTKM_EXEC void RaiseToTest(FunctionType function, ComponentType exponent) const
439 {
440 Lists table;
441 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS - NUM_COMPONENTS + 1; index++)
442 {
443 VectorType original;
444 VectorType raiseresult;
445 for (vtkm::IdComponent componentIndex = 0; componentIndex < NUM_COMPONENTS; componentIndex++)
446 {
447 ComponentType x = static_cast<ComponentType>(table.NumberList(componentIndex + index));
448 Traits::SetComponent(original, componentIndex, x);
449 Traits::SetComponent(raiseresult, componentIndex, vtkm::Pow(x, exponent));
450 }
451
452 VectorType mathresult = function(original);
453
454 VTKM_MATH_ASSERT(test_equal(mathresult, raiseresult), "Exponent functions do not agree.");
455 }
456 }
457
458 struct SqrtFunctor
459 {
460 VTKM_EXEC
operatorScalarVectorFieldTests::SqrtFunctor461 VectorType operator()(VectorType x) const { return vtkm::Sqrt(x); }
462 };
463 VTKM_EXEC
TestSqrtScalarVectorFieldTests464 void TestSqrt() const
465 {
466 // std::cout << " Testing Sqrt" << std::endl;
467 RaiseToTest(SqrtFunctor(), 0.5);
468 }
469
470 struct RSqrtFunctor
471 {
472 VTKM_EXEC
operatorScalarVectorFieldTests::RSqrtFunctor473 VectorType operator()(VectorType x) const { return vtkm::RSqrt(x); }
474 };
475 VTKM_EXEC
TestRSqrtScalarVectorFieldTests476 void TestRSqrt() const
477 {
478 // std::cout << " Testing RSqrt"<< std::endl;
479 RaiseToTest(RSqrtFunctor(), -0.5);
480 }
481
482 struct CbrtFunctor
483 {
484 VTKM_EXEC
operatorScalarVectorFieldTests::CbrtFunctor485 VectorType operator()(VectorType x) const { return vtkm::Cbrt(x); }
486 };
487 VTKM_EXEC
TestCbrtScalarVectorFieldTests488 void TestCbrt() const
489 {
490 // std::cout << " Testing Cbrt" << std::endl;
491 RaiseToTest(CbrtFunctor(), vtkm::Float32(1.0 / 3.0));
492 }
493
494 struct RCbrtFunctor
495 {
496 VTKM_EXEC
operatorScalarVectorFieldTests::RCbrtFunctor497 VectorType operator()(VectorType x) const { return vtkm::RCbrt(x); }
498 };
499 VTKM_EXEC
TestRCbrtScalarVectorFieldTests500 void TestRCbrt() const
501 {
502 // std::cout << " Testing RCbrt" << std::endl;
503 RaiseToTest(RCbrtFunctor(), vtkm::Float32(-1.0 / 3.0));
504 }
505
506 template <typename FunctionType>
507 VTKM_EXEC void RaiseByTest(FunctionType function,
508 ComponentType base,
509 ComponentType exponentbias = 0.0,
510 ComponentType resultbias = 0.0) const
511 {
512 Lists table;
513 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS - NUM_COMPONENTS + 1; index++)
514 {
515 VectorType original;
516 VectorType raiseresult;
517 for (vtkm::IdComponent componentIndex = 0; componentIndex < NUM_COMPONENTS; componentIndex++)
518 {
519 ComponentType x = static_cast<ComponentType>(table.NumberList(componentIndex + index));
520 Traits::SetComponent(original, componentIndex, x);
521 Traits::SetComponent(
522 raiseresult, componentIndex, vtkm::Pow(base, x + exponentbias) + resultbias);
523 }
524
525 VectorType mathresult = function(original);
526
527 VTKM_MATH_ASSERT(test_equal(mathresult, raiseresult), "Exponent functions do not agree.");
528 }
529 }
530
531 struct ExpFunctor
532 {
533 VTKM_EXEC
operatorScalarVectorFieldTests::ExpFunctor534 VectorType operator()(VectorType x) const { return vtkm::Exp(x); }
535 };
536 VTKM_EXEC
TestExpScalarVectorFieldTests537 void TestExp() const
538 {
539 // std::cout << " Testing Exp" << std::endl;
540 RaiseByTest(ExpFunctor(), vtkm::Float32(2.71828183));
541 }
542
543 struct Exp2Functor
544 {
545 VTKM_EXEC
operatorScalarVectorFieldTests::Exp2Functor546 VectorType operator()(VectorType x) const { return vtkm::Exp2(x); }
547 };
548 VTKM_EXEC
TestExp2ScalarVectorFieldTests549 void TestExp2() const
550 {
551 // std::cout << " Testing Exp2" << std::endl;
552 RaiseByTest(Exp2Functor(), 2.0);
553 }
554
555 struct ExpM1Functor
556 {
557 VTKM_EXEC
operatorScalarVectorFieldTests::ExpM1Functor558 VectorType operator()(VectorType x) const { return vtkm::ExpM1(x); }
559 };
560 VTKM_EXEC
TestExpM1ScalarVectorFieldTests561 void TestExpM1() const
562 {
563 // std::cout << " Testing ExpM1" << std::endl;
564 RaiseByTest(ExpM1Functor(), ComponentType(2.71828183), 0.0, -1.0);
565 }
566
567 struct Exp10Functor
568 {
569 VTKM_EXEC
operatorScalarVectorFieldTests::Exp10Functor570 VectorType operator()(VectorType x) const { return vtkm::Exp10(x); }
571 };
572 VTKM_EXEC
TestExp10ScalarVectorFieldTests573 void TestExp10() const
574 {
575 // std::cout << " Testing Exp10" << std::endl;
576 RaiseByTest(Exp10Functor(), 10.0);
577 }
578
579 template <typename FunctionType>
580 VTKM_EXEC void LogBaseTest(FunctionType function,
581 ComponentType base,
582 ComponentType bias = 0.0) const
583 {
584 Lists table;
585 for (vtkm::IdComponent index = 0; index < Lists::NUM_NUMBERS - NUM_COMPONENTS + 1; index++)
586 {
587 VectorType basevector(base);
588 VectorType original;
589 VectorType biased;
590 for (vtkm::IdComponent componentIndex = 0; componentIndex < NUM_COMPONENTS; componentIndex++)
591 {
592 ComponentType x = static_cast<ComponentType>(table.NumberList(componentIndex + index));
593 Traits::SetComponent(original, componentIndex, x);
594 Traits::SetComponent(biased, componentIndex, x + bias);
595 }
596
597 VectorType logresult = vtkm::Log2(biased) / vtkm::Log2(basevector);
598
599 VectorType mathresult = function(original);
600
601 VTKM_MATH_ASSERT(test_equal(mathresult, logresult), "Exponent functions do not agree.");
602 }
603 }
604
605 struct LogFunctor
606 {
607 VTKM_EXEC
operatorScalarVectorFieldTests::LogFunctor608 VectorType operator()(VectorType x) const { return vtkm::Log(x); }
609 };
610 VTKM_EXEC
TestLogScalarVectorFieldTests611 void TestLog() const
612 {
613 // std::cout << " Testing Log" << std::endl;
614 LogBaseTest(LogFunctor(), vtkm::Float32(2.71828183));
615 }
616
617 struct Log10Functor
618 {
619 VTKM_EXEC
operatorScalarVectorFieldTests::Log10Functor620 VectorType operator()(VectorType x) const { return vtkm::Log10(x); }
621 };
622 VTKM_EXEC
TestLog10ScalarVectorFieldTests623 void TestLog10() const
624 {
625 // std::cout << " Testing Log10" << std::endl;
626 LogBaseTest(Log10Functor(), 10.0);
627 }
628
629 struct Log1PFunctor
630 {
631 VTKM_EXEC
operatorScalarVectorFieldTests::Log1PFunctor632 VectorType operator()(VectorType x) const { return vtkm::Log1P(x); }
633 };
634 VTKM_EXEC
TestLog1PScalarVectorFieldTests635 void TestLog1P() const
636 {
637 // std::cout << " Testing Log1P" << std::endl;
638 LogBaseTest(Log1PFunctor(), ComponentType(2.71828183), 1.0);
639 }
640
641 VTKM_EXEC
TestCopySignScalarVectorFieldTests642 void TestCopySign() const
643 {
644 // std::cout << "Testing CopySign." << std::endl;
645 // Assuming all TestValues positive.
646 VectorType positive1 = TestValue(1, VectorType());
647 VectorType positive2 = TestValue(2, VectorType());
648 VectorType negative1 = -positive1;
649 VectorType negative2 = -positive2;
650
651 VTKM_MATH_ASSERT(test_equal(vtkm::CopySign(positive1, positive2), positive1),
652 "CopySign failed.");
653 VTKM_MATH_ASSERT(test_equal(vtkm::CopySign(negative1, positive2), positive1),
654 "CopySign failed.");
655 VTKM_MATH_ASSERT(test_equal(vtkm::CopySign(positive1, negative2), negative1),
656 "CopySign failed.");
657 VTKM_MATH_ASSERT(test_equal(vtkm::CopySign(negative1, negative2), negative1),
658 "CopySign failed.");
659 }
660
661 VTKM_EXEC
operatorScalarVectorFieldTests662 void operator()(vtkm::Id) const
663 {
664 this->TestTriangleTrig();
665 this->TestHyperbolicTrig();
666 this->TestSqrt();
667 this->TestRSqrt();
668 this->TestCbrt();
669 this->TestRCbrt();
670 this->TestExp();
671 this->TestExp2();
672 this->TestExpM1();
673 this->TestExp10();
674 this->TestLog();
675 this->TestLog10();
676 this->TestLog1P();
677 this->TestCopySign();
678 }
679 };
680
681 template <typename Device>
682 struct TryScalarVectorFieldTests
683 {
684 template <typename VectorType>
operatorTryScalarVectorFieldTests685 void operator()(const VectorType&) const
686 {
687 vtkm::cont::DeviceAdapterAlgorithm<Device>::Schedule(ScalarVectorFieldTests<VectorType>(), 1);
688 }
689 };
690
691 //-----------------------------------------------------------------------------
692 template <typename T>
693 struct AllTypesTests : public vtkm::exec::FunctorBase
694 {
695 VTKM_EXEC
TestMinMaxAllTypesTests696 void TestMinMax() const
697 {
698 T low = TestValue(2, T());
699 T high = TestValue(10, T());
700 // std::cout << "Testing min/max " << low << " " << high << std::endl;
701 VTKM_MATH_ASSERT(test_equal(vtkm::Min(low, high), low), "Wrong min.");
702 VTKM_MATH_ASSERT(test_equal(vtkm::Min(high, low), low), "Wrong min.");
703 VTKM_MATH_ASSERT(test_equal(vtkm::Max(low, high), high), "Wrong max.");
704 VTKM_MATH_ASSERT(test_equal(vtkm::Max(high, low), high), "Wrong max.");
705
706 using Traits = vtkm::VecTraits<T>;
707 T mixed1 = low;
708 T mixed2 = high;
709 Traits::SetComponent(mixed1, 0, Traits::GetComponent(high, 0));
710 Traits::SetComponent(mixed2, 0, Traits::GetComponent(low, 0));
711 VTKM_MATH_ASSERT(test_equal(vtkm::Min(mixed1, mixed2), low), "Wrong min.");
712 VTKM_MATH_ASSERT(test_equal(vtkm::Min(mixed2, mixed1), low), "Wrong min.");
713 VTKM_MATH_ASSERT(test_equal(vtkm::Max(mixed1, mixed2), high), "Wrong max.");
714 VTKM_MATH_ASSERT(test_equal(vtkm::Max(mixed2, mixed1), high), "Wrong max.");
715 }
716
717 VTKM_EXEC
operatorAllTypesTests718 void operator()(vtkm::Id) const { this->TestMinMax(); }
719 };
720
721 template <typename Device>
722 struct TryAllTypesTests
723 {
724 template <typename T>
operatorTryAllTypesTests725 void operator()(const T&) const
726 {
727 vtkm::cont::DeviceAdapterAlgorithm<Device>::Schedule(AllTypesTests<T>(), 1);
728 }
729 };
730
731 //-----------------------------------------------------------------------------
732 template <typename T>
733 struct AbsTests : public vtkm::exec::FunctorBase
734 {
735 VTKM_EXEC
operatorAbsTests736 void operator()(vtkm::Id index) const
737 {
738 // std::cout << "Testing Abs." << std::endl;
739 T positive = TestValue(index, T()); // Assuming all TestValues positive.
740 T negative = -positive;
741
742 VTKM_MATH_ASSERT(test_equal(vtkm::Abs(positive), positive), "Abs returned wrong value.");
743 VTKM_MATH_ASSERT(test_equal(vtkm::Abs(negative), positive), "Abs returned wrong value.");
744 }
745 };
746
747 template <typename Device>
748 struct TryAbsTests
749 {
750 template <typename T>
operatorTryAbsTests751 void operator()(const T&) const
752 {
753 vtkm::cont::DeviceAdapterAlgorithm<Device>::Schedule(AbsTests<T>(), 10);
754 }
755 };
756
757 struct TypeListTagAbs
758 : vtkm::ListTagJoin<
759 vtkm::ListTagJoin<vtkm::ListTagBase<vtkm::Int32, vtkm::Int64>, vtkm::TypeListTagIndex>,
760 vtkm::TypeListTagField>
761 {
762 };
763
764 //-----------------------------------------------------------------------------
765 template <typename Device>
RunMathTests()766 void RunMathTests()
767 {
768 std::cout << "Tests for scalar types." << std::endl;
769 vtkm::testing::Testing::TryTypes(TryScalarFieldTests<Device>(), vtkm::TypeListTagFieldScalar());
770 std::cout << "Test for scalar and vector types." << std::endl;
771 vtkm::testing::Testing::TryTypes(TryScalarVectorFieldTests<Device>(), vtkm::TypeListTagField());
772 std::cout << "Test for exemplar types." << std::endl;
773 vtkm::testing::Testing::TryTypes(TryAllTypesTests<Device>());
774 std::cout << "Test all Abs types" << std::endl;
775 vtkm::testing::Testing::TryTypes(TryAbsTests<Device>(), TypeListTagAbs());
776 }
777
778 } // namespace UnitTestMathNamespace
779
780 #endif //vtk_m_testing_TestingMath_h
781