aboutsummaryrefslogtreecommitdiffhomepage
path: root/unsupported/test/cxx11_tensor_of_const_values.cpp
blob: 344d678ef9f76c8c811672cd1cf7a40e91b80c59 (plain)
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
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2014 Benoit Steiner <benoit.steiner.goog@gmail.com>
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.

#include "main.h"

#include <Eigen/CXX11/Tensor>

using Eigen::Tensor;
using Eigen::RowMajor;

static void test_assign()
{
  float data1[6];
  TensorMap<Tensor<const float, 2>> mat1(data1, 2, 3);
  float data2[6];
  const TensorMap<Tensor<float, 2>> mat2(data2, 2, 3);

  for (int i = 0; i < 6; ++i) {
    data1[i] = i;
    data2[i] = -i;
  }

  Tensor<float, 2> rslt1;
  rslt1 = mat1;
  Tensor<float, 2> rslt2;
  rslt2 = mat2;

  Tensor<float, 2> rslt3 = mat1;
  Tensor<float, 2> rslt4 = mat2;

  Tensor<float, 2> rslt5(mat1);
  Tensor<float, 2> rslt6(mat2);

  for (int i = 0; i < 2; ++i) {
    for (int j = 0; j < 3; ++j) {
      VERIFY_IS_APPROX(rslt1(i,j), static_cast<float>(i + 2*j));
      VERIFY_IS_APPROX(rslt2(i,j), static_cast<float>(-i - 2*j));
      VERIFY_IS_APPROX(rslt3(i,j), static_cast<float>(i + 2*j));
      VERIFY_IS_APPROX(rslt4(i,j), static_cast<float>(-i - 2*j));
      VERIFY_IS_APPROX(rslt5(i,j), static_cast<float>(i + 2*j));
      VERIFY_IS_APPROX(rslt6(i,j), static_cast<float>(-i - 2*j));
    }
  }
}


static void test_plus()
{
  float data1[6];
  TensorMap<Tensor<const float, 2>> mat1(data1, 2, 3);
  float data2[6];
  TensorMap<Tensor<float, 2>> mat2(data2, 2, 3);

  for (int i = 0; i < 6; ++i) {
    data1[i] = i;
    data2[i] = -i;
  }

  Tensor<float, 2> sum1;
  sum1 = mat1 + mat2;
  Tensor<float, 2> sum2;
  sum2 = mat2 + mat1;

  for (int i = 0; i < 2; ++i) {
    for (int j = 0; j < 3; ++j) {
      VERIFY_IS_APPROX(sum1(i,j), 0.0f);
      VERIFY_IS_APPROX(sum2(i,j), 0.0f);
    }
  }
}


static void test_plus_equal()
{
  float data1[6];
  TensorMap<Tensor<const float, 2>> mat1(data1, 2, 3);
  float data2[6];
  TensorMap<Tensor<float, 2>> mat2(data2, 2, 3);

  for (int i = 0; i < 6; ++i) {
    data1[i] = i;
    data2[i] = -i;
  }
  mat2 += mat1;

  for (int i = 0; i < 2; ++i) {
    for (int j = 0; j < 3; ++j) {
      VERIFY_IS_APPROX(mat2(i,j), 0.0f);
    }
  }
}


EIGEN_DECLARE_TEST(cxx11_tensor_of_const_values)
{
  CALL_SUBTEST(test_assign());
  CALL_SUBTEST(test_plus());
  CALL_SUBTEST(test_plus_equal());
}