testSMatrixSVector.cc

Go to the documentation of this file.

#include <iostream>

#include <utl/SVector.h>
#include <utl/SMatrix.h>

#include <tst/Verify.h>
#include <cppunit/extensions/HelperMacros.h>

using namespace std;
using namespace utl;
using namespace tst;

#define ASSERT_CLOSE(x, y) CPPUNIT_ASSERT(Verify<CloseTo>(x, y))
#define ASSERT_CLOSE_EPS(x, y, eps) CPPUNIT_ASSERT(Verify<CloseTo>(x, y, eps))
#define ASSERT_EQUAL(x, y) CPPUNIT_ASSERT(Verify<Equal>(x, y))


class Foo {
public:
  Foo(const int off = 0) : fOffset(off) { }

  int operator[](const unsigned int i) const { return 2*int(i) + fOffset; }

private:
  const int fOffset;
};


class Bar {
public:
  Foo operator[](const unsigned int i) const { return Foo(i); }
};


class TestSMatrixSVector : public CppUnit::TestFixture {

  CPPUNIT_TEST_SUITE(TestSMatrixSVector);
  CPPUNIT_TEST(TestSVectorCtor);
  CPPUNIT_TEST(TestSVectorOp);
  CPPUNIT_TEST(TestSVectorAngle);
  CPPUNIT_TEST(TestSVectorListAssignment);
  CPPUNIT_TEST(TestSVectorIterators);
  CPPUNIT_TEST(TestSMatrixCtor);
  CPPUNIT_TEST(TestSMatrixOp);
  //CPPUNIT_TEST(TestSMatrixSVectorOp);
  CPPUNIT_TEST(TestSMatrixListAssignment);
  CPPUNIT_TEST(TestSVectorNontrivialTypes);
  CPPUNIT_TEST_SUITE_END();

public:
  void
  TestSVectorCtor()
  {
    // default ctor
    SVector<10> v1;

    ASSERT_EQUAL((unsigned int)(SVector<10>::Size), 10U);
    ASSERT_EQUAL(v1.GetSize(), std::size_t(10));

    // clear
    v1.Clear();

    for (unsigned int i = 0; i < v1.GetSize(); ++i)
      ASSERT_EQUAL(v1[i], 0.);

    ASSERT_EQUAL(v1.Get<0>(), 0.);
    ASSERT_EQUAL(v1.Get<1>(), 0.);
    ASSERT_EQUAL(v1.Get<2>(), 0.);
    ASSERT_EQUAL(v1.Get<3>(), 0.);
    ASSERT_EQUAL(v1.Get<4>(), 0.);

    const int v1value = 1;
    v1.Clear(v1value);

    for (unsigned int i = 0; i < v1.GetSize(); ++i)
      ASSERT_EQUAL(v1[i], double(v1value));

    // scalar initializer
    const int v2value = 13;
    SVector<100, int> v2(v2value);

    for (unsigned int i = 0; i < v2.GetSize(); ++i)
      ASSERT_EQUAL(v2[i], 13);

    ASSERT_EQUAL(v2.Get<0>(), 13);
    ASSERT_EQUAL(v2.Get<1>(), 13);
    ASSERT_EQUAL(v2.Get<10>(), 13);
    ASSERT_EQUAL(v2.Get<50>(), 13);
    ASSERT_EQUAL(v2.Get<99>(), 13);

    // functor initializer
    Foo foo;
    SVector<5> v3(foo);

    for (unsigned int i = 0; i < v3.GetSize(); ++i)
      ASSERT_EQUAL(v3[i], double(foo[i]));

    // c-array initializer
    const int a[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 };

    SVector<5, int> v4(a);

    for (unsigned int i = 0; i < v4.GetSize(); ++i)
      ASSERT_EQUAL(v4[i], int(i+1));

    SVector<7> v41(a);

    for (unsigned int i = 0; i < v41.GetSize(); ++i)
      ASSERT_EQUAL(v41[i], i+1.);

    // pointer initializer
    const int* const pa = &(a[0]);

    SVector<5, int> v5(pa);

    for (unsigned int i = 0; i < v5.GetSize(); ++i)
      ASSERT_EQUAL(v5[i], int(i+1));

    // copy ctor
    SVector<5, int> v6(v5);

    for (unsigned int i = 0; i < v6.GetSize(); ++i)
      ASSERT_EQUAL(v6[i], v5[i]);

    // non-identical copy
    SVector<100> v7(v2);

    for (unsigned int i = 0; i < v7.GetSize(); ++i)
      ASSERT_EQUAL(v7[i], double(v2[i]));
  }

  void
  TestSVectorOp()
  {
    SVector<10> v1(1.);

    v1 *= 2;

    for (unsigned int i = 0; i < v1.GetSize(); ++i)
      ASSERT_EQUAL(v1[i], 2.);

    SVector<10> v2(0.);

    v2 += v1;

    for (unsigned int i = 0; i < v2.GetSize(); ++i)
      ASSERT_EQUAL(v2[i], 2.);

    v2 += v1;
    v2 *= 0.5;
    v2 -= v1;

    for (unsigned int i = 0; i < v2.GetSize(); ++i)
      ASSERT_EQUAL(v2[i], 0.);

    SVector<10> v3;
    SVector<10> v4;

    for (unsigned int i = 0; i < v3.GetSize(); ++i) {
      v3[i] = i;
      v4[i] = -2 * int(i);
    }

    const double dot = v3 * v4;

    ASSERT_EQUAL(dot, -570.);

    const SVector<10> v5 = v3 + v4;

    for (unsigned int i = 0; i < v5.GetSize(); ++i)
      ASSERT_EQUAL(v5[i], -double(i));

    const SVector<10> v6 = v3 - v4;

    for (unsigned int i = 0; i < v5.GetSize(); ++i)
      ASSERT_EQUAL(v6[i], 3.*i);
  }

  void
  TestSVectorAngle()
  {
    const double eps = 1e-8;
    SVector<3> v1;
    v1 = 1., 1., 1.;
    SVector<3> v2;
    v2 = 1., 1., 1.+eps;
    ASSERT_CLOSE(CosAngle(v1, v2), 1.);
    ASSERT_CLOSE(Angle(v1, v2), 0.);
    SVector<3> v3;
    v3 = -1., -1., -1.+eps;
    ASSERT_CLOSE(CosAngle(v1, v3), -1.);
    ASSERT_CLOSE(Angle(v1, v3), M_PI);
  }

  void
  TestSVectorListAssignment()
  {
    {
      SVector<5> v;
      v = 1., 2, 3, 4, 5;
      double count = 0;
      for (unsigned int i = 0; i < v.GetSize(); ++i)
        ASSERT_EQUAL(v[i], ++count);
    }

    {
      SVector<10, int> v(0);
      v = 1, 2, 3;

      for (unsigned int i = 0; i < 3; ++i)
        ASSERT_EQUAL(v[i], int(i+1));
      for (unsigned int i = 3; i < v.GetSize(); ++i)
        ASSERT_EQUAL(v[i], 0);
    }
    {
      SVector<10> v(0.);
      v = 1., 2, 3;

      for (unsigned int i = 0; i < 3; ++i)
        ASSERT_EQUAL(v[i], double(i+1));
      for (unsigned int i = 3; i < v.GetSize(); ++i)
        ASSERT_EQUAL(v[i], 0.);
    }
    {
      SVector<10> v(0.);
      v = 1., 2, 3;

      SVector<10, SVector<10> > m(SVector<10>(0.));
      m = v, v, v;

      for (unsigned int i = 0; i < 3; ++i)
        for (unsigned int j = 0; j < m[0].GetSize(); ++j)
          ASSERT_EQUAL(m[i][j], v[j]);
      for (unsigned int i = 3; i < v.GetSize(); ++i)
        for (unsigned int j = 0; j < m[0].GetSize(); ++j)
        ASSERT_EQUAL(m[i][j], 0.);
    }
  }

  void
  TestSVectorIterators()
  {
    const int n = 10;
    typedef SVector<n, int> SVNI;
    SVNI v;
    for (int i = 0; i < n; ++i)
      v[i] = i;
    int i = 0;
    for (SVNI::Iterator it = v.Begin(); it != v.End(); ++it) {
      ASSERT_EQUAL(*it, i);
      ++i;
    }
    i = 0;
    for (SVNI::ConstIterator it = v.Begin(); it != v.End(); ++it) {
      ASSERT_EQUAL(*it, i);
      ++i;
    }
  }

  void
  TestSMatrixCtor()
  {
    {
      // default ctor
      SMatrix<10, 100> m;

      ASSERT_EQUAL((SMatrix<10, 100>::NRows), 10U);
      ASSERT_EQUAL((SMatrix<10, 100>::NColumns), 100U);
      ASSERT_EQUAL(m.GetNRows(), 10U);
      ASSERT_EQUAL(m.GetNColumns(), 100U);

      // clear
      m.Clear();

      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], 0.);

      ASSERT_EQUAL((m.Get<0, 0>()), 0.);
      ASSERT_EQUAL((m.Get<9, 0>()), 0.);
      ASSERT_EQUAL((m.Get<0, 99>()), 0.);
      ASSERT_EQUAL((m.Get<9, 99>()), 0.);
      ASSERT_EQUAL((m.Get<5, 50>()), 0.);

      const int value = 1;
      m.Clear(value);

      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], double(value));
    }

    {
      // scalar initializer
      const int value = 13;
      SMatrix<100, 5, int> m(value);

      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], 13);

      ASSERT_EQUAL((m.Get<0, 0>()), 13);
      ASSERT_EQUAL((m.Get<99, 0>()), 13);
      ASSERT_EQUAL((m.Get<0, 4>()), 13);
      ASSERT_EQUAL((m.Get<99, 4>()), 13);
      ASSERT_EQUAL((m.Get<50, 2>()), 13);

      // non-identical copy
      SMatrix<100, 5> md(m);

      for (unsigned int i = 0; i < md.GetNRows(); ++i)
        for (unsigned int j = 0; j < md.GetNColumns(); ++j)
          ASSERT_EQUAL(md[i][j], double(m[i][j]));
    }

    {
      // functor initializer
      Bar bar;
      SMatrix<5, 5> m(bar);

      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], double(bar[i][j]));
    }

    {
      // c-array initializer
      const int a[][3] = {
        { 1, 2, 3 },
        { 4, 5, 6 },
        { 7, 8, 9 },
        { 10, 11, 12 }
      };

      SMatrix<4, 3, int> m(a);

      {
        int count = 0;
        for (unsigned int i = 0; i < m.GetNRows(); ++i)
          for (unsigned int j = 0; j < m.GetNColumns(); ++j)
            ASSERT_EQUAL(m[i][j], ++count);
      }

      // nonidentical copy
      SMatrix<4, 3> md(a);

      {
        double count = 0;
        for (unsigned int i = 0; i < md.GetNRows(); ++i)
          for (unsigned int j = 0; j < md.GetNColumns(); ++j)
            ASSERT_EQUAL(md[i][j], ++count);
      }

      // copy ctor
      SMatrix<4, 3, int> m2(m);

      for (unsigned int i = 0; i < m2.GetNRows(); ++i)
        for (unsigned int j = 0; j < m2.GetNColumns(); ++j)
          ASSERT_EQUAL(m2[i][j], m[i][j]);
    }

    {
      int** ma = new int*[10];
      int count = 0;
      for (unsigned int i = 0; i < 10; ++i) {
        ma[i] = new int[100];
        for (unsigned int j = 0; j < 100; ++j)
          ma[i][j] = ++count;
      }

      // ** ctor
      SMatrix<10, 100> m(ma);

      for (unsigned int i = 0; i < 10; ++i)
        for (unsigned int j = 0; j < 100; j++)
          ASSERT_EQUAL(m[i][j], double(ma[i][j]));

      for (unsigned int i = 0; i < 10; ++i)
        delete[] ma[i];
      delete[] ma;
    }
  }

  void
  TestSMatrixOp()
  {
    SMatrix<10, 5> m1(1.);

    m1 *= 2;

    for (unsigned int i = 0; i < m1.GetNRows(); ++i)
      for (unsigned int j = 0; j < m1.GetNColumns(); ++j)
        ASSERT_EQUAL(m1[i][j], 2.);

    SMatrix<10, 5> m2(0.);

    m2 += m1;

    for (unsigned int i = 0; i < m2.GetNRows(); ++i)
      for (unsigned int j = 0; j < m2.GetNColumns(); ++j)
        ASSERT_EQUAL(m2[i][j], 2.);

    m2 += m1;
    m2 *= 0.5;
    m2 -= m1;

    for (unsigned int i = 0; i < m2.GetNRows(); ++i)
      for (unsigned int j = 0; j < m2.GetNColumns(); ++j)
        ASSERT_EQUAL(m2[i][j], 0.);

    {
      const unsigned int n = 100;
      const unsigned int m = 50;
      SMatrix<n, m> m3;
      SMatrix<m, n> m4;

      for (unsigned int i = 0; i < m3.GetNRows(); ++i)
        for (unsigned int j = 0; j < m3.GetNColumns(); ++j) {
          m3[i][j] = i + j;
          m4[j][i] = -2 * int(i) + int(j);
      }

      const SMatrix<10, 10> prod(m3 * m4);

      {
        const int m = m3.GetNColumns();
        ASSERT_EQUAL(m, int(m4.GetNRows()));
        for (unsigned int i = 0; i < prod.GetNRows(); ++i)
          for (unsigned int j = 0; j < prod.GetNColumns(); ++j) {
            const int ii = i;
            const int jj = j;
            ASSERT_EQUAL(
              prod[i][j],
              double(-2*m*ii*jj+(ii-2*jj)*(m-1)*m/2+(m-1)*m*(2*m-1)/6)
            );
          }
      }
    }
  }

  /*void
  TestSMatrixSVectorOp()
  {
    const unsigned int n = 10;
    const unsigned int m = 5;

    SMatrix<n, m> a;

    for (unsigned int i = 0; i < a.GetNRows(); ++i)
      for (unsigned int j = 0; j < a.GetNColumns(); ++j)
        a[i][j] = int(i) - 3*int(j);

    SVector<m> x;
    SVector<m, int> xi;

    for (unsigned int i = 0; i < x.GetSize(); ++i) {
      x[i] = i;
      xi[i] = i;
    }

    const SVector<n> y = a * x;
    const SVector<n> yi = a * xi;

    {
      const int mm = m;
      for (unsigned int i = 0; i < x.GetSize(); ++i) {
        const int ii = i;
        ASSERT_EQUAL(y[i],
                     double(ii*(mm-1)*mm/2-(mm-1)*mm*(2*mm-1)/2));
        ASSERT_EQUAL(yi[i],
                     double(ii*(mm-1)*mm/2-(mm-1)*mm*(2*mm-1)/2));
      }
    }
  }*/

  void
  TestSMatrixListAssignment()
  {
    {
      SMatrix<3, 4> m;
      m =
        1., 2., 3., 4.,
        5., 6., 7., 8.,
        9., 10, 11, 12;
      double count = 0;
      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], ++count);
    }

    {
      SMatrix<21, 13, int> m;
      m =
        1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,
        14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26,
        27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39,
        40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52,
        53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65,
        66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78,
        79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91,
        92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104,
        105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117,
        118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130,
        131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143,
        144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156,
        157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169,
        170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 182,
        183, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195,
        196, 197, 198, 199, 200, 201, 202, 203, 204, 205, 206, 207, 208,
        209, 210, 211, 212, 213, 214, 215, 216, 217, 218, 219, 220, 221,
        222, 223, 224, 225, 226, 227, 228, 229, 230, 231, 232, 233, 234,
        235, 236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247,
        248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260,
        261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273;
      int count = 0;
      for (unsigned int i = 0; i < m.GetNRows(); ++i)
        for (unsigned int j = 0; j < m.GetNColumns(); ++j)
          ASSERT_EQUAL(m[i][j], ++count);
    }
  }

  void
  TestSVectorNontrivialTypes()
  {
    {
      // dot product of two vectors of vectors is a scalar
      // (recursive dot product)

      SVector<100, SVector<10> > v1(SVector<10>(1.));
      SVector<100, SVector<10> > v2(SVector<10>(2.));

      const double partialDot = v1[0] * v2[0];

      ASSERT_EQUAL(partialDot, 2.*v1[0].GetSize());

      const double dot = v1 * v2;

      ASSERT_EQUAL(dot, 2.*v1.GetSize()*v1[0].GetSize());
    }

    {
      // dot product of two vectors of matrices is a matrix

      SVector<10, SMatrix<20, 30> > v1(SMatrix<20, 30>(1.));
      SVector<10, SMatrix<30, 20> > v2(SMatrix<30, 20>(2.));

      const SMatrix<20, 20> partialDot = v1[0] * v2[0];

      const SMatrix<20, 20> partialRes(2.*v1[0].GetNColumns());

      ASSERT_EQUAL(partialDot, partialRes);

      const SMatrix<20, 20> dot = v1 * v2;

      const SMatrix<20, 20> dotRes(2.*v1[0].GetNColumns()*v1.GetSize());

      ASSERT_EQUAL(dot, dotRes);
    }

    {
      // dot product of two vectors of matrices of different type
      // is a matrix

      SVector<10, SMatrix<20, 30, int> > v1(SMatrix<20, 30, int>(1));
      SVector<10, SMatrix<30, 20, double> > v2(SMatrix<30, 20, double>(2.));

      const SMatrix<20, 20, double> partialDot = v1[0] * v2[0];

      const SMatrix<20, 20, double> partialRes(2.*v1[0].GetNColumns());

      ASSERT_EQUAL(partialDot, partialRes);

      const SMatrix<20, 20, double> dot = v1 * v2;

      const SMatrix<20, 20, double> dotRes(2.*v1[0].GetNColumns()*v1.GetSize());

      ASSERT_EQUAL(dot, dotRes);
    }
  }

};


CPPUNIT_TEST_SUITE_REGISTRATION(TestSMatrixSVector);