src/ImathTest/testMatrix.cpp

   1 //
   2 // SPDX-License-Identifier: BSD-3-Clause
   3 // Copyright Contributors to the OpenEXR Project.
   4 //
   5
   6 #ifdef NDEBUG
   7 #    undef NDEBUG
   8 #endif
   9
  10 #include <ImathInt64.h>
  11 #include <ImathMath.h>
  12 #include <ImathMatrix.h>
  13 #include <ImathMatrixAlgo.h>
  14 #include <ImathRandom.h>
  15 #include <ImathVec.h>
  16 #include <assert.h>
  17 #include <iostream>
  18 #include "testMatrix.h"
  19 #include <sstream>
  20
  21 // Include ImathForward *after* other headers to validate forward declarations
  22 #include <ImathForward.h>
  23
  24 using namespace std;
  25 using IMATH_INTERNAL_NAMESPACE::Int64;
  26
  27 //
  28 // This file is not currently intended to exhaustively test
  29 // the Imath Matrix33<T> and Matrix44<T> classes.  We leave
  30 // that to PyImathTest.
  31 //
  32 // Instead, in this file we test only those aspects of the
  33 // Imath Matrix33<T> and Matrix44<T> classes that must be
  34 // or are more convenient to test from C++.
  35 //
  36
  37 void
  38 testMatrix22ArrayConstructor(const float a[2][2])
  39 {
  40     IMATH_INTERNAL_NAMESPACE::M22f m(a);
  41     assert(m == IMATH_INTERNAL_NAMESPACE::M22f());
  42 }
  43
  44 void
  45 testMatrix33ArrayConstructor(const float a[3][3])
  46 {
  47     IMATH_INTERNAL_NAMESPACE::M33f m(a);
  48     assert(m == IMATH_INTERNAL_NAMESPACE::M33f());
  49 }
  50
  51 void
  52 testMatrix44ArrayConstructor(const float a[4][4])
  53 {
  54     IMATH_INTERNAL_NAMESPACE::M44f m(a);
  55     assert(m == IMATH_INTERNAL_NAMESPACE::M44f());
  56 }
  57
  58
  59 void
  60 testMatrix ()
  61 {
  62     cout << "Testing functions in ImathMatrix.h" << endl;
  63
  64     union
  65     {
  66         float f;
  67         int i;
  68     } nanf;
  69     nanf.i = 0x7f800001; //  NAN
  70
  71     union
  72     {
  73         double d;
  74         uint64_t i;
  75     } nand;
  76     nand.i = 0x7ff0000000000001ULL; //  NAN
  77
  78     {
  79         cout << "Imath::M22f constructors and equality operators" << endl;
  80
  81         IMATH_INTERNAL_NAMESPACE::M22f m1;
  82         m1[0][0] = 99.0f;
  83         m1[1][1] = 101.0f;
  84
  85         const IMATH_INTERNAL_NAMESPACE::M22f test (m1);
  86         assert (test == m1);
  87
  88         IMATH_INTERNAL_NAMESPACE::M22f test2;
  89         assert (test != test2);
  90
  91         IMATH_INTERNAL_NAMESPACE::M22f test3;
  92         test3.makeIdentity();
  93         assert (test2 == test3);
  94
  95         const float a[2][2] = {
  96             { 1.0f, 0.0f },
  97             { 0.0f, 1.0f }
  98         };
  99         testMatrix22ArrayConstructor(a);
 100
 101         m1 = 42;
 102         assert (m1[0][0] == 42 && m1[0][1] == 42 && m1[1][0] == 42 && m1[1][1] == 42);
 103
 104         const float* i1 = test.getValue();
 105         assert(i1[0] == 99.0f);
 106         assert(i1[3] == 101.0f);
 107
 108         float* i2 = m1.getValue();
 109         assert(i2[0] == 42.0f);
 110         assert(i2[1] == 42.0f);
 111         assert(i2[2] == 42.0f);
 112         assert(i2[3] == 42.0f);
 113
 114         IMATH_INTERNAL_NAMESPACE::M22f test4;
 115         test.getValue(test4);
 116         assert (test == test4);
 117
 118         test4.setTheMatrix(test3);
 119         assert(test4 == test3);
 120     }
 121
 122     {
 123         cout << "M22d constructors and equality operators" << endl;
 124
 125         IMATH_INTERNAL_NAMESPACE::M22d m2;
 126         m2[0][0] = 99.0f;
 127         m2[1][1] = 101.0f;
 128
 129         IMATH_INTERNAL_NAMESPACE::M22d test (m2);
 130         assert (test == m2);
 131
 132         IMATH_INTERNAL_NAMESPACE::M22d test2;
 133         assert (test != test2);
 134
 135         IMATH_INTERNAL_NAMESPACE::M22d test3;
 136         test3.makeIdentity();
 137         assert (test2 == test3);
 138
 139         IMATH_INTERNAL_NAMESPACE::M22f test4 (1.0f, 2.0f, 3.0f, 4.0f);
 140
 141         IMATH_INTERNAL_NAMESPACE::M22d test5 = IMATH_INTERNAL_NAMESPACE::M22d (test4);
 142
 143         assert (test5[0][0] == 1.0);
 144         assert (test5[0][1] == 2.0);
 145
 146         assert (test5[1][0] == 3.0);
 147         assert (test5[1][1] == 4.0);
 148
 149         const float a[3][3] = {
 150             { 1.0f, 0.0f, 0.0f },
 151             { 0.0f, 1.0f, 0.0f },
 152             { 0.0f, 0.0f, 1.0f }
 153         };
 154         testMatrix33ArrayConstructor(a);
 155     }
 156
 157     {
 158         cout << "M22f inversion operators" << endl;
 159
 160         IMATH_INTERNAL_NAMESPACE::M22f m1 (3.0f, 3.0f, 5.0f, 5.0f);
 161         IMATH_INTERNAL_NAMESPACE::M22f m2 = m1;
 162         assert(m1.inverse(false) == m1.inverse());
 163         m2.invert(false);
 164         m1.invert();
 165         assert(m1 == m2);
 166
 167         IMATH_INTERNAL_NAMESPACE::M22f m3 (4.0f, 7.0f, 2.0f, 6.0f);
 168         m2 = m3;
 169         assert(m2.inverse(true) == m2.inverse());
 170         m3.invert(true);
 171         m2.invert();
 172         assert(m3 == m2);
 173     }
 174
 175     {
 176         cout << "Imath::M33f shear functions" << endl;
 177
 178         IMATH_INTERNAL_NAMESPACE::M33f m1, m2;
 179         m1.setShear (2.0f);
 180         assert (m1[0][0] == 1.0f && m1[0][1] == 0.0f && m1[0][2] == 0.0f && m1[1][0] == 2.0f &&
 181                 m1[1][1] == 1.0f && m1[1][2] == 0.0f && m1[2][0] == 0.0f && m1[2][1] == 0.0f &&
 182                 m1[2][2] == 1.0f);
 183
 184         m2.setShear (IMATH_INTERNAL_NAMESPACE::V2f (3.0f, 4.0f));
 185         assert (m2[0][0] == 1.0f && m2[0][1] == 4.0f && m2[0][2] == 0.0f && m2[1][0] == 3.0f &&
 186                 m2[1][1] == 1.0f && m2[1][2] == 0.0f && m2[2][0] == 0.0f && m2[2][1] == 0.0f &&
 187                 m2[2][2] == 1.0f);
 188
 189         m1.shear (IMATH_INTERNAL_NAMESPACE::V2f (5.0f, 6.0f));
 190         assert (m1[0][0] == 13.0f && m1[0][1] == 6.0f && m1[0][2] == 0.0f && m1[1][0] == 7.0f &&
 191                 m1[1][1] == 1.0f && m1[1][2] == 0.0f && m1[2][0] == 0.0f && m1[2][1] == 0.0f &&
 192                 m1[2][2] == 1.0f);
 193
 194         m2.shear (7.0f);
 195         assert (m2[0][0] == 1.0f && m2[0][1] == 4.0f && m2[0][2] == 0.0f && m2[1][0] == 10.0f &&
 196                 m2[1][1] == 29.0f && m2[1][2] == 0.0f && m2[2][0] == 0.0f && m2[2][1] == 0.0f &&
 197                 m2[2][2] == 1.0f);
 198
 199         cout << "M33f constructors and equality operators" << endl;
 200
 201         const IMATH_INTERNAL_NAMESPACE::M33f test (m2);
 202         assert (test == m2);
 203
 204         IMATH_INTERNAL_NAMESPACE::M33f test2;
 205         assert (test != test2);
 206
 207         IMATH_INTERNAL_NAMESPACE::M33f test3;
 208         test3.makeIdentity();
 209         assert (test2 == test3);
 210
 211         m1 = 42;
 212         assert (m1[0][0] == 42 && m1[0][1] == 42 && m1[0][2] == 42 &&
 213                 m1[1][0] == 42 && m1[1][1] == 42 && m1[1][2] == 42 &&
 214                 m1[2][0] == 42 && m1[2][1] == 42 && m1[2][2] == 42);
 215
 216         const float* i1 = test.getValue();
 217         assert (i1[0] == 1.0f  && i1[1] == 4.0f  && i1[2] == 0.0f &&
 218                 i1[3] == 10.0f && i1[4] == 29.0f && i1[5] == 0.0f &&
 219                 i1[6] == 0.0f  && i1[7] == 0.0f  && i1[8] == 1.0f);
 220
 221         float* i2 = m1.getValue();
 222         assert(i2[0] == 42.0f);
 223         assert(i2[1] == 42.0f);
 224         assert(i2[2] == 42.0f);
 225         assert(i2[3] == 42.0f);
 226
 227         IMATH_INTERNAL_NAMESPACE::M33f test4;
 228         test.getValue(test4);
 229         assert (test == test4);
 230
 231         test4.setTheMatrix(test3);
 232         assert(test4 == test3);
 233
 234         IMATH_INTERNAL_NAMESPACE::V3f v(2.0f);
 235         IMATH_INTERNAL_NAMESPACE::M33f m(2.0f);
 236         v *= m;
 237         assert (IMATH_INTERNAL_NAMESPACE::equal(v[0], 12.0f, 0.0001f));
 238         assert (IMATH_INTERNAL_NAMESPACE::equal(v[1], 12.0f, 0.0001f));
 239         assert (IMATH_INTERNAL_NAMESPACE::equal(v[2], 12.0f, 0.0001f));
 240     }
 241
 242     {
 243       cout << "M33f inversion operators" << endl;
 244
 245       IMATH_INTERNAL_NAMESPACE::M33f m1 (0.0f, 2.0f, -1.0f, 3.0f, -2.0f, 1.0f, 3.0f, 2.0f, -1.0f);
 246       IMATH_INTERNAL_NAMESPACE::M33f m2 = m1;
 247       assert(m1.inverse(false) == m1.inverse());
 248       m2.invert(false);
 249       m1.invert();
 250       assert(m1 == m2);
 251
 252       IMATH_INTERNAL_NAMESPACE::M33f m3 (1.0f, 0.0f, 5.0f, 2.0f, 1.0f, 6.0f, 3.0f, 4.0f, 0.0f);
 253       m2 = m3;
 254       assert(m3.inverse(true) == m3.inverse());
 255       m3.invert(true);
 256       m2.invert();
 257       assert(m3 == m2);
 258
 259       IMATH_INTERNAL_NAMESPACE::M33f m4 (0.0f, 2.0f, -1.0f, 3.0f, -2.0f, 1.0f, 3.0f, 2.0f, -1.0f);
 260       m2 = m4;
 261       assert(m4.gjInverse(false) == m4.gjInverse());
 262       m2.gjInvert(false);
 263       m4.gjInvert();
 264       assert(m4 == m2);
 265
 266       IMATH_INTERNAL_NAMESPACE::M33f m5 (1.0f, 0.0f, 5.0f, 2.0f, 1.0f, 6.0f, 3.0f, 4.0f, 0.0f);
 267       m2 = m5;
 268       assert(m5.gjInverse(true) == m5.gjInverse());
 269       m5.gjInvert(true);
 270       m2.gjInvert();
 271       assert(m5 == m2);
 272     }
 273
 274     {
 275         cout << "M33d constructors and equality operators" << endl;
 276
 277         IMATH_INTERNAL_NAMESPACE::M33d m2;
 278         m2[0][0] = 99.0f;
 279         m2[1][2] = 101.0f;
 280
 281         IMATH_INTERNAL_NAMESPACE::M33d test (m2);
 282         assert (test == m2);
 283
 284         IMATH_INTERNAL_NAMESPACE::M33d test2;
 285         assert (test != test2);
 286
 287         IMATH_INTERNAL_NAMESPACE::M33d test3;
 288         test3.makeIdentity();
 289         assert (test2 == test3);
 290
 291         IMATH_INTERNAL_NAMESPACE::M33f test4 (1.0f, 2.0f, 3.0f, 4.0f, 5.0f, 6.0f, 7.0f, 8.0f, 9.0f);
 292
 293         IMATH_INTERNAL_NAMESPACE::M33d test5 = IMATH_INTERNAL_NAMESPACE::M33d (test4);
 294
 295         assert (test5[0][0] == 1.0);
 296         assert (test5[0][1] == 2.0);
 297         assert (test5[0][2] == 3.0);
 298
 299         assert (test5[1][0] == 4.0);
 300         assert (test5[1][1] == 5.0);
 301         assert (test5[1][2] == 6.0);
 302
 303         assert (test5[2][0] == 7.0);
 304         assert (test5[2][1] == 8.0);
 305         assert (test5[2][2] == 9.0);
 306     }
 307
 308     {
 309         IMATH_INTERNAL_NAMESPACE::M44f m2;
 310         m2[0][0] = 99.0f;
 311         m2[1][2] = 101.0f;
 312
 313         cout << "M44f constructors and equality operators" << endl;
 314
 315         IMATH_INTERNAL_NAMESPACE::M44f test (m2);
 316         assert (test == m2);
 317
 318         IMATH_INTERNAL_NAMESPACE::M44f test2;
 319         assert (test != test2);
 320
 321         IMATH_INTERNAL_NAMESPACE::M44f test3;
 322         test3.makeIdentity();
 323         assert (test2 == test3);
 324
 325         //
 326         // Test non-equality when a NAN is in the same
 327         // place in two identical matrices
 328         //
 329
 330         test2[0][0] = nanf.f;
 331         test3       = test2;
 332         assert (test2 != test3);
 333     }
 334
 335     {
 336         IMATH_INTERNAL_NAMESPACE::M44d m2;
 337         m2[0][0] = 99.0f;
 338         m2[1][2] = 101.0f;
 339
 340         cout << "M44d constructors and equality operators" << endl;
 341
 342         IMATH_INTERNAL_NAMESPACE::M44d test (m2);
 343         assert (test == m2);
 344
 345         IMATH_INTERNAL_NAMESPACE::M44d test2;
 346         assert (test != test2);
 347
 348         IMATH_INTERNAL_NAMESPACE::M44d test3;
 349         test3.makeIdentity();
 350         assert (test2 == test3);
 351
 352         const float a[4][4] = {
 353             { 1.0f, 0.0f, 0.0f, 0.0f },
 354             { 0.0f, 1.0f, 0.0f, 0.0f },
 355             { 0.0f, 0.0f, 1.0f, 0.0f },
 356             { 0.0f, 0.0f, 0.0f, 1.0f }
 357         };
 358         testMatrix44ArrayConstructor(a);
 359
 360         //
 361         // Test non-equality when a NAN is in the same
 362         // place in two identical matrices
 363         //
 364
 365         test2[0][0] = nand.d;
 366         test3       = test2;
 367         assert (test2 != test3);
 368
 369         IMATH_INTERNAL_NAMESPACE::M44f test4 (1.0f,
 370                                               2.0f,
 371                                               3.0f,
 372                                               4.0f,
 373                                               5.0f,
 374                                               6.0f,
 375                                               7.0f,
 376                                               8.0f,
 377                                               9.0f,
 378                                               10.0f,
 379                                               11.0f,
 380                                               12.0f,
 381                                               13.0f,
 382                                               14.0f,
 383                                               15.0f,
 384                                               16.0f);
 385
 386         IMATH_INTERNAL_NAMESPACE::M44d test5 = IMATH_INTERNAL_NAMESPACE::M44d (test4);
 387
 388         assert (test5[0][0] == 1.0);
 389         assert (test5[0][1] == 2.0);
 390         assert (test5[0][2] == 3.0);
 391         assert (test5[0][3] == 4.0);
 392
 393         assert (test5[1][0] == 5.0);
 394         assert (test5[1][1] == 6.0);
 395         assert (test5[1][2] == 7.0);
 396         assert (test5[1][3] == 8.0);
 397
 398         assert (test5[2][0] == 9.0);
 399         assert (test5[2][1] == 10.0);
 400         assert (test5[2][2] == 11.0);
 401         assert (test5[2][3] == 12.0);
 402
 403         assert (test5[3][0] == 13.0);
 404         assert (test5[3][1] == 14.0);
 405         assert (test5[3][2] == 15.0);
 406         assert (test5[3][3] == 16.0);
 407     }
 408
 409     {
 410         cout << "M44f *= operators" << endl;
 411
 412         IMATH_INTERNAL_NAMESPACE::V3f v(2.0f);
 413         IMATH_INTERNAL_NAMESPACE::M44f m(2.0f);
 414         m.setScale(2.0f);
 415         v *= m;
 416         assert (IMATH_INTERNAL_NAMESPACE::equal(v[0], 4.0f, 0.0001f));
 417         assert (IMATH_INTERNAL_NAMESPACE::equal(v[1], 4.0f, 0.0001f));
 418         assert (IMATH_INTERNAL_NAMESPACE::equal(v[2], 4.0f, 0.0001f));
 419
 420         IMATH_INTERNAL_NAMESPACE::V4f v4f(2.0f);
 421         IMATH_INTERNAL_NAMESPACE::V4f v4f2 = v4f * m;
 422
 423         assert (IMATH_INTERNAL_NAMESPACE::equal(v4f2[0], 4.0f, 0.0001f));
 424         assert (IMATH_INTERNAL_NAMESPACE::equal(v4f2[1], 4.0f, 0.0001f));
 425         assert (IMATH_INTERNAL_NAMESPACE::equal(v4f2[2], 4.0f, 0.0001f));
 426         assert (IMATH_INTERNAL_NAMESPACE::equal(v4f2[3], 2.0f, 0.0001f));
 427
 428         v4f *= m;
 429         assert (v4f == v4f2);
 430
 431         IMATH_INTERNAL_NAMESPACE::M44f a(2.0f);
 432         IMATH_INTERNAL_NAMESPACE::M44f b(3.0f);
 433         IMATH_INTERNAL_NAMESPACE::M44f c;
 434
 435         IMATH_INTERNAL_NAMESPACE::M44f::multiply(a, b, c);
 436         assert (IMATH_INTERNAL_NAMESPACE::equal(c[0][0], 24.0f, 0.0001f));
 437         assert (IMATH_INTERNAL_NAMESPACE::equal(c[1][1], 24.0f, 0.0001f));
 438         assert (IMATH_INTERNAL_NAMESPACE::equal(c[2][2], 24.0f, 0.0001f));
 439         assert (IMATH_INTERNAL_NAMESPACE::equal(c[3][3], 24.0f, 0.0001f));
 440     }
 441
 442     cout << "Matrix << operators" << endl;
 443
 444     {
 445         std::stringstream s;
 446         s << IMATH_INTERNAL_NAMESPACE::identity22f;
 447         const char v[] = "(  1.000000e+00   0.000000e+00\n   0.000000e+00   1.000000e+00)\n";
 448         assert (s.str() == v);
 449     }
 450
 451     {
 452         std::stringstream s;
 453         s << IMATH_INTERNAL_NAMESPACE::identity33f;
 454         const char v[] = "(  1.000000e+00   0.000000e+00   0.000000e+00\n   0.000000e+00   1.000000e+00   0.000000e+00\n   0.000000e+00   0.000000e+00   1.000000e+00)\n";
 455         assert (s.str() == v);
 456     }
 457
 458     {
 459         std::stringstream s;
 460         s << IMATH_INTERNAL_NAMESPACE::identity44f;
 461         const char v[] = "(  1.000000e+00   0.000000e+00   0.000000e+00   0.000000e+00\n   0.000000e+00   1.000000e+00   0.000000e+00   0.000000e+00\n   0.000000e+00   0.000000e+00   1.000000e+00   0.000000e+00\n   0.000000e+00   0.000000e+00   0.000000e+00   1.000000e+00)\n";
 462         assert (s.str() == v);
 463     }
 464
 465     {
 466         cout << "M44f inversion operators" << endl;
 467
 468       IMATH_INTERNAL_NAMESPACE::M44f m1 (1.0f,
 469                                          0.0f,
 470                                          0.0f,
 471                                          0.0f,
 472                                          0.0f,
 473                                          1.0f,
 474                                          0.0f,
 475                                          0.0f,
 476                                          0.0f,
 477                                          0.0f,
 478                                          1.0f,
 479                                          0.0f,
 480                                          0.0f,
 481                                          0.0f,
 482                                          0.0f,
 483                                          0.0f);
 484       IMATH_INTERNAL_NAMESPACE::M44f m2 = m1;
 485       assert(m1.inverse(false) == m1.inverse());
 486       m2.invert(false);
 487       m1.invert();
 488       assert(m1 == m2);
 489
 490       IMATH_INTERNAL_NAMESPACE::M44f m3 (5.0f,
 491                                          6.0f,
 492                                          6.0f,
 493                                          8.0f,
 494                                          2.0f,
 495                                          2.0f,
 496                                          2.0f,
 497                                          8.0f,
 498                                          6.0f,
 499                                          6.0f,
 500                                          2.0f,
 501                                          8.0f,
 502                                          2.0f,
 503                                          3.0f,
 504                                          6.0f,
 505                                          7.0f);
 506       m2 = m3;
 507       assert(m3.inverse(true) == m3.inverse());
 508       m3.invert(true);
 509       m2.invert();
 510       assert(m3 == m2);
 511
 512       IMATH_INTERNAL_NAMESPACE::M44f m4 (1.0f,
 513                                          0.0f,
 514                                          0.0f,
 515                                          0.0f,
 516                                          0.0f,
 517                                          1.0f,
 518                                          0.0f,
 519                                          0.0f,
 520                                          0.0f,
 521                                          0.0f,
 522                                          1.0f,
 523                                          0.0f,
 524                                          0.0f,
 525                                          0.0f,
 526                                          0.0f,
 527                                          0.0f);
 528       m2 = m4;
 529       assert(m4.gjInverse(false) == m4.gjInverse());
 530       m2.gjInvert(false);
 531       m4.gjInvert();
 532       assert(m4 == m2);
 533
 534       IMATH_INTERNAL_NAMESPACE::M44f m5 (5.0f,
 535                                          6.0f,
 536                                          6.0f,
 537                                          8.0f,
 538                                          2.0f,
 539                                          2.0f,
 540                                          2.0f,
 541                                          8.0f,
 542                                          6.0f,
 543                                          6.0f,
 544                                          2.0f,
 545                                          8.0f,
 546                                          2.0f,
 547                                          3.0f,
 548                                          6.0f,
 549                                          7.0f);
 550       m2 = m5;
 551       assert(m5.gjInverse(true) == m5.gjInverse());
 552       m5.gjInvert(true);
 553       m2.gjInvert();
 554       assert(m5 == m2);
 555     }
 556
 557     {
 558         cout << "Converting between M33 and M44" << endl;
 559
 560         IMATH_INTERNAL_NAMESPACE::M44d m1;
 561         m1[0][0] = 99;
 562         IMATH_INTERNAL_NAMESPACE::M44f m2;
 563         m2.setValue (m1);
 564         assert (m2[0][0] == (float) m1[0][0]);
 565         m1[0][0] = 101;
 566         m1.setValue (m2);
 567         assert (m2[0][0] == (float) m1[0][0]);
 568     }
 569
 570     // Matrix minors
 571     {
 572         cout << "3x3 Matrix minors" << endl;
 573
 574         IMATH_INTERNAL_NAMESPACE::M33f a (1, 2, 3, 4, 5, 6, 7, 8, 9);
 575
 576         assert (a.minorOf (0, 0) == a.fastMinor (1, 2, 1, 2));
 577         assert (a.minorOf (0, 1) == a.fastMinor (1, 2, 0, 2));
 578         assert (a.minorOf (0, 2) == a.fastMinor (1, 2, 0, 1));
 579         assert (a.minorOf (1, 0) == a.fastMinor (0, 2, 1, 2));
 580         assert (a.minorOf (1, 1) == a.fastMinor (0, 2, 0, 2));
 581         assert (a.minorOf (1, 2) == a.fastMinor (0, 2, 0, 1));
 582         assert (a.minorOf (2, 0) == a.fastMinor (0, 1, 1, 2));
 583         assert (a.minorOf (2, 1) == a.fastMinor (0, 1, 0, 2));
 584         assert (a.minorOf (2, 2) == a.fastMinor (0, 1, 0, 1));
 585     }
 586     {
 587         IMATH_INTERNAL_NAMESPACE::M33d a (1, 2, 3, 4, 5, 6, 7, 8, 9);
 588
 589         assert (a.minorOf (0, 0) == a.fastMinor (1, 2, 1, 2));
 590         assert (a.minorOf (0, 1) == a.fastMinor (1, 2, 0, 2));
 591         assert (a.minorOf (0, 2) == a.fastMinor (1, 2, 0, 1));
 592         assert (a.minorOf (1, 0) == a.fastMinor (0, 2, 1, 2));
 593         assert (a.minorOf (1, 1) == a.fastMinor (0, 2, 0, 2));
 594         assert (a.minorOf (1, 2) == a.fastMinor (0, 2, 0, 1));
 595         assert (a.minorOf (2, 0) == a.fastMinor (0, 1, 1, 2));
 596         assert (a.minorOf (2, 1) == a.fastMinor (0, 1, 0, 2));
 597         assert (a.minorOf (2, 2) == a.fastMinor (0, 1, 0, 1));
 598     }
 599
 600     // Determinants (by building a random singular value decomposition)
 601
 602     {
 603         cout << "2x2 determinant" << endl;
 604
 605         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 606
 607         IMATH_INTERNAL_NAMESPACE::M22f u;
 608         IMATH_INTERNAL_NAMESPACE::M22f v;
 609         IMATH_INTERNAL_NAMESPACE::M22f s;
 610
 611         u.setRotation (random.nextf());
 612         v.setRotation (random.nextf());
 613         s[0][0] = random.nextf();
 614         s[1][1] = random.nextf();
 615
 616         IMATH_INTERNAL_NAMESPACE::M22f c = u * s * v.transpose();
 617         assert (fabsf (c.determinant() - s[0][0] * s[1][1]) <= u.baseTypeEpsilon());
 618     }
 619     {
 620         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 621
 622         IMATH_INTERNAL_NAMESPACE::M22d u;
 623         IMATH_INTERNAL_NAMESPACE::M22d v;
 624         IMATH_INTERNAL_NAMESPACE::M22d s;
 625
 626         u.setRotation ((double) random.nextf());
 627         v.setRotation ((double) random.nextf());
 628         s[0][0] = (double) random.nextf();
 629         s[1][1] = (double) random.nextf();
 630
 631         IMATH_INTERNAL_NAMESPACE::M22d c = u * s * v.transpose();
 632         assert (fabs (c.determinant() - s[0][0] * s[1][1]) <= u.baseTypeEpsilon());
 633     }
 634
 635     {
 636         cout << "3x3 determinant" << endl;
 637
 638         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 639
 640         IMATH_INTERNAL_NAMESPACE::M33f u;
 641         IMATH_INTERNAL_NAMESPACE::M33f v;
 642         IMATH_INTERNAL_NAMESPACE::M33f s;
 643
 644         u.setRotation (random.nextf());
 645         v.setRotation (random.nextf());
 646         s[0][0] = random.nextf();
 647         s[1][1] = random.nextf();
 648         s[2][2] = random.nextf();
 649
 650         IMATH_INTERNAL_NAMESPACE::M33f c = u * s * v.transpose();
 651         assert (fabsf (c.determinant() - s[0][0] * s[1][1] * s[2][2]) <= u.baseTypeEpsilon());
 652     }
 653     {
 654         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 655
 656         IMATH_INTERNAL_NAMESPACE::M33d u;
 657         IMATH_INTERNAL_NAMESPACE::M33d v;
 658         IMATH_INTERNAL_NAMESPACE::M33d s;
 659
 660         u.setRotation ((double) random.nextf());
 661         v.setRotation ((double) random.nextf());
 662         s[0][0] = (double) random.nextf();
 663         s[1][1] = (double) random.nextf();
 664         s[2][2] = (double) random.nextf();
 665
 666         IMATH_INTERNAL_NAMESPACE::M33d c = u * s * v.transpose();
 667         assert (fabs (c.determinant() - s[0][0] * s[1][1] * s[2][2]) <= u.baseTypeEpsilon());
 668     }
 669
 670     // Outer product of two 3D vectors
 671     {
 672         cout << "Outer product of two 3D vectors" << endl;
 673
 674         IMATH_INTERNAL_NAMESPACE::V3f a (1, 2, 3);
 675         IMATH_INTERNAL_NAMESPACE::V3f b (4, 5, 6);
 676         IMATH_INTERNAL_NAMESPACE::M33f p = IMATH_INTERNAL_NAMESPACE::outerProduct (a, b);
 677
 678         for (int i = 0; i < 3; i++)
 679         {
 680             for (int j = 0; j < 3; j++)
 681             {
 682                 assert (p[i][j] == a[i] * b[j]);
 683             }
 684         }
 685     }
 686     {
 687         IMATH_INTERNAL_NAMESPACE::V3d a (1, 2, 3);
 688         IMATH_INTERNAL_NAMESPACE::V3d b (4, 5, 6);
 689         IMATH_INTERNAL_NAMESPACE::M33d p = IMATH_INTERNAL_NAMESPACE::outerProduct (a, b);
 690
 691         for (int i = 0; i < 3; i++)
 692         {
 693             for (int j = 0; j < 3; j++)
 694             {
 695                 assert (p[i][j] == a[i] * b[j]);
 696             }
 697         }
 698     }
 699
 700     // Determinants (by building a random singular value decomposition)
 701     {
 702         cout << "4x4 determinants" << endl;
 703
 704         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 705
 706         IMATH_INTERNAL_NAMESPACE::M44f u = IMATH_INTERNAL_NAMESPACE::rotationMatrix (
 707             IMATH_INTERNAL_NAMESPACE::V3f (random.nextf(), random.nextf(), random.nextf())
 708                 .normalize(),
 709             IMATH_INTERNAL_NAMESPACE::V3f (random.nextf(), random.nextf(), random.nextf())
 710                 .normalize());
 711         IMATH_INTERNAL_NAMESPACE::M44f v = IMATH_INTERNAL_NAMESPACE::rotationMatrix (
 712             IMATH_INTERNAL_NAMESPACE::V3f (random.nextf(), random.nextf(), random.nextf())
 713                 .normalize(),
 714             IMATH_INTERNAL_NAMESPACE::V3f (random.nextf(), random.nextf(), random.nextf())
 715                 .normalize());
 716         IMATH_INTERNAL_NAMESPACE::M44f s;
 717
 718         s[0][0] = random.nextf();
 719         s[1][1] = random.nextf();
 720         s[2][2] = random.nextf();
 721         s[3][3] = random.nextf();
 722
 723         IMATH_INTERNAL_NAMESPACE::M44f c = u * s * v.transpose();
 724         assert (fabsf (c.determinant() - s[0][0] * s[1][1] * s[2][2] * s[3][3]) <=
 725                 u.baseTypeEpsilon());
 726     }
 727     {
 728         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 729
 730         IMATH_INTERNAL_NAMESPACE::M44d u = IMATH_INTERNAL_NAMESPACE::rotationMatrix (
 731             IMATH_INTERNAL_NAMESPACE::V3d (random.nextf(), random.nextf(), random.nextf())
 732                 .normalize(),
 733             IMATH_INTERNAL_NAMESPACE::V3d (random.nextf(), random.nextf(), random.nextf())
 734                 .normalize());
 735         IMATH_INTERNAL_NAMESPACE::M44d v = IMATH_INTERNAL_NAMESPACE::rotationMatrix (
 736             IMATH_INTERNAL_NAMESPACE::V3d (random.nextf(), random.nextf(), random.nextf())
 737                 .normalize(),
 738             IMATH_INTERNAL_NAMESPACE::V3d (random.nextf(), random.nextf(), random.nextf())
 739                 .normalize());
 740         IMATH_INTERNAL_NAMESPACE::M44d s;
 741
 742         s[0][0] = random.nextf();
 743         s[1][1] = random.nextf();
 744         s[2][2] = random.nextf();
 745         s[3][3] = random.nextf();
 746
 747         IMATH_INTERNAL_NAMESPACE::M44d c = u * s * v.transpose();
 748         assert (fabs (c.determinant() - s[0][0] * s[1][1] * s[2][2] * s[3][3]) <=
 749                 u.baseTypeEpsilon());
 750     }
 751
 752     // Trace
 753     {
 754         cout << "2x2 trace" << endl;
 755
 756         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 757
 758         IMATH_INTERNAL_NAMESPACE::M22f u;
 759         float                          trace = 0;
 760         for (int i = 0; i < 2; i++)
 761         {
 762             for (int j = 0; j < 2; j++)
 763             {
 764                 const float randomNum = random.nextf ();
 765                 u[i][j]               = randomNum;
 766                 if (i == j) { trace += randomNum; }
 767             }
 768         }
 769
 770         assert (fabsf (u.trace () - trace) <= u.baseTypeEpsilon ());
 771     }
 772     {
 773         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 774
 775         IMATH_INTERNAL_NAMESPACE::M22d u;
 776         double                         trace = 0;
 777         for (int i = 0; i < 2; i++)
 778         {
 779             for (int j = 0; j < 2; j++)
 780             {
 781                 const double randomNum = random.nextf ();
 782                 u[i][j]                = randomNum;
 783                 if (i == j) { trace += randomNum; }
 784             }
 785         }
 786
 787         assert (fabs (u.trace () - trace) <= u.baseTypeEpsilon ());
 788     }
 789     {
 790         cout << "3x3 trace" << endl;
 791
 792         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 793
 794         IMATH_INTERNAL_NAMESPACE::M33f u;
 795         float                          trace = 0;
 796         for (int i = 0; i < 3; i++)
 797         {
 798             for (int j = 0; j < 3; j++)
 799             {
 800                 const float randomNum = random.nextf ();
 801                 u[i][j]               = randomNum;
 802                 if (i == j) { trace += randomNum; }
 803             }
 804         }
 805
 806         assert (fabsf (u.trace () - trace) <= u.baseTypeEpsilon ());
 807     }
 808     {
 809         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 810
 811         IMATH_INTERNAL_NAMESPACE::M33d u;
 812         double                         trace = 0;
 813         for (int i = 0; i < 3; i++)
 814         {
 815             for (int j = 0; j < 3; j++)
 816             {
 817                 const double randomNum = random.nextf ();
 818                 u[i][j]                = randomNum;
 819                 if (i == j) { trace += randomNum; }
 820             }
 821         }
 822
 823         assert (fabs (u.trace () - trace) <= u.baseTypeEpsilon ());
 824     }
 825     {
 826         cout << "4x4 trace" << endl;
 827
 828         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 829
 830         IMATH_INTERNAL_NAMESPACE::M44f u;
 831         float                          trace = 0;
 832         for (int i = 0; i < 4; i++)
 833         {
 834             for (int j = 0; j < 4; j++)
 835             {
 836                 const float randomNum = random.nextf ();
 837                 u[i][j]               = randomNum;
 838                 if (i == j) { trace += randomNum; }
 839             }
 840         }
 841
 842         assert (fabsf (u.trace () - trace) <= u.baseTypeEpsilon ());
 843     }
 844     {
 845         IMATH_INTERNAL_NAMESPACE::Rand32 random;
 846
 847         IMATH_INTERNAL_NAMESPACE::M44d u;
 848         double                         trace = 0;
 849         for (int i = 0; i < 4; i++)
 850         {
 851             for (int j = 0; j < 4; j++)
 852             {
 853                 const double randomNum = random.nextf ();
 854                 u[i][j]                = randomNum;
 855                 if (i == j) { trace += randomNum; }
 856             }
 857         }
 858
 859         assert (fabs (u.trace () - trace) <= u.baseTypeEpsilon ());
 860     }
 861
 862     // Matrix minors
 863     {
 864         cout << "4x4 matrix minors" << endl;
 865
 866         IMATH_INTERNAL_NAMESPACE::M44d a (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16);
 867
 868         assert (a.minorOf (0, 0) == a.fastMinor (1, 2, 3, 1, 2, 3));
 869         assert (a.minorOf (0, 1) == a.fastMinor (1, 2, 3, 0, 2, 3));
 870         assert (a.minorOf (0, 2) == a.fastMinor (1, 2, 3, 0, 1, 3));
 871         assert (a.minorOf (0, 3) == a.fastMinor (1, 2, 3, 0, 1, 2));
 872         assert (a.minorOf (1, 0) == a.fastMinor (0, 2, 3, 1, 2, 3));
 873         assert (a.minorOf (1, 1) == a.fastMinor (0, 2, 3, 0, 2, 3));
 874         assert (a.minorOf (1, 2) == a.fastMinor (0, 2, 3, 0, 1, 3));
 875         assert (a.minorOf (1, 3) == a.fastMinor (0, 2, 3, 0, 1, 2));
 876         assert (a.minorOf (2, 0) == a.fastMinor (0, 1, 3, 1, 2, 3));
 877         assert (a.minorOf (2, 1) == a.fastMinor (0, 1, 3, 0, 2, 3));
 878         assert (a.minorOf (2, 2) == a.fastMinor (0, 1, 3, 0, 1, 3));
 879         assert (a.minorOf (2, 3) == a.fastMinor (0, 1, 3, 0, 1, 2));
 880         assert (a.minorOf (3, 0) == a.fastMinor (0, 1, 2, 1, 2, 3));
 881         assert (a.minorOf (3, 1) == a.fastMinor (0, 1, 2, 0, 2, 3));
 882         assert (a.minorOf (3, 2) == a.fastMinor (0, 1, 2, 0, 1, 3));
 883         assert (a.minorOf (3, 3) == a.fastMinor (0, 1, 2, 0, 1, 2));
 884     }
 885     {
 886         IMATH_INTERNAL_NAMESPACE::M44f a (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16);
 887
 888         assert (a.minorOf (0, 0) == a.fastMinor (1, 2, 3, 1, 2, 3));
 889         assert (a.minorOf (0, 1) == a.fastMinor (1, 2, 3, 0, 2, 3));
 890         assert (a.minorOf (0, 2) == a.fastMinor (1, 2, 3, 0, 1, 3));
 891         assert (a.minorOf (0, 3) == a.fastMinor (1, 2, 3, 0, 1, 2));
 892         assert (a.minorOf (1, 0) == a.fastMinor (0, 2, 3, 1, 2, 3));
 893         assert (a.minorOf (1, 1) == a.fastMinor (0, 2, 3, 0, 2, 3));
 894         assert (a.minorOf (1, 2) == a.fastMinor (0, 2, 3, 0, 1, 3));
 895         assert (a.minorOf (1, 3) == a.fastMinor (0, 2, 3, 0, 1, 2));
 896         assert (a.minorOf (2, 0) == a.fastMinor (0, 1, 3, 1, 2, 3));
 897         assert (a.minorOf (2, 1) == a.fastMinor (0, 1, 3, 0, 2, 3));
 898         assert (a.minorOf (2, 2) == a.fastMinor (0, 1, 3, 0, 1, 3));
 899         assert (a.minorOf (2, 3) == a.fastMinor (0, 1, 3, 0, 1, 2));
 900         assert (a.minorOf (3, 0) == a.fastMinor (0, 1, 2, 1, 2, 3));
 901         assert (a.minorOf (3, 1) == a.fastMinor (0, 1, 2, 0, 2, 3));
 902         assert (a.minorOf (3, 2) == a.fastMinor (0, 1, 2, 0, 1, 3));
 903         assert (a.minorOf (3, 3) == a.fastMinor (0, 1, 2, 0, 1, 2));
 904     }
 905
 906     // VC 2005 64 bits compiler has a bug with __restrict keword.
 907     // Pointers with __restrict should not alias the same symbol.
 908     // But, with optimization on, VC removes intermediate temp variable
 909     // and ignores __restrict.
 910     {
 911         cout << "M44 multiplicaftion test" << endl;
 912         IMATH_INTERNAL_NAMESPACE::M44f M (1.0f,
 913                                           2.0f,
 914                                           3.0f,
 915                                           4.0f,
 916                                           5.0f,
 917                                           6.0f,
 918                                           7.0f,
 919                                           8.0f,
 920                                           9.0f,
 921                                           10.0f,
 922                                           11.0f,
 923                                           12.0f,
 924                                           13.0f,
 925                                           14.0f,
 926                                           15.0f,
 927                                           16.0f);
 928
 929         IMATH_INTERNAL_NAMESPACE::M44f N;
 930         N.makeIdentity();
 931
 932         // N should be equal to M
 933         // This typical test fails
 934         // when __restrict is used for pointers in "multiply" function.
 935         N = N * M;
 936
 937         assert (N == M);
 938
 939         if (N != M)
 940         {
 941             cout << "M44 multiplication test has failed, error." << endl
 942                  << "M" << endl
 943                  << M << endl
 944                  << "N" << endl
 945                  << N << endl;
 946         }
 947     }
 948
 949     cout << "ok\n" << endl;
 950 }