testSparseMatrix.cc

Go to the documentation of this file.

#include <iostream>
#include <numeric>

#include <utl/SparseMatrix.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))
#define ASSERT_NOT_EQUAL(x, y) CPPUNIT_ASSERT(Verify<Not<Equal> >(x, y))
#define ASSERT_ITERATOR_EQUAL(it, e, r, c, v) \
  CPPUNIT_ASSERT(it != e); \
  ASSERT_EQUAL(it.GetRow(), r); \
  ASSERT_EQUAL(it.GetColumn(), c); \
  ASSERT_EQUAL(it(), v)


class TestSparseMatrix : public CppUnit::TestFixture {

  CPPUNIT_TEST_SUITE(TestSparseMatrix);
  CPPUNIT_TEST(TestAccess);
  CPPUNIT_TEST(TestIterators);
  CPPUNIT_TEST(TestOperators);
  CPPUNIT_TEST(TestMultiplication);
  CPPUNIT_TEST(TestNonSquareMultiplication);
  CPPUNIT_TEST(TestNonSquareMultiplication2);
  CPPUNIT_TEST(TestTranspose);
  CPPUNIT_TEST_SUITE_END();

public:
  void
  TestAccess()
  {
    // assume 10 x 10
    const int rows = 5;
    const int cols = 5;
    SparseMatrix<> m;

    m(0, 0) = 1;
    m(2, 4) = 2;
    m(4, 4) = 3;
    m(1, 3) = 4;
    m(3, 1) = 5;

    stringstream ss;

    OutputSparse(ss, m, rows, cols);
    const string test1 =
      "1 . . . .\n"
      ". . . 4 .\n"
      ". . . . 2\n"
      ". 5 . . .\n"
      ". . . . 3\n";
    ASSERT_EQUAL(ss.str(), test1);

    for (int i = 0; i < rows; ++i)
      for (int j = 0; j < cols; ++j)
        m[i][j];

    ss.str("");
    OutputSparse(ss, m, rows, cols);
    ASSERT_EQUAL(ss.str(), test1);

    for (int i = 0; i < rows; ++i)
      for (int j = 0; j < cols; ++j)
        m(i, j);

    ss.str("");
    OutputSparse(ss, m, rows, cols);
    const string test2 =
      "1 0 0 0 0\n"
      "0 0 0 4 0\n"
      "0 0 0 0 2\n"
      "0 5 0 0 0\n"
      "0 0 0 0 3\n";
    ASSERT_EQUAL(ss.str(), test2);

    m.Reclaim();

    ss.str("");
    OutputSparse(ss, m, rows, cols);
    ASSERT_EQUAL(ss.str(), test1);
  }

  void
  TestIterators()
  {
    // assume 10 x 20
    SparseMatrix<int> m;

    m(1, 1) = 1;
    m(5, 7) = 2;
    m(7, 15) = 3;
    m(9, 19) = 4;
    m(0, 0) = 5;
    m(9, 0) = 6;

    SparseMatrix<int>::Iterator it = m.SparseBegin();
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 0, 0, 5);
    ++it;
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 1, 1, 1);
    ++it;
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 5, 7, 2);
    ++it;
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 7, 15, 3);
    ++it;
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 9, 0, 6);
    ++it;
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 9, 19, 4);
    ++it;
    CPPUNIT_ASSERT(it == m.SparseEnd());

    it = m.Find(0, 0);
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 0, 0, 5);

    it = m.Find(7, 15);
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 7, 15, 3);

    it = m.Find(9, 19);
    ASSERT_ITERATOR_EQUAL(it, m.SparseEnd(), 9, 19, 4);
  }

  void
  TestOperators()
  {
    SparseMatrix<double> m;

    m(0, 0) = 1;
    m(1, 1) = 2;
    m(1, 4) = 3;
    m(2, 0) = 4;
    m(2, 4) = 5;

    SparseMatrix<double> mm(m);

    m *= 2;

    for (SparseMatrix<double>::Iterator it = m.SparseBegin();
         it != m.SparseEnd(); ++it)
      ASSERT_EQUAL(it(), 2 * mm[it.GetRow()][it.GetColumn()]);

    m += mm;

    for (SparseMatrix<double>::Iterator it = m.SparseBegin();
         it != m.SparseEnd(); ++it)
      ASSERT_EQUAL(it(), 3 * mm[it.GetRow()][it.GetColumn()]);

    m -= mm;

    for (SparseMatrix<double>::Iterator it = m.SparseBegin();
         it != m.SparseEnd(); ++it)
      ASSERT_EQUAL(it(), 2 * mm[it.GetRow()][it.GetColumn()]);

    m -= m;
    
    m.Reclaim();

    ASSERT_EQUAL(m.GetNonEmptySize(), 0U);
  }

  void
  TestMultiplication()
  {
    // 4 x 4
    SparseMatrix<double> a;

    a(0, 0) = 1; a(0, 2) = 2;
    a(1, 1) = 3; a(1, 3) = 4;
    a(2, 0) = 5; a(2, 2) = 6;
    a(3, 1) = 7; a(3, 3) = 8;

    SparseMatrix<double> b;

    b(0, 1) = 1;  b(0, 3) = -2;
    b(1, 0) = 4;  b(1, 2) = -8;
    b(2, 1) = 16; b(2, 3) = -32;
    b(3, 0) = 64; b(3, 2) = -128;

    SparseMatrix<double> r = a * b;

    SparseMatrix<double> c;

    c(0, 1) = 33;  c(0, 3) = -66;
    c(1, 0) = 268; c(1, 2) = -536;
    c(2, 1) = 101; c(2, 3) = -202;
    c(3, 0) = 540; c(3, 2) = -1080;

    CPPUNIT_ASSERT(r == c);
  }

  void
  TestNonSquareMultiplication()
  {
    // 2 x 10
    SparseMatrix<int> a;

    a(0, 0) = 1; a(0, 9) = 1;
    a(1, 0) = 1; a(1, 9) = 1;

    // 10 x 2
    SparseMatrix<float> b;
    b(0, 0) = 1; b(0, 1) = 1;
    b(9, 0) = 1; b(9, 1) = 1;

    // 2 x 2
    SparseMatrix<float> c = a * b;

    SparseMatrix<float> t;

    t(0, 0) = 2; t(0, 1) = 2;
    t(1, 0) = 2; t(1, 1) = 2;

    ASSERT_EQUAL(t, c);

    SparseMatrix<double> td(t);

    CPPUNIT_ASSERT(td == c);
  }

  void
  TestNonSquareMultiplication2()
  {
    const int n = 50;
    const int m = 25;

    SparseMatrix<> a;
    SparseMatrix<> b;

    for (int i = 0; i < n; ++i)
      for (int j = 0; j < m; ++j) {
        a(i, j) = i + j;
        b(j, i) = -2 * i + j;
    }

    const SparseMatrix<> c = a * b;

    for (int i = 0; i < n; ++i)
      for (int j = 0; j < n; ++j)
        ASSERT_EQUAL(c[i][j],
                     double(-2*m*i*j+(i-2*j)*(m-1)*m/2+(m-1)*m*(2*m-1)/6));
  }

  void
  TestTranspose()
  {
    // 3 x 10
    SparseMatrix<double> a;

    a(0, 0) = 1; a(0, 2) = 2; a(0, 9) = 1;
    a(1, 1) = 3; a(1, 3) = 4;
    a(2, 0) = 5;

    SparseMatrix<double> at(a);
    at.Transpose();

    SparseMatrix<double> t;

    t(0, 0) = 1; t(2, 0) = 2; t(9, 0) = 1;
    t(1, 1) = 3; t(3, 1) = 4;
    t(0, 2) = 5;

    CPPUNIT_ASSERT(t == at);

    at.Transpose();

    CPPUNIT_ASSERT(a == at);
  }

};


CPPUNIT_TEST_SUITE_REGISTRATION(TestSparseMatrix);