Nektar++
TestFullMatrixOperations.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: TestFullMatrixOperations.cpp
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description:
32//
33///////////////////////////////////////////////////////////////////////////////
34
35#include <boost/test/unit_test.hpp>
36
37#include <LibUtilities/LinearAlgebra/NekMatrix.hpp>
38#include <UnitTests/CountedObject.h>
39
40#include <boost/core/ignore_unused.hpp>
41#include <boost/test/tools/floating_point_comparison.hpp>
42#include <boost/test/unit_test.hpp>
43
44namespace Nektar
45{
46namespace FullMatrixOperationsUnitTests
47{
48template <typename DataType, typename MatrixType>
49int foo(NekMatrix<DataType, MatrixType> &d)
50{
51 boost::ignore_unused(d);
52 return 1;
53}
54
55template <typename DataType> int foo(NekMatrix<DataType, BlockMatrixTag> &d)
56{
57 boost::ignore_unused(d);
58 return 2;
59}
60
61BOOST_AUTO_TEST_CASE(TestDoubleSquareFullVectorMultiplication)
62{
63 NekMatrix<double> m1;
64 NekMatrix<NekMatrix<double>, BlockMatrixTag> m2(2, 2, 2, 2);
65 foo(m1);
66 foo(m2);
67 double m_buf[] = {1, 2, 3, 4, 5, 6, 7, 8, 9};
68 double v_buf[] = {4, 5, 6};
69
70 NekMatrix<double> m(3, 3, m_buf);
71 NekVector<double> v(3, v_buf);
72
73 double expected_result_buf[] = {66, 81, 96};
74 NekVector<double> expected_result(3, expected_result_buf);
75
76 NekVector<double> result = m * v;
77
78 BOOST_CHECK_EQUAL(expected_result, result);
79}
80
81BOOST_AUTO_TEST_CASE(TestDoubleSquareFullVectorMultiplicationWithAliasing)
82{
83 double m_buf[] = {1, 2, 3, 4, 5, 6, 7, 8, 9};
84 double v_buf[] = {4, 5, 6};
85
86 NekMatrix<double> m(3, 3, m_buf);
87 NekVector<double> v(3, v_buf);
88
89 double expected_result_buf[] = {66, 81, 96};
90 NekVector<double> expected_result(3, expected_result_buf);
91
92 v = m * v;
93
94 BOOST_CHECK_EQUAL(expected_result, v);
95}
96
97BOOST_AUTO_TEST_CASE(
98 TestDoubleSquareFullVectorMultiplicationWithSharedArrayAliasing)
99{
100 double m_buf[] = {1, 2, 3, 4, 5, 6, 7, 8, 9};
101 double v_buf[] = {4, 5, 6};
102
103 Array<OneD, double> vector_array(3, v_buf);
104
105 NekMatrix<double> m(3, 3, m_buf);
106 NekVector<double> v(3, vector_array, eWrapper);
107 NekVector<double> result(3, vector_array, eWrapper);
108
109 double expected_result_buf[] = {66, 81, 96};
110 NekVector<double> expected_result(3, expected_result_buf);
111
112 result = m * v;
113
114 BOOST_CHECK_EQUAL(expected_result, result);
115 BOOST_CHECK_EQUAL(vector_array[0], 66);
116 BOOST_CHECK_EQUAL(vector_array[1], 81);
117 BOOST_CHECK_EQUAL(vector_array[2], 96);
118}
119
120BOOST_AUTO_TEST_CASE(
121 TestDoubleSquareFullVectorMultiplicationWithSharedArrayAliasingAndMatrixTranspose)
122{
123 double m_buf[] = {1, 4, 7, 2, 5, 8, 3, 6, 9};
124 double v_buf[] = {4, 5, 6};
125
126 Array<OneD, double> vector_array(3, v_buf);
127
128 NekMatrix<double> m(3, 3, m_buf);
129 NekVector<double> v(3, vector_array, eWrapper);
130 NekVector<double> result(3, vector_array, eWrapper);
131
132 double expected_result_buf[] = {66, 81, 96};
133 NekVector<double> expected_result(3, expected_result_buf);
134
135 result = Transpose(m) * v;
136
137 BOOST_CHECK_EQUAL(expected_result, result);
138 BOOST_CHECK_EQUAL(vector_array[0], 66);
139 BOOST_CHECK_EQUAL(vector_array[1], 81);
140 BOOST_CHECK_EQUAL(vector_array[2], 96);
141}
142
143BOOST_AUTO_TEST_CASE(
144 TestDoubleSquareFullVectorMultiplicationWithSharedArrayAliasingAndOverlap)
145{
146 double m_buf[] = {1, 2, 3, 4, 5, 6, 7, 8, 9};
147 double v_buf[] = {4, 5, 6, 0};
148
149 Array<OneD, double> vector_array(4, v_buf);
150 Array<OneD, double> tmp;
151
152 NekMatrix<double> m(3, 3, m_buf);
153 NekVector<double> v(3, vector_array, eWrapper);
154 NekVector<double> result(3, tmp = vector_array + 1, eWrapper);
155
156 double expected_result_buf[] = {66, 81, 96};
157 NekVector<double> expected_result(3, expected_result_buf);
158
159 result = m * v;
160
161 BOOST_CHECK_EQUAL(expected_result, result);
162 BOOST_CHECK_EQUAL(vector_array[1], 66);
163 BOOST_CHECK_EQUAL(vector_array[2], 81);
164 BOOST_CHECK_EQUAL(vector_array[3], 96);
165}
166
167BOOST_AUTO_TEST_CASE(TestScaledDoubleSquareFullVectorMultiplication)
168{
169 double m_buf[] = {1, 2, 3, 4, 5, 6, 7, 8, 9};
170 double v_buf[] = {4, 5, 6};
171
172 std::shared_ptr<NekMatrix<double>> inner(
173 new NekMatrix<double>(3, 3, m_buf));
174 NekMatrix<NekMatrix<double>, ScaledMatrixTag> m(7.0, inner);
175 NekVector<double> v(3, v_buf);
176
177 double expected_result_buf[] = {462, 567, 672};
178 NekVector<double> expected_result(3, expected_result_buf);
179
180 NekVector<double> result = m * v;
181
182 BOOST_CHECK_EQUAL(expected_result, result);
183}
184
185BOOST_AUTO_TEST_CASE(TestThreeMatrixMultiplication)
186{
187 double buf1[] = {1.0, 2.0, 3.0, 4.0};
188 double buf2[] = {5.0, 6.0, 7.0, 8.0};
189 double buf3[] = {9.0, 10.0, 11.0, 12.0};
190
191 NekMatrix<double> m1(2, 2, buf1);
192 NekMatrix<double> m2(2, 2, buf2);
193 NekMatrix<double> m3(2, 2, buf3);
194
195 double expected_result_buf[] = {517, 766, 625, 926};
196 NekMatrix<double> expected_result(2, 2, expected_result_buf);
197
198 NekMatrix<double> result1 = m1 * m2 * m3;
199 NekMatrix<double> result2 = m1 * (m2 * m3);
200 BOOST_CHECK_EQUAL(expected_result, m1 * m2 * m3);
201
202 BOOST_CHECK_EQUAL(expected_result, m1 * (m2 * m3));
203}
204
205BOOST_AUTO_TEST_CASE(TestThreeMatrixMultiplicationWithTranspose)
206{
207 double buf1[] = {1.0, 2.0, 3.0, 4.0};
208 double buf2[] = {5.0, 6.0, 7.0, 8.0};
209 double buf3[] = {9.0, 11.0, 10.0, 12.0};
210
211 NekMatrix<double> m1(2, 2, buf1);
212 NekMatrix<double> m2(2, 2, buf2);
213 NekMatrix<double> m3(2, 2, buf3);
214
215 double expected_result_buf[] = {517, 766, 625, 926};
216 NekMatrix<double> expected_result(2, 2, expected_result_buf);
217
218 BOOST_CHECK_EQUAL(expected_result, m1 * m2 * Transpose(m3));
219}
220
221BOOST_AUTO_TEST_CASE(TestThreeWrappedMatrixMultiplicationWithTranspose)
222{
223 double buf1[] = {1.0, 2.0, 3.0, 4.0};
224 double buf2[] = {5.0, 6.0, 7.0, 8.0};
225 double buf3[] = {9.0, 11.0, 10.0, 12.0};
226 double out_buf[] = {0.0, 0.0, 0.0, 0.0};
227
228 Array<OneD, double> array_buf1(4, buf1);
229 Array<OneD, double> array_buf2(4, buf2);
230 Array<OneD, double> array_buf3(4, buf3);
231 Array<OneD, double> array_out_buf(4, out_buf);
232
233 NekMatrix<double> m1(2, 2, array_buf1, eWrapper);
234 NekMatrix<double> m2(2, 2, array_buf2, eWrapper);
235 NekMatrix<double> m3(2, 2, array_buf3, eWrapper);
236 NekMatrix<double> result(2, 2, array_out_buf, eWrapper);
237
238 double expected_result_buf[] = {517, 766, 625, 926};
239 NekMatrix<double> expected_result(2, 2, expected_result_buf);
240
241 result = m1 * m2 * Transpose(m3);
242 BOOST_CHECK_EQUAL(expected_result, result);
243 BOOST_CHECK_EQUAL(array_out_buf[0], 517.0);
244 BOOST_CHECK_EQUAL(array_out_buf[1], 766.0);
245 BOOST_CHECK_EQUAL(array_out_buf[2], 625.0);
246 BOOST_CHECK_EQUAL(array_out_buf[3], 926.0);
247 BOOST_CHECK_EQUAL(result(0, 0), 517.0);
248 BOOST_CHECK_EQUAL(result(1, 0), 766.0);
249 BOOST_CHECK_EQUAL(result(0, 1), 625.0);
250 BOOST_CHECK_EQUAL(result(1, 1), 926.0);
251}
252} // namespace FullMatrixOperationsUnitTests
253} // namespace Nektar
Nektar::FullMatrixOperationsUnitTests::BOOST_AUTO_TEST_CASE
BOOST_AUTO_TEST_CASE(TestDoubleSquareFullVectorMultiplication)
Definition: TestFullMatrixOperations.cpp:61
Nektar::FullMatrixOperationsUnitTests::foo
int foo(NekMatrix< DataType, MatrixType > &d)
Definition: TestFullMatrixOperations.cpp:49
Nektar::UnitTests::d
std::vector< double > d(NPUPPER *NPUPPER)
Nektar
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:2
Nektar::Transpose
NekMatrix< InnerMatrixType, BlockMatrixTag > Transpose(NekMatrix< InnerMatrixType, BlockMatrixTag > &rhs)
Definition: BlockMatrix.hpp:235
Generated on Wed Jan 31 2024 08:26:49 for Nektar++ by doxygen 1.9.4