ViennaCL - The Vienna Computing Library  1.7.1
Free open-source GPU-accelerated linear algebra and solver library.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
matrix_float_double.hpp
Go to the documentation of this file.
1 /* =========================================================================
2  Copyright (c) 2010-2016, Institute for Microelectronics,
3  Institute for Analysis and Scientific Computing,
4  TU Wien.
5  Portions of this software are copyright by UChicago Argonne, LLC.
6 
7  -----------------
8  ViennaCL - The Vienna Computing Library
9  -----------------
10 
11  Project Head: Karl Rupp rupp@iue.tuwien.ac.at
12 
13  (A list of authors and contributors can be found in the PDF manual)
14 
15  License: MIT (X11), see file LICENSE in the base directory
16 ============================================================================= */
17 
18 #include <utility>
19 #include <iostream>
20 #include <fstream>
21 #include <string>
22 #include <cmath>
23 #include <algorithm>
24 #include <cstdio>
25 #include <ctime>
26 
27 #include "viennacl/scalar.hpp"
28 #include "viennacl/matrix.hpp"
29 #include "viennacl/linalg/prod.hpp"
30 #include "viennacl/linalg/norm_1.hpp"
31 #include "viennacl/linalg/norm_inf.hpp"
32 #include "viennacl/linalg/norm_frobenius.hpp"
33 #include "viennacl/matrix_proxy.hpp"
34 #include "viennacl/vector_proxy.hpp"
35 #include "viennacl/linalg/norm_1.hpp"
36 #include "viennacl/linalg/norm_2.hpp"
37 #include "viennacl/linalg/norm_inf.hpp"
38 #include "viennacl/linalg/norm_frobenius.hpp"
39 
40 template<typename NumericT, typename VCLMatrixType>
41 bool check_for_equality(std::vector<std::vector<NumericT> > const & std_A, VCLMatrixType const & vcl_A, double epsilon)
42 {
43  std::vector<std::vector<NumericT> > vcl_A_cpu(vcl_A.size1(), std::vector<NumericT>(vcl_A.size2()));
44  viennacl::backend::finish(); //workaround for a bug in APP SDK 2.7 on Trinity APUs (with Catalyst 12.8)
45  viennacl::copy(vcl_A, vcl_A_cpu);
46 
47  for (std::size_t i=0; i<std_A.size(); ++i)
48  {
49  for (std::size_t j=0; j<std_A[i].size(); ++j)
50  {
51  if (std::fabs(std_A[i][j] - vcl_A_cpu[i][j]) > 0)
52  {
53  if ( (std::abs(std_A[i][j] - vcl_A_cpu[i][j]) / std::max(std::fabs(std_A[i][j]), std::fabs(vcl_A_cpu[i][j])) > epsilon) || std::fabs(vcl_A_cpu[i][j] - vcl_A_cpu[i][j]) > 0 )
54  {
55  std::cout << "Error at index (" << i << ", " << j << "): " << std_A[i][j] << " vs " << vcl_A_cpu[i][j] << std::endl;
56  std::cout << std::endl << "TEST failed!" << std::endl;
57  return false;
58  }
59  }
60  }
61  }
62 
63  std::cout << "PASSED!" << std::endl;
64  return true;
65 }
66 
67 
68 
69 
70 template<typename STLMatrixType,
71  typename ViennaCLMatrixType1, typename ViennaCLMatrixType2, typename ViennaCLMatrixType3>
72 int run_test(double epsilon,
73  STLMatrixType & std_A, STLMatrixType & std_B, STLMatrixType & std_C,
74  ViennaCLMatrixType1 & vcl_A, ViennaCLMatrixType2 & vcl_B, ViennaCLMatrixType3 vcl_C)
75 {
76 
77  typedef typename viennacl::result_of::cpu_value_type<typename ViennaCLMatrixType1::value_type>::type cpu_value_type;
78 
79  cpu_value_type alpha = cpu_value_type(3.1415);
80  viennacl::scalar<cpu_value_type> gpu_alpha = alpha;
81 
82  cpu_value_type beta = cpu_value_type(2.7182);
83  viennacl::scalar<cpu_value_type> gpu_beta = beta;
84 
85 
86  //
87  // Initializer:
88  //
89  std::cout << "Checking for zero_matrix initializer..." << std::endl;
90  std_A = std::vector<std::vector<cpu_value_type> >(std_A.size(), std::vector<cpu_value_type>(std_A[0].size()));
91  vcl_A = viennacl::zero_matrix<cpu_value_type>(vcl_A.size1(), vcl_A.size2());
92  if (!check_for_equality(std_A, vcl_A, epsilon))
93  return EXIT_FAILURE;
94 
95  std::cout << "Checking for scalar_matrix initializer..." << std::endl;
96  std_A = std::vector<std::vector<cpu_value_type> >(std_A.size(), std::vector<cpu_value_type>(std_A[0].size(), alpha));
97  vcl_A = viennacl::scalar_matrix<cpu_value_type>(vcl_A.size1(), vcl_A.size2(), alpha);
98  if (!check_for_equality(std_A, vcl_A, epsilon))
99  return EXIT_FAILURE;
100 
101  std_A = std::vector<std::vector<cpu_value_type> >(std_A.size(), std::vector<cpu_value_type>(std_A[0].size(), gpu_beta));
102  vcl_A = viennacl::scalar_matrix<cpu_value_type>( vcl_A.size1(), vcl_A.size2(), gpu_beta);
103  if (!check_for_equality(std_A, vcl_A, epsilon))
104  return EXIT_FAILURE;
105 
106  /*std::cout << "Checking for identity initializer..." << std::endl;
107  std_A = ...;
108  vcl_A = viennacl::identity_matrix<cpu_value_type>(vcl_A.size1());
109  if (!check_for_equality(std_A, vcl_A, epsilon))
110  return EXIT_FAILURE;*/
111 
112 
113  std::cout << std::endl;
114  //std::cout << "//" << std::endl;
115  //std::cout << "////////// Test: Assignments //////////" << std::endl;
116  //std::cout << "//" << std::endl;
117 
118  if (!check_for_equality(std_B, vcl_B, epsilon))
119  return EXIT_FAILURE;
120 
121  std::cout << "Testing matrix assignment... ";
122  //std::cout << std_B(0,0) << " vs. " << vcl_B(0,0) << std::endl;
123  std_A = std_B;
124  vcl_A = vcl_B;
125  if (!check_for_equality(std_A, vcl_A, epsilon))
126  return EXIT_FAILURE;
127 
128  if (std_A.size() == std_A[0].size())
129  {
130  std::cout << "Testing matrix assignment (transposed)... ";
131  for (std::size_t i=0; i<std_A.size(); ++i)
132  for (std::size_t j=0; j<std_A[i].size(); ++j)
133  std_A[i][j] = std_B[j][i];
134  vcl_A = viennacl::trans(vcl_B);
135 
136  if (!check_for_equality(std_A, vcl_A, epsilon))
137  return EXIT_FAILURE;
138  }
139 
140 
141 
142  //std::cout << std::endl;
143  //std::cout << "//" << std::endl;
144  //std::cout << "////////// Test 1: Copy to GPU //////////" << std::endl;
145  //std::cout << "//" << std::endl;
146 
147  std_A = std_B;
148  viennacl::copy(std_B, vcl_A);
149  std::cout << "Testing upper left copy to GPU... ";
150  if (!check_for_equality(std_A, vcl_A, epsilon))
151  return EXIT_FAILURE;
152 
153 
154  std_C = std_B;
155  viennacl::copy(std_B, vcl_C);
156  std::cout << "Testing lower right copy to GPU... ";
157  if (!check_for_equality(std_C, vcl_C, epsilon))
158  return EXIT_FAILURE;
159 
160 
161  //std::cout << std::endl;
162  //std::cout << "//" << std::endl;
163  //std::cout << "////////// Test 2: Copy from GPU //////////" << std::endl;
164  //std::cout << "//" << std::endl;
165 
166  std::cout << "Testing upper left copy to A... ";
167  if (!check_for_equality(std_A, vcl_A, epsilon))
168  return EXIT_FAILURE;
169 
170  std::cout << "Testing lower right copy to C... ";
171  if (!check_for_equality(std_C, vcl_C, epsilon))
172  return EXIT_FAILURE;
173 
174 
175 
176  //std::cout << "//" << std::endl;
177  //std::cout << "////////// Test 3: Addition //////////" << std::endl;
178  //std::cout << "//" << std::endl;
179  viennacl::copy(std_C, vcl_C);
180 
181  std::cout << "Inplace add: ";
182  for (std::size_t i=0; i<std_C.size(); ++i)
183  for (std::size_t j=0; j<std_C[i].size(); ++j)
184  std_C[i][j] += std_C[i][j];
185  vcl_C += vcl_C;
186 
187  if (!check_for_equality(std_C, vcl_C, epsilon))
188  return EXIT_FAILURE;
189 
190  if (std_C.size() == std_C[0].size())
191  {
192  std::cout << "Inplace add (transposed): ";
193  STLMatrixType std_C2(std_C);
194  for (std::size_t i=0; i<std_C.size(); ++i)
195  for (std::size_t j=0; j<std_C[i].size(); ++j)
196  std_C[i][j] += std_C2[j][i];
197  vcl_C += viennacl::trans(vcl_C);
198 
199  if (!check_for_equality(std_C, vcl_C, epsilon))
200  return EXIT_FAILURE;
201 
202  std::cout << "Inplace add (transposed, expression): ";
203  std_C2 = std_C;
204  for (std::size_t i=0; i<std_C.size(); ++i)
205  for (std::size_t j=0; j<std_C[i].size(); ++j)
206  std_C[i][j] += std_C2[j][i] + std_C2[j][i];
207  vcl_C += viennacl::trans(vcl_C + vcl_C);
208 
209  if (!check_for_equality(std_C, vcl_C, epsilon))
210  return EXIT_FAILURE;
211  }
212 
213  std::cout << "Scaled inplace add: ";
214  for (std::size_t i=0; i<std_C.size(); ++i)
215  for (std::size_t j=0; j<std_C[i].size(); ++j)
216  std_C[i][j] += beta * std_A[i][j];
217  vcl_C += gpu_beta * vcl_A;
218 
219  if (!check_for_equality(std_C, vcl_C, epsilon))
220  return EXIT_FAILURE;
221 
222  std::cout << "Add: ";
223  for (std::size_t i=0; i<std_C.size(); ++i)
224  for (std::size_t j=0; j<std_C[i].size(); ++j)
225  std_C[i][j] = std_A[i][j] + std_B[i][j];
226  vcl_C = vcl_A + vcl_B;
227 
228  if (!check_for_equality(std_C, vcl_C, epsilon))
229  return EXIT_FAILURE;
230 
231  std::cout << "Add with flipsign: ";
232  for (std::size_t i=0; i<std_C.size(); ++i)
233  for (std::size_t j=0; j<std_C[i].size(); ++j)
234  std_C[i][j] = -std_A[i][j] + std_B[i][j];
235  vcl_C = - vcl_A + vcl_B;
236 
237  if (!check_for_equality(std_C, vcl_C, epsilon))
238  return EXIT_FAILURE;
239 
240 
241  std::cout << "Scaled add (left): ";
242  for (std::size_t i=0; i<std_C.size(); ++i)
243  for (std::size_t j=0; j<std_C[i].size(); ++j)
244  std_C[i][j] = cpu_value_type(long(alpha)) * std_A[i][j] + std_B[i][j];
245  vcl_C = long(alpha) * vcl_A + vcl_B;
246 
247  if (!check_for_equality(std_C, vcl_C, epsilon))
248  return EXIT_FAILURE;
249 
250  for (std::size_t i=0; i<std_C.size(); ++i)
251  for (std::size_t j=0; j<std_C[i].size(); ++j)
252  std_C[i][j] = cpu_value_type(float(alpha)) * std_A[i][j] + std_B[i][j];
253  vcl_C = float(alpha) * vcl_A + vcl_B;
254 
255  if (!check_for_equality(std_C, vcl_C, epsilon))
256  return EXIT_FAILURE;
257 
258  for (std::size_t i=0; i<std_C.size(); ++i)
259  for (std::size_t j=0; j<std_C[i].size(); ++j)
260  std_C[i][j] = cpu_value_type(double(alpha)) * std_A[i][j] + std_B[i][j];
261  vcl_C = double(alpha) * vcl_A + vcl_B;
262 
263  if (!check_for_equality(std_C, vcl_C, epsilon))
264  return EXIT_FAILURE;
265 
266  std::cout << "Scaled add (left): ";
267  vcl_C = gpu_alpha * vcl_A + vcl_B;
268  if (!check_for_equality(std_C, vcl_C, epsilon))
269  return EXIT_FAILURE;
270 
271 
272  std::cout << "Scaled add (right): ";
273  for (std::size_t i=0; i<std_C.size(); ++i)
274  for (std::size_t j=0; j<std_C[i].size(); ++j)
275  std_C[i][j] = std_A[i][j] + std_B[i][j] * cpu_value_type(long(beta));
276  vcl_C = vcl_A + vcl_B * long(beta);
277 
278  if (!check_for_equality(std_C, vcl_C, epsilon))
279  return EXIT_FAILURE;
280 
281  for (std::size_t i=0; i<std_C.size(); ++i)
282  for (std::size_t j=0; j<std_C[i].size(); ++j)
283  std_C[i][j] = std_A[i][j] + std_B[i][j] * cpu_value_type(float(beta));
284  vcl_C = vcl_A + vcl_B * float(beta);
285 
286  if (!check_for_equality(std_C, vcl_C, epsilon))
287  return EXIT_FAILURE;
288 
289  for (std::size_t i=0; i<std_C.size(); ++i)
290  for (std::size_t j=0; j<std_C[i].size(); ++j)
291  std_C[i][j] = std_A[i][j] + std_B[i][j] * cpu_value_type(double(beta));
292  vcl_C = vcl_A + vcl_B * double(beta);
293 
294  if (!check_for_equality(std_C, vcl_C, epsilon))
295  return EXIT_FAILURE;
296 
297  std::cout << "Scaled add (right): ";
298  vcl_C = vcl_A + gpu_beta * vcl_B;
299  if (!check_for_equality(std_C, vcl_C, epsilon))
300  return EXIT_FAILURE;
301 
302  std::cout << "Scaled add (right, with division): ";
303  for (std::size_t i=0; i<std_C.size(); ++i)
304  for (std::size_t j=0; j<std_C[i].size(); ++j)
305  std_C[i][j] = std_A[i][j] + std_B[i][j] / cpu_value_type(long(beta));
306  vcl_C = vcl_A + vcl_B / long(beta);
307 
308  if (!check_for_equality(std_C, vcl_C, epsilon))
309  return EXIT_FAILURE;
310 
311  for (std::size_t i=0; i<std_C.size(); ++i)
312  for (std::size_t j=0; j<std_C[i].size(); ++j)
313  std_C[i][j] = std_A[i][j] + std_B[i][j] / cpu_value_type(float(beta));
314  vcl_C = vcl_A + vcl_B / float(beta);
315 
316  if (!check_for_equality(std_C, vcl_C, epsilon))
317  return EXIT_FAILURE;
318 
319  for (std::size_t i=0; i<std_C.size(); ++i)
320  for (std::size_t j=0; j<std_C[i].size(); ++j)
321  std_C[i][j] = std_A[i][j] + std_B[i][j] / cpu_value_type(double(beta));
322  vcl_C = vcl_A + vcl_B / double(beta);
323 
324  if (!check_for_equality(std_C, vcl_C, epsilon))
325  return EXIT_FAILURE;
326 
327 
328  std::cout << "Scaled add (both): ";
329  for (std::size_t i=0; i<std_C.size(); ++i)
330  for (std::size_t j=0; j<std_C[i].size(); ++j)
331  std_C[i][j] = alpha * std_A[i][j] + beta * std_B[i][j];
332  vcl_C = alpha * vcl_A + beta * vcl_B;
333 
334  if (!check_for_equality(std_C, vcl_C, epsilon))
335  return EXIT_FAILURE;
336 
337  std::cout << "Scaled add (both): ";
338  vcl_C = gpu_alpha * vcl_A + gpu_beta * vcl_B;
339  if (!check_for_equality(std_C, vcl_C, epsilon))
340  return EXIT_FAILURE;
341 
342  //std::cout << "//" << std::endl;
343  //std::cout << "////////// Test 4: Subtraction //////////" << std::endl;
344  //std::cout << "//" << std::endl;
345  viennacl::copy(std_C, vcl_C);
346 
347  std::cout << "Inplace sub: ";
348  for (std::size_t i=0; i<std_C.size(); ++i)
349  for (std::size_t j=0; j<std_C[i].size(); ++j)
350  std_C[i][j] -= std_B[i][j];
351  vcl_C -= vcl_B;
352 
353  if (!check_for_equality(std_C, vcl_C, epsilon))
354  return EXIT_FAILURE;
355 
356  if (std_C.size() == std_C[0].size())
357  {
358  std::cout << "Inplace sub (transposed): ";
359  STLMatrixType std_C2(std_C);
360  for (std::size_t i=0; i<std_C.size(); ++i)
361  for (std::size_t j=0; j<std_C[i].size(); ++j)
362  std_C[i][j] -= cpu_value_type(2) * std_C2[j][i];
363  vcl_C -= cpu_value_type(2) * viennacl::trans(vcl_C);
364 
365  if (!check_for_equality(std_C, vcl_C, epsilon))
366  return EXIT_FAILURE;
367 
368  std::cout << "Inplace sub (transposed, expression): ";
369  std_C2 = std_C;
370  for (std::size_t i=0; i<std_C.size(); ++i)
371  for (std::size_t j=0; j<std_C[i].size(); ++j)
372  std_C[i][j] -= std_C2[j][i] + std_C2[j][i];
373  vcl_C -= viennacl::trans(vcl_C + vcl_C);
374 
375  if (!check_for_equality(std_C, vcl_C, epsilon))
376  return EXIT_FAILURE;
377  }
378 
379  std::cout << "Scaled Inplace sub: ";
380  for (std::size_t i=0; i<std_C.size(); ++i)
381  for (std::size_t j=0; j<std_C[i].size(); ++j)
382  std_C[i][j] -= alpha * std_B[i][j];
383  vcl_C -= alpha * vcl_B;
384 
385  if (!check_for_equality(std_C, vcl_C, epsilon))
386  return EXIT_FAILURE;
387 
388 
389 
390 
391  std::cout << "Sub: ";
392  for (std::size_t i=0; i<std_C.size(); ++i)
393  for (std::size_t j=0; j<std_C[i].size(); ++j)
394  std_C[i][j] = std_A[i][j] - std_B[i][j];
395  vcl_C = vcl_A - vcl_B;
396 
397  if (!check_for_equality(std_C, vcl_C, epsilon))
398  return EXIT_FAILURE;
399 
400  std::cout << "Scaled sub (left): ";
401  for (std::size_t i=0; i<std_C.size(); ++i)
402  for (std::size_t j=0; j<std_C[i].size(); ++j)
403  std_B[i][j] = alpha * std_A[i][j] - std_C[i][j];
404  vcl_B = alpha * vcl_A - vcl_C;
405 
406  if (!check_for_equality(std_B, vcl_B, epsilon))
407  return EXIT_FAILURE;
408 
409  std::cout << "Scaled sub (left): ";
410  vcl_B = gpu_alpha * vcl_A - vcl_C;
411  if (!check_for_equality(std_B, vcl_B, epsilon))
412  return EXIT_FAILURE;
413 
414 
415  std::cout << "Scaled sub (right): ";
416  for (std::size_t i=0; i<std_C.size(); ++i)
417  for (std::size_t j=0; j<std_C[i].size(); ++j)
418  std_B[i][j] = std_A[i][j] - beta * std_C[i][j];
419  vcl_B = vcl_A - vcl_C * beta;
420 
421  if (!check_for_equality(std_B, vcl_B, epsilon))
422  return EXIT_FAILURE;
423 
424  std::cout << "Scaled sub (right): ";
425  vcl_B = vcl_A - vcl_C * gpu_beta;
426  if (!check_for_equality(std_B, vcl_B, epsilon))
427  return EXIT_FAILURE;
428 
429 
430  std::cout << "Scaled sub (both): ";
431  for (std::size_t i=0; i<std_C.size(); ++i)
432  for (std::size_t j=0; j<std_C[i].size(); ++j)
433  std_B[i][j] = alpha * std_A[i][j] - beta * std_C[i][j];
434  vcl_B = alpha * vcl_A - vcl_C * beta;
435 
436  if (!check_for_equality(std_B, vcl_B, epsilon))
437  return EXIT_FAILURE;
438 
439  std::cout << "Scaled sub (both): ";
440  vcl_B = gpu_alpha * vcl_A - vcl_C * gpu_beta;
441  if (!check_for_equality(std_B, vcl_B, epsilon))
442  return EXIT_FAILURE;
443 
444 
445  std::cout << "Unary operator-: ";
446  for (std::size_t i=0; i<std_C.size(); ++i)
447  for (std::size_t j=0; j<std_C[i].size(); ++j)
448  std_C[i][j] = - std_A[i][j];
449  vcl_C = -vcl_A;
450 
451  if (!check_for_equality(std_C, vcl_C, epsilon))
452  return EXIT_FAILURE;
453 
454 
455 
456  //std::cout << "//" << std::endl;
457  //std::cout << "////////// Test 5: Scaling //////////" << std::endl;
458  //std::cout << "//" << std::endl;
459  viennacl::copy(std_A, vcl_A);
460 
461  std::cout << "Multiplication with CPU scalar: ";
462  for (std::size_t i=0; i<std_A.size(); ++i)
463  for (std::size_t j=0; j<std_A[i].size(); ++j)
464  std_A[i][j] *= cpu_value_type(long(alpha));
465  vcl_A *= long(alpha);
466 
467  if (!check_for_equality(std_A, vcl_A, epsilon))
468  return EXIT_FAILURE;
469 
470  for (std::size_t i=0; i<std_A.size(); ++i)
471  for (std::size_t j=0; j<std_A[i].size(); ++j)
472  std_A[i][j] *= cpu_value_type(float(alpha));
473  vcl_A *= float(alpha);
474 
475  if (!check_for_equality(std_A, vcl_A, epsilon))
476  return EXIT_FAILURE;
477 
478  for (std::size_t i=0; i<std_A.size(); ++i)
479  for (std::size_t j=0; j<std_A[i].size(); ++j)
480  std_A[i][j] *= cpu_value_type(double(alpha));
481  vcl_A *= double(alpha);
482 
483  if (!check_for_equality(std_A, vcl_A, epsilon))
484  return EXIT_FAILURE;
485 
486  std::cout << "Multiplication with GPU scalar: ";
487  for (std::size_t i=0; i<std_A.size(); ++i)
488  for (std::size_t j=0; j<std_A[i].size(); ++j)
489  std_A[i][j] *= beta;
490  vcl_A *= gpu_beta;
491 
492  if (!check_for_equality(std_A, vcl_A, epsilon))
493  return EXIT_FAILURE;
494 
495 
496  std::cout << "Division with CPU scalar: ";
497  for (std::size_t i=0; i<std_A.size(); ++i)
498  for (std::size_t j=0; j<std_A[i].size(); ++j)
499  std_A[i][j] /= cpu_value_type(long(alpha));
500  vcl_A /= long(alpha);
501 
502  if (!check_for_equality(std_A, vcl_A, epsilon))
503  return EXIT_FAILURE;
504 
505  for (std::size_t i=0; i<std_A.size(); ++i)
506  for (std::size_t j=0; j<std_A[i].size(); ++j)
507  std_A[i][j] /= cpu_value_type(float(alpha));
508  vcl_A /= float(alpha);
509 
510  if (!check_for_equality(std_A, vcl_A, epsilon))
511  return EXIT_FAILURE;
512 
513  for (std::size_t i=0; i<std_A.size(); ++i)
514  for (std::size_t j=0; j<std_A[i].size(); ++j)
515  std_A[i][j] /= cpu_value_type(double(alpha));
516  vcl_A /= double(alpha);
517 
518  if (!check_for_equality(std_A, vcl_A, epsilon))
519  return EXIT_FAILURE;
520 
521  std::cout << "Division with GPU scalar: ";
522  for (std::size_t i=0; i<std_A.size(); ++i)
523  for (std::size_t j=0; j<std_A[i].size(); ++j)
524  std_A[i][j] /= beta;
525  vcl_A /= gpu_beta;
526 
527  if (!check_for_equality(std_A, vcl_A, epsilon))
528  return EXIT_FAILURE;
529 
530 
531 
532  std::cout << "Testing elementwise multiplication..." << std::endl;
533  std_B = std::vector<std::vector<cpu_value_type> >(std_B.size(), std::vector<cpu_value_type>(std_B[0].size(), cpu_value_type(1.4142)));
534  for (std::size_t i=0; i<std_A.size(); ++i)
535  for (std::size_t j=0; j<std_A[i].size(); ++j)
536  std_A[i][j] = cpu_value_type(3.1415) * std_B[i][j];
537  viennacl::copy(std_A, vcl_A);
538  viennacl::copy(std_B, vcl_B);
539  viennacl::copy(std_B, vcl_B);
540  for (std::size_t i=0; i<std_A.size(); ++i)
541  for (std::size_t j=0; j<std_A[i].size(); ++j)
542  std_A[i][j] = std_A[i][j] * std_B[i][j];
543  vcl_A = viennacl::linalg::element_prod(vcl_A, vcl_B);
544 
545  if (!check_for_equality(std_A, vcl_A, epsilon))
546  return EXIT_FAILURE;
547 
548  for (std::size_t i=0; i<std_A.size(); ++i)
549  for (std::size_t j=0; j<std_A[i].size(); ++j)
550  std_A[i][j] += std_A[i][j] * std_B[i][j];
551  vcl_A += viennacl::linalg::element_prod(vcl_A, vcl_B);
552 
553  if (!check_for_equality(std_A, vcl_A, epsilon))
554  return EXIT_FAILURE;
555 
556  for (std::size_t i=0; i<std_A.size(); ++i)
557  for (std::size_t j=0; j<std_A[i].size(); ++j)
558  std_A[i][j] -= std_A[i][j] * std_B[i][j];
559  vcl_A -= viennacl::linalg::element_prod(vcl_A, vcl_B);
560 
561  if (!check_for_equality(std_A, vcl_A, epsilon))
562  return EXIT_FAILURE;
563 
565  for (std::size_t i=0; i<std_A.size(); ++i)
566  for (std::size_t j=0; j<std_A[i].size(); ++j)
567  std_A[i][j] = (std_A[i][j] + std_B[i][j]) * std_B[i][j];
568  vcl_A = viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B);
569 
570  if (!check_for_equality(std_A, vcl_A, epsilon))
571  return EXIT_FAILURE;
572 
573  for (std::size_t i=0; i<std_A.size(); ++i)
574  for (std::size_t j=0; j<std_A[i].size(); ++j)
575  std_A[i][j] += (std_A[i][j] + std_B[i][j]) * std_B[i][j];
576  vcl_A += viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B);
577 
578  if (!check_for_equality(std_A, vcl_A, epsilon))
579  return EXIT_FAILURE;
580 
581  for (std::size_t i=0; i<std_A.size(); ++i)
582  for (std::size_t j=0; j<std_A[i].size(); ++j)
583  std_A[i][j] -= (std_A[i][j] + std_B[i][j]) * std_B[i][j];
584  vcl_A -= viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B);
585 
586  if (!check_for_equality(std_A, vcl_A, epsilon))
587  return EXIT_FAILURE;
588 
590  for (std::size_t i=0; i<std_A.size(); ++i)
591  for (std::size_t j=0; j<std_A[i].size(); ++j)
592  std_A[i][j] = std_A[i][j] * (std_B[i][j] + std_A[i][j]);
593  vcl_A = viennacl::linalg::element_prod(vcl_A, vcl_B + vcl_A);
594 
595  if (!check_for_equality(std_A, vcl_A, epsilon))
596  return EXIT_FAILURE;
597 
598  for (std::size_t i=0; i<std_A.size(); ++i)
599  for (std::size_t j=0; j<std_A[i].size(); ++j)
600  std_A[i][j] += std_A[i][j] * (std_B[i][j] + std_A[i][j]);
601  vcl_A += viennacl::linalg::element_prod(vcl_A, vcl_B + vcl_A);
602 
603  if (!check_for_equality(std_A, vcl_A, epsilon))
604  return EXIT_FAILURE;
605 
606  for (std::size_t i=0; i<std_A.size(); ++i)
607  for (std::size_t j=0; j<std_A[i].size(); ++j)
608  std_A[i][j] -= std_A[i][j] * (std_B[i][j] + std_A[i][j]);
609  vcl_A -= viennacl::linalg::element_prod(vcl_A, vcl_B + vcl_A);
610 
611  if (!check_for_equality(std_A, vcl_A, epsilon))
612  return EXIT_FAILURE;
613 
615  for (std::size_t i=0; i<std_A.size(); ++i)
616  for (std::size_t j=0; j<std_A[i].size(); ++j)
617  std_A[i][j] = (std_A[i][j] + std_B[i][j]) * (std_B[i][j] + std_A[i][j]);
618  vcl_A = viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B + vcl_A);
619 
620  if (!check_for_equality(std_A, vcl_A, epsilon))
621  return EXIT_FAILURE;
622 
623  for (std::size_t i=0; i<std_A.size(); ++i)
624  for (std::size_t j=0; j<std_A[i].size(); ++j)
625  std_A[i][j] += (std_A[i][j] + std_B[i][j]) * (std_B[i][j] + std_A[i][j]);
626  vcl_A += viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B + vcl_A);
627 
628  if (!check_for_equality(std_A, vcl_A, epsilon))
629  return EXIT_FAILURE;
630 
631  for (std::size_t i=0; i<std_A.size(); ++i)
632  for (std::size_t j=0; j<std_A[i].size(); ++j)
633  std_A[i][j] -= (std_A[i][j] + std_B[i][j]) * (std_B[i][j] + std_A[i][j]);
634  vcl_A -= viennacl::linalg::element_prod(vcl_A + vcl_B, vcl_B + vcl_A);
635 
636  if (!check_for_equality(std_A, vcl_A, epsilon))
637  return EXIT_FAILURE;
638 
639 
640  std_B = std::vector<std::vector<cpu_value_type> >(std_B.size(), std::vector<cpu_value_type>(std_B[0].size(), cpu_value_type(1.4142)));
641  for (std::size_t i=0; i<std_A.size(); ++i)
642  for (std::size_t j=0; j<std_A[i].size(); ++j)
643  std_A[i][j] = cpu_value_type(3.1415) * std_B[i][j];
644  viennacl::copy(std_A, vcl_A);
645  viennacl::copy(std_B, vcl_B);
646 
647  for (std::size_t i=0; i<std_A.size(); ++i)
648  for (std::size_t j=0; j<std_A[i].size(); ++j)
649  std_A[i][j] = std_A[i][j] / std_B[i][j];
650  vcl_A = viennacl::linalg::element_div(vcl_A, vcl_B);
651 
652  if (!check_for_equality(std_A, vcl_A, epsilon))
653  return EXIT_FAILURE;
654 
655  for (std::size_t i=0; i<std_A.size(); ++i)
656  for (std::size_t j=0; j<std_A[i].size(); ++j)
657  std_A[i][j] += std_A[i][j] / std_B[i][j];
658  vcl_A += viennacl::linalg::element_div(vcl_A, vcl_B);
659 
660  if (!check_for_equality(std_A, vcl_A, epsilon))
661  return EXIT_FAILURE;
662 
663  for (std::size_t i=0; i<std_A.size(); ++i)
664  for (std::size_t j=0; j<std_A[i].size(); ++j)
665  std_A[i][j] -= std_A[i][j] / std_B[i][j];
666  vcl_A -= viennacl::linalg::element_div(vcl_A, vcl_B);
667 
668  if (!check_for_equality(std_A, vcl_A, epsilon))
669  return EXIT_FAILURE;
670 
672  for (std::size_t i=0; i<std_A.size(); ++i)
673  for (std::size_t j=0; j<std_A[i].size(); ++j)
674  std_A[i][j] = (std_A[i][j] + std_B[i][j]) / std_B[i][j];
675  vcl_A = viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B);
676 
677  if (!check_for_equality(std_A, vcl_A, epsilon))
678  return EXIT_FAILURE;
679 
680  for (std::size_t i=0; i<std_A.size(); ++i)
681  for (std::size_t j=0; j<std_A[i].size(); ++j)
682  std_A[i][j] += (std_A[i][j] + std_B[i][j]) / std_B[i][j];
683  vcl_A += viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B);
684 
685  if (!check_for_equality(std_A, vcl_A, epsilon))
686  return EXIT_FAILURE;
687 
688  for (std::size_t i=0; i<std_A.size(); ++i)
689  for (std::size_t j=0; j<std_A[i].size(); ++j)
690  std_A[i][j] -= (std_A[i][j] + std_B[i][j]) / std_B[i][j];
691  vcl_A -= viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B);
692 
693  if (!check_for_equality(std_A, vcl_A, epsilon))
694  return EXIT_FAILURE;
695 
697  for (std::size_t i=0; i<std_A.size(); ++i)
698  for (std::size_t j=0; j<std_A[i].size(); ++j)
699  std_A[i][j] = std_A[i][j] / (std_B[i][j] + std_A[i][j]);
700  vcl_A = viennacl::linalg::element_div(vcl_A, vcl_B + vcl_A);
701 
702  if (!check_for_equality(std_A, vcl_A, epsilon))
703  return EXIT_FAILURE;
704 
705  for (std::size_t i=0; i<std_A.size(); ++i)
706  for (std::size_t j=0; j<std_A[i].size(); ++j)
707  std_A[i][j] += std_A[i][j] / (std_B[i][j] + std_A[i][j]);
708  vcl_A += viennacl::linalg::element_div(vcl_A, vcl_B + vcl_A);
709 
710  if (!check_for_equality(std_A, vcl_A, epsilon))
711  return EXIT_FAILURE;
712 
713  for (std::size_t i=0; i<std_A.size(); ++i)
714  for (std::size_t j=0; j<std_A[i].size(); ++j)
715  std_A[i][j] -= std_A[i][j] / (std_B[i][j] + std_A[i][j]);
716  vcl_A -= viennacl::linalg::element_div(vcl_A, vcl_B + vcl_A);
717 
718  if (!check_for_equality(std_A, vcl_A, epsilon))
719  return EXIT_FAILURE;
720 
722  for (std::size_t i=0; i<std_A.size(); ++i)
723  for (std::size_t j=0; j<std_A[i].size(); ++j)
724  std_A[i][j] = (std_A[i][j] + std_B[i][j]) / (std_B[i][j] + std_A[i][j]);
725  vcl_A = viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B + vcl_A);
726 
727  if (!check_for_equality(std_A, vcl_A, epsilon))
728  return EXIT_FAILURE;
729 
730  for (std::size_t i=0; i<std_A.size(); ++i)
731  for (std::size_t j=0; j<std_A[i].size(); ++j)
732  std_A[i][j] += (std_A[i][j] + std_B[i][j]) / (std_B[i][j] + std_A[i][j]);
733  vcl_A += viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B + vcl_A);
734 
735  if (!check_for_equality(std_A, vcl_A, epsilon))
736  return EXIT_FAILURE;
737 
738  for (std::size_t i=0; i<std_A.size(); ++i)
739  for (std::size_t j=0; j<std_A[i].size(); ++j)
740  std_A[i][j] -= (std_A[i][j] + std_B[i][j]) / (std_B[i][j] + std_A[i][j]);
741  vcl_A -= viennacl::linalg::element_div(vcl_A + vcl_B, vcl_B + vcl_A);
742 
743  if (!check_for_equality(std_A, vcl_A, epsilon))
744  return EXIT_FAILURE;
745 
746  // element_pow
747  std::cout << "Testing unary element_pow()..." << std::endl;
748 
749  std_B = std::vector<std::vector<cpu_value_type> >(std_B.size(), std::vector<cpu_value_type>(std_B[0].size(), cpu_value_type(1.4142)));
750  for (std::size_t i=0; i<std_A.size(); ++i)
751  for (std::size_t j=0; j<std_A[i].size(); ++j)
752  std_A[i][j] = cpu_value_type(3.1415) * std_B[i][j];
753  viennacl::copy(std_A, vcl_A);
754  viennacl::copy(std_B, vcl_B);
755 
756  for (std::size_t i=0; i<std_C.size(); ++i)
757  for (std::size_t j=0; j<std_C[i].size(); ++j)
758  std_C[i][j] = std::pow(std_A[i][j], std_B[i][j]);
759  vcl_C = viennacl::linalg::element_pow(vcl_A, vcl_B);
760 
761  if (!check_for_equality(std_C, vcl_C, epsilon))
762  return EXIT_FAILURE;
763 
764  for (std::size_t i=0; i<std_C.size(); ++i)
765  for (std::size_t j=0; j<std_C[i].size(); ++j)
766  std_C[i][j] += std::pow(std_A[i][j], std_B[i][j]);
767  vcl_C += viennacl::linalg::element_pow(vcl_A, vcl_B);
768 
769  if (!check_for_equality(std_C, vcl_C, epsilon))
770  return EXIT_FAILURE;
771 
772  for (std::size_t i=0; i<std_C.size(); ++i)
773  for (std::size_t j=0; j<std_C[i].size(); ++j)
774  std_C[i][j] -= std::pow(std_A[i][j], std_B[i][j]);
775  vcl_C -= viennacl::linalg::element_pow(vcl_A, vcl_B);
776 
777  if (!check_for_equality(std_C, vcl_C, epsilon))
778  return EXIT_FAILURE;
779 
781  for (std::size_t i=0; i<std_C.size(); ++i)
782  for (std::size_t j=0; j<std_C[i].size(); ++j)
783  std_C[i][j] = std::pow(std_A[i][j] + std_B[i][j], std_B[i][j]);
784  vcl_C = viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B);
785 
786  if (!check_for_equality(std_C, vcl_C, epsilon))
787  return EXIT_FAILURE;
788 
789  for (std::size_t i=0; i<std_C.size(); ++i)
790  for (std::size_t j=0; j<std_C[i].size(); ++j)
791  std_C[i][j] += std::pow(std_A[i][j] + std_B[i][j], std_B[i][j]);
792  vcl_C += viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B);
793 
794  if (!check_for_equality(std_C, vcl_C, epsilon))
795  return EXIT_FAILURE;
796 
797  for (std::size_t i=0; i<std_C.size(); ++i)
798  for (std::size_t j=0; j<std_C[i].size(); ++j)
799  std_C[i][j] -= std::pow(std_A[i][j] + std_B[i][j], std_B[i][j]);
800  vcl_C -= viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B);
801 
802  if (!check_for_equality(std_C, vcl_C, epsilon))
803  return EXIT_FAILURE;
804 
806  for (std::size_t i=0; i<std_C.size(); ++i)
807  for (std::size_t j=0; j<std_C[i].size(); ++j)
808  std_C[i][j] = std::pow(std_A[i][j], std_B[i][j] + std_A[i][j]);
809  vcl_C = viennacl::linalg::element_pow(vcl_A, vcl_B + vcl_A);
810 
811  if (!check_for_equality(std_C, vcl_C, epsilon))
812  return EXIT_FAILURE;
813 
814  for (std::size_t i=0; i<std_C.size(); ++i)
815  for (std::size_t j=0; j<std_C[i].size(); ++j)
816  std_C[i][j] += std::pow(std_A[i][j], std_B[i][j] + std_A[i][j]);
817  vcl_C += viennacl::linalg::element_pow(vcl_A, vcl_B + vcl_A);
818 
819  if (!check_for_equality(std_C, vcl_C, epsilon))
820  return EXIT_FAILURE;
821 
822  for (std::size_t i=0; i<std_C.size(); ++i)
823  for (std::size_t j=0; j<std_C[i].size(); ++j)
824  std_C[i][j] -= std::pow(std_A[i][j], std_B[i][j] + std_A[i][j]);
825  vcl_C -= viennacl::linalg::element_pow(vcl_A, vcl_B + vcl_A);
826 
827  if (!check_for_equality(std_C, vcl_C, epsilon))
828  return EXIT_FAILURE;
829 
831  for (std::size_t i=0; i<std_C.size(); ++i)
832  for (std::size_t j=0; j<std_C[i].size(); ++j)
833  std_C[i][j] = std::pow(std_A[i][j] + std_B[i][j], std_B[i][j] + std_A[i][j]);
834  vcl_C = viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B + vcl_A);
835 
836  if (!check_for_equality(std_C, vcl_C, epsilon))
837  return EXIT_FAILURE;
838 
839  for (std::size_t i=0; i<std_C.size(); ++i)
840  for (std::size_t j=0; j<std_C[i].size(); ++j)
841  std_C[i][j] += std::pow(std_A[i][j] + std_B[i][j], std_B[i][j] + std_A[i][j]);
842  vcl_C += viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B + vcl_A);
843 
844  if (!check_for_equality(std_C, vcl_C, epsilon))
845  return EXIT_FAILURE;
846 
847  for (std::size_t i=0; i<std_C.size(); ++i)
848  for (std::size_t j=0; j<std_C[i].size(); ++j)
849  std_C[i][j] -= std::pow(std_A[i][j] + std_B[i][j], std_B[i][j] + std_A[i][j]);
850  vcl_C -= viennacl::linalg::element_pow(vcl_A + vcl_B, vcl_B + vcl_A);
851 
852  if (!check_for_equality(std_C, vcl_C, epsilon))
853  return EXIT_FAILURE;
854 
855 
856  std::cout << "Testing unary elementwise operations..." << std::endl;
857 
858 #define GENERATE_UNARY_OP_TEST(FUNCNAME) \
859  std_B = std::vector<std::vector<cpu_value_type> >(std_B.size(), std::vector<cpu_value_type>(std_B[0].size(), cpu_value_type(1.4142))); \
860  for (std::size_t i=0; i<std_C.size(); ++i) \
861  for (std::size_t j=0; j<std_C[i].size(); ++j) \
862  std_A[i][j] = cpu_value_type(3.1415) * std_B[i][j]; \
863  for (std::size_t i=0; i<std_C.size(); ++i) \
864  for (std::size_t j=0; j<std_C[i].size(); ++j) \
865  std_C[i][j] = cpu_value_type(2.7172) * std_A[i][j]; \
866  viennacl::copy(std_A, vcl_A); \
867  viennacl::copy(std_B, vcl_B); \
868  viennacl::copy(std_C, vcl_C); \
869  viennacl::copy(std_B, vcl_B); \
870  \
871  for (std::size_t i=0; i<std_C.size(); ++i) \
872  for (std::size_t j=0; j<std_C[i].size(); ++j) \
873  std_C[i][j] = std::FUNCNAME(std_A[i][j]); \
874  vcl_C = viennacl::linalg::element_##FUNCNAME(vcl_A); \
875  \
876  if (!check_for_equality(std_C, vcl_C, epsilon)) \
877  { \
878  std::cout << "Failure at C = " << #FUNCNAME << "(A)" << std::endl; \
879  return EXIT_FAILURE; \
880  } \
881  \
882  for (std::size_t i=0; i<std_C.size(); ++i) \
883  for (std::size_t j=0; j<std_C[i].size(); ++j) \
884  std_C[i][j] = std::FUNCNAME(std_A[i][j] + std_B[i][j]); \
885  vcl_C = viennacl::linalg::element_##FUNCNAME(vcl_A + vcl_B); \
886  \
887  if (!check_for_equality(std_C, vcl_C, epsilon)) \
888  { \
889  std::cout << "Failure at C = " << #FUNCNAME << "(A + B)" << std::endl; \
890  return EXIT_FAILURE; \
891  } \
892  \
893  for (std::size_t i=0; i<std_C.size(); ++i) \
894  for (std::size_t j=0; j<std_C[i].size(); ++j) \
895  std_C[i][j] += std::FUNCNAME(std_A[i][j]); \
896  vcl_C += viennacl::linalg::element_##FUNCNAME(vcl_A); \
897  \
898  if (!check_for_equality(std_C, vcl_C, epsilon)) \
899  { \
900  std::cout << "Failure at C += " << #FUNCNAME << "(A)" << std::endl; \
901  return EXIT_FAILURE; \
902  } \
903  \
904  for (std::size_t i=0; i<std_C.size(); ++i) \
905  for (std::size_t j=0; j<std_C[i].size(); ++j) \
906  std_C[i][j] += std::FUNCNAME(std_A[i][j] + std_B[i][j]); \
907  vcl_C += viennacl::linalg::element_##FUNCNAME(vcl_A + vcl_B); \
908  \
909  if (!check_for_equality(std_C, vcl_C, epsilon)) \
910  { \
911  std::cout << "Failure at C += " << #FUNCNAME << "(A + B)" << std::endl; \
912  return EXIT_FAILURE; \
913  } \
914  \
915  for (std::size_t i=0; i<std_C.size(); ++i) \
916  for (std::size_t j=0; j<std_C[i].size(); ++j) \
917  std_C[i][j] -= std::FUNCNAME(std_A[i][j]); \
918  vcl_C -= viennacl::linalg::element_##FUNCNAME(vcl_A); \
919  \
920  if (!check_for_equality(std_C, vcl_C, epsilon)) \
921  { \
922  std::cout << "Failure at C -= " << #FUNCNAME << "(A)" << std::endl; \
923  return EXIT_FAILURE; \
924  } \
925  \
926  for (std::size_t i=0; i<std_C.size(); ++i) \
927  for (std::size_t j=0; j<std_C[i].size(); ++j) \
928  std_C[i][j] -= std::FUNCNAME(std_A[i][j] + std_B[i][j]); \
929  vcl_C -= viennacl::linalg::element_##FUNCNAME(vcl_A + vcl_B); \
930  \
931  if (!check_for_equality(std_C, vcl_C, epsilon)) \
932  { \
933  std::cout << "Failure at C -= " << #FUNCNAME << "(A + B)" << std::endl; \
934  return EXIT_FAILURE; \
935  } \
936  \
937 
938  GENERATE_UNARY_OP_TEST(cos);
939  GENERATE_UNARY_OP_TEST(cosh);
940  GENERATE_UNARY_OP_TEST(exp);
941  GENERATE_UNARY_OP_TEST(floor);
942  GENERATE_UNARY_OP_TEST(fabs);
943  GENERATE_UNARY_OP_TEST(log);
944  GENERATE_UNARY_OP_TEST(log10);
945  GENERATE_UNARY_OP_TEST(sin);
946  GENERATE_UNARY_OP_TEST(sinh);
947  GENERATE_UNARY_OP_TEST(fabs);
948  GENERATE_UNARY_OP_TEST(sqrt);
949  GENERATE_UNARY_OP_TEST(tan);
950  GENERATE_UNARY_OP_TEST(tanh);
951 
952  std::cout << "Complicated expressions: ";
953  //std::cout << "std_A: " << std_A << std::endl;
954  //std::cout << "std_B: " << std_B << std::endl;
955  //std::cout << "std_C: " << std_C << std::endl;
956  for (std::size_t i=0; i<std_C.size(); ++i)
957  for (std::size_t j=0; j<std_C[i].size(); ++j)
958  std_B[i][j] += alpha * (- std_A[i][j] - beta * std_C[i][j] + std_A[i][j]);
959  vcl_B += gpu_alpha * (- vcl_A - vcl_C * beta + vcl_A);
960 
961  if (!check_for_equality(std_B, vcl_B, epsilon))
962  return EXIT_FAILURE;
963 
964  for (std::size_t i=0; i<std_C.size(); ++i)
965  for (std::size_t j=0; j<std_C[i].size(); ++j)
966  std_B[i][j] += (- std_A[i][j] - beta * std_C[i][j] + std_A[i][j] * beta) / gpu_alpha;
967  vcl_B += (- vcl_A - vcl_C * beta + gpu_beta * vcl_A) / gpu_alpha;
968 
969  if (!check_for_equality(std_B, vcl_B, epsilon))
970  return EXIT_FAILURE;
971 
972 
973  for (std::size_t i=0; i<std_C.size(); ++i)
974  for (std::size_t j=0; j<std_C[i].size(); ++j)
975  std_B[i][j] -= alpha * (- std_A[i][j] - beta * std_C[i][j] - std_A[i][j]);
976  vcl_B -= gpu_alpha * (- vcl_A - vcl_C * beta - vcl_A);
977 
978  if (!check_for_equality(std_B, vcl_B, epsilon))
979  return EXIT_FAILURE;
980 
981  for (std::size_t i=0; i<std_C.size(); ++i)
982  for (std::size_t j=0; j<std_C[i].size(); ++j)
983  std_B[i][j] -= (- std_A[i][j] - beta * std_C[i][j] - std_A[i][j] * beta) / alpha;
984  vcl_B -= (- vcl_A - vcl_C * beta - gpu_beta * vcl_A) / gpu_alpha;
985 
986  if (!check_for_equality(std_B, vcl_B, epsilon))
987  return EXIT_FAILURE;
988 
989  std::cout << std::endl;
990  std::cout << "----------------------------------------------" << std::endl;
991  std::cout << std::endl;
992 
993 
994  return EXIT_SUCCESS;
995 }
996 
997 
998 
999 
1000 template<typename T, typename ScalarType>
1001 int run_test(double epsilon)
1002 {
1003 
1004  typedef viennacl::matrix<ScalarType, T> VCLMatrixType;
1005 
1006  std::size_t dim_rows = 131;
1007  std::size_t dim_cols = 33;
1008  //std::size_t dim_rows = 5;
1009  //std::size_t dim_cols = 3;
1010 
1011  //setup std objects:
1012  std::vector<std::vector<ScalarType> > std_A(dim_rows, std::vector<ScalarType>(dim_cols));
1013  std::vector<std::vector<ScalarType> > std_B(dim_rows, std::vector<ScalarType>(dim_cols));
1014  std::vector<std::vector<ScalarType> > std_C(dim_rows, std::vector<ScalarType>(dim_cols));
1015 
1016  for (std::size_t i=0; i<std_A.size(); ++i)
1017  for (std::size_t j=0; j<std_A[i].size(); ++j)
1018  {
1019  std_A[i][j] = ScalarType((i+2) + (j+1)*(i+2));
1020  std_B[i][j] = ScalarType((j+2) + (j+1)*(j+2));
1021  std_C[i][j] = ScalarType((i+1) + (i+1)*(i+2));
1022  }
1023 
1024  std::vector<std::vector<ScalarType> > std_A_large(4 * dim_rows, std::vector<ScalarType>(4 * dim_cols));
1025  for (std::size_t i=0; i<std_A_large.size(); ++i)
1026  for (std::size_t j=0; j<std_A_large[i].size(); ++j)
1027  std_A_large[i][j] = ScalarType(i * std_A_large[i].size() + j);
1028 
1029  //Setup ViennaCL objects
1030  VCLMatrixType vcl_A_full(4 * dim_rows, 4 * dim_cols);
1031  VCLMatrixType vcl_B_full(4 * dim_rows, 4 * dim_cols);
1032  VCLMatrixType vcl_C_full(4 * dim_rows, 4 * dim_cols);
1033 
1034  viennacl::copy(std_A_large, vcl_A_full);
1035  viennacl::copy(std_A_large, vcl_B_full);
1036  viennacl::copy(std_A_large, vcl_C_full);
1037 
1038  //
1039  // Create A
1040  //
1041  VCLMatrixType vcl_A(dim_rows, dim_cols);
1042 
1043  viennacl::range vcl_A_r1(2 * dim_rows, 3 * dim_rows);
1044  viennacl::range vcl_A_r2(dim_cols, 2 * dim_cols);
1045  viennacl::matrix_range<VCLMatrixType> vcl_range_A(vcl_A_full, vcl_A_r1, vcl_A_r2);
1046 
1047  viennacl::slice vcl_A_s1(2, 3, dim_rows);
1048  viennacl::slice vcl_A_s2(2 * dim_cols, 2, dim_cols);
1049  viennacl::matrix_slice<VCLMatrixType> vcl_slice_A(vcl_A_full, vcl_A_s1, vcl_A_s2);
1050 
1051 
1052  //
1053  // Create B
1054  //
1055  VCLMatrixType vcl_B(dim_rows, dim_cols);
1056 
1057  viennacl::range vcl_B_r1(dim_rows, 2 * dim_rows);
1058  viennacl::range vcl_B_r2(2 * dim_cols, 3 * dim_cols);
1059  viennacl::matrix_range<VCLMatrixType> vcl_range_B(vcl_B_full, vcl_B_r1, vcl_B_r2);
1060 
1061  viennacl::slice vcl_B_s1(2 * dim_rows, 2, dim_rows);
1062  viennacl::slice vcl_B_s2(dim_cols, 3, dim_cols);
1063  viennacl::matrix_slice<VCLMatrixType> vcl_slice_B(vcl_B_full, vcl_B_s1, vcl_B_s2);
1064 
1065 
1066  //
1067  // Create C
1068  //
1069  VCLMatrixType vcl_C(dim_rows, dim_cols);
1070 
1071  viennacl::range vcl_C_r1(2 * dim_rows, 3 * dim_rows);
1072  viennacl::range vcl_C_r2(3 * dim_cols, 4 * dim_cols);
1073  viennacl::matrix_range<VCLMatrixType> vcl_range_C(vcl_C_full, vcl_C_r1, vcl_C_r2);
1074 
1075  viennacl::slice vcl_C_s1(dim_rows, 2, dim_rows);
1076  viennacl::slice vcl_C_s2(0, 3, dim_cols);
1077  viennacl::matrix_slice<VCLMatrixType> vcl_slice_C(vcl_C_full, vcl_C_s1, vcl_C_s2);
1078 
1079  viennacl::copy(std_A, vcl_A);
1080  viennacl::copy(std_A, vcl_range_A);
1081  viennacl::copy(std_A, vcl_slice_A);
1082 
1083  viennacl::copy(std_B, vcl_B);
1084  viennacl::copy(std_B, vcl_range_B);
1085  viennacl::copy(std_B, vcl_slice_B);
1086 
1087  viennacl::copy(std_C, vcl_C);
1088  viennacl::copy(std_C, vcl_range_C);
1089  viennacl::copy(std_C, vcl_slice_C);
1090 
1091 
1092  std::cout << std::endl;
1093  std::cout << "//" << std::endl;
1094  std::cout << "////////// Test: Copy CTOR //////////" << std::endl;
1095  std::cout << "//" << std::endl;
1096 
1097  {
1098  std::cout << "Testing matrix created from range... ";
1099  VCLMatrixType vcl_temp = vcl_range_A;
1100  if (check_for_equality(std_A, vcl_temp, epsilon))
1101  std::cout << "PASSED!" << std::endl;
1102  else
1103  {
1104  std::cout << "vcl_temp: " << vcl_temp << std::endl;
1105  std::cout << "vcl_range_A: " << vcl_range_A << std::endl;
1106  std::cout << "vcl_A: " << vcl_A << std::endl;
1107  std::cout << std::endl << "TEST failed!" << std::endl;
1108  return EXIT_FAILURE;
1109  }
1110 
1111  std::cout << "Testing matrix created from slice... ";
1112  VCLMatrixType vcl_temp2 = vcl_range_B;
1113  if (check_for_equality(std_B, vcl_temp2, epsilon))
1114  std::cout << "PASSED!" << std::endl;
1115  else
1116  {
1117  std::cout << std::endl << "TEST failed!" << std::endl;
1118  return EXIT_FAILURE;
1119  }
1120  }
1121 
1122  std::cout << "//" << std::endl;
1123  std::cout << "////////// Test: Initializer for matrix type //////////" << std::endl;
1124  std::cout << "//" << std::endl;
1125 
1126  {
1127  std::vector<std::vector<ScalarType> > std_dummy1(std_A.size(), std::vector<ScalarType>(std_A.size()));
1128  for (std::size_t i=0; i<std_A.size(); ++i) std_dummy1[i][i] = ScalarType(1);
1129  std::vector<std::vector<ScalarType> > std_dummy2(std_A.size(), std::vector<ScalarType>(std_A.size(), ScalarType(3)));
1130  std::vector<std::vector<ScalarType> > std_dummy3(std_A.size(), std::vector<ScalarType>(std_A.size()));
1131 
1132  viennacl::matrix<ScalarType> vcl_dummy1 = viennacl::identity_matrix<ScalarType>(std_A.size());
1133  viennacl::matrix<ScalarType> vcl_dummy2 = viennacl::scalar_matrix<ScalarType>(std_A.size(), std_A.size(), 3.0);
1134  viennacl::matrix<ScalarType> vcl_dummy3 = viennacl::zero_matrix<ScalarType>(std_A.size(), std_A.size());
1135 
1136  std::cout << "Testing initializer CTOR... ";
1137  if ( check_for_equality(std_dummy1, vcl_dummy1, epsilon)
1138  && check_for_equality(std_dummy2, vcl_dummy2, epsilon)
1139  && check_for_equality(std_dummy3, vcl_dummy3, epsilon)
1140  )
1141  std::cout << "PASSED!" << std::endl;
1142  else
1143  {
1144  std::cout << std::endl << "TEST failed!" << std::endl;
1145  return EXIT_FAILURE;
1146  }
1147 
1148  std_dummy1 = std::vector<std::vector<ScalarType> >(std_A.size(), std::vector<ScalarType>(std_A.size()));
1149  std_dummy2 = std::vector<std::vector<ScalarType> >(std_A.size(), std::vector<ScalarType>(std_A.size()));
1150  for (std::size_t i=0; i<std_A.size(); ++i) std_dummy2[i][i] = ScalarType(1);
1151  std_dummy3 = std::vector<std::vector<ScalarType> >(std_A.size(), std::vector<ScalarType>(std_A.size(), 3));
1152 
1153  vcl_dummy1 = viennacl::zero_matrix<ScalarType>(std_A.size(), std_A.size());
1154  vcl_dummy2 = viennacl::identity_matrix<ScalarType>(std_A.size());
1155  vcl_dummy3 = viennacl::scalar_matrix<ScalarType>(std_A.size(), std_A.size(), 3.0);
1156 
1157  std::cout << "Testing initializer assignment... ";
1158  if ( check_for_equality(std_dummy1, vcl_dummy1, epsilon)
1159  && check_for_equality(std_dummy2, vcl_dummy2, epsilon)
1160  && check_for_equality(std_dummy3, vcl_dummy3, epsilon)
1161  )
1162  std::cout << "PASSED!" << std::endl;
1163  else
1164  {
1165  std::cout << std::endl << "TEST failed!" << std::endl;
1166  return EXIT_FAILURE;
1167  }
1168  }
1169 
1170  std::cout << "//" << std::endl;
1171  std::cout << "////////// Test: Norms //////////" << std::endl;
1172  std::cout << "//" << std::endl;
1173 
1174  /*ScalarType std_norm_1 = viennacl::linalg::norm_1(std_C);
1175  ScalarType vcl_norm_1 = viennacl::linalg::norm_1(vcl_C);
1176  if ( std::fabs(std_norm_1 - vcl_norm_1) / std_norm_1 > epsilon)
1177  {
1178  std::cerr << "Failure at norm_1(): " << std::fabs(std_norm_1 - vcl_norm_1) / std_norm_1 << std::endl;
1179  return EXIT_FAILURE;
1180  }
1181 
1182  ScalarType std_norm_inf = boost::numeric::std::norm_inf(std_C);
1183  ScalarType vcl_norm_inf = viennacl::linalg::norm_inf(vcl_C);
1184  if ( std::fabs(std_norm_inf - vcl_norm_inf) / std_norm_inf > epsilon)
1185  {
1186  std::cerr << "Failure at norm_inf(): " << std::fabs(std_norm_inf - vcl_norm_inf) / std_norm_inf << std::endl;
1187  return EXIT_FAILURE;
1188  }*/
1189 
1190  ScalarType std_norm_frobenius = 0;
1191  for (std::size_t i=0; i<std_C.size(); ++i)
1192  for (std::size_t j=0; j<std_C[i].size(); ++j)
1193  std_norm_frobenius += std_C[i][j] * std_C[i][j];
1194  std_norm_frobenius = std::sqrt(std_norm_frobenius);
1195  ScalarType vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_C);
1196  if ( std::fabs(std_norm_frobenius - vcl_norm_frobenius) / std_norm_frobenius > epsilon)
1197  {
1198  std::cerr << "Failure at norm_frobenius()" << std::endl;
1199  return EXIT_FAILURE;
1200  }
1201 
1202  viennacl::scalar<ScalarType> device_std_norm_frobenius = std_norm_frobenius;
1203  viennacl::scalar<ScalarType> device_vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_C);
1204  if ( std::fabs(device_std_norm_frobenius - device_vcl_norm_frobenius) / device_std_norm_frobenius > epsilon)
1205  {
1206  std::cerr << "Failure at norm_frobenius()" << std::endl;
1207  return EXIT_FAILURE;
1208  }
1209 
1210  vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_range_C);
1211  if ( std::fabs(std_norm_frobenius - vcl_norm_frobenius) / std_norm_frobenius > epsilon)
1212  {
1213  std::cerr << "Failure at norm_frobenius() with range" << std::endl;
1214  return EXIT_FAILURE;
1215  }
1216 
1217  device_vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_range_C);
1218  if ( std::fabs(device_std_norm_frobenius - device_vcl_norm_frobenius) / device_std_norm_frobenius > epsilon)
1219  {
1220  std::cerr << "Failure at norm_frobenius() with range" << std::endl;
1221  return EXIT_FAILURE;
1222  }
1223 
1224  vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_slice_C);
1225  if ( std::fabs(std_norm_frobenius - vcl_norm_frobenius) / std_norm_frobenius > epsilon)
1226  {
1227  std::cerr << "Failure at norm_frobenius() with slice" << std::endl;
1228  return EXIT_FAILURE;
1229  }
1230 
1231  device_vcl_norm_frobenius = viennacl::linalg::norm_frobenius(vcl_slice_C);
1232  if ( std::fabs(device_std_norm_frobenius - device_vcl_norm_frobenius) / device_std_norm_frobenius > epsilon)
1233  {
1234  std::cerr << "Failure at norm_frobenius() with slice" << std::endl;
1235  return EXIT_FAILURE;
1236  }
1237 
1238  std::cout << "PASSED!" << std::endl;
1239 
1240 
1241  //
1242  // run operation tests:
1243  //
1244 
1245  std::cout << "//" << std::endl;
1246  std::cout << "////////// Test: Operations //////////" << std::endl;
1247  std::cout << "//" << std::endl;
1248 
1249 
1251  std::cout << "Testing A=matrix, B=matrix, C=matrix ..." << std::endl;
1252  viennacl::copy(std_A, vcl_A);
1253  viennacl::copy(std_B, vcl_B);
1254  viennacl::copy(std_C, vcl_C);
1255  if (run_test(epsilon,
1256  std_A, std_B, std_C,
1257  vcl_A, vcl_B, vcl_C) != EXIT_SUCCESS)
1258  {
1259  return EXIT_FAILURE;
1260  }
1261 
1262  std::cout << "Testing A=matrix, B=matrix, C=range ..." << std::endl;
1263  viennacl::copy(std_A, vcl_A);
1264  viennacl::copy(std_B, vcl_B);
1265  viennacl::copy(std_C, vcl_range_C);
1266  if (run_test(epsilon,
1267  std_A, std_B, std_C,
1268  vcl_A, vcl_B, vcl_range_C) != EXIT_SUCCESS)
1269  {
1270  return EXIT_FAILURE;
1271  }
1272 
1273  std::cout << "Testing A=matrix, B=matrix, C=slice ..." << std::endl;
1274  viennacl::copy(std_A, vcl_A);
1275  viennacl::copy(std_B, vcl_B);
1276  viennacl::copy(std_C, vcl_slice_C);
1277  if (run_test(epsilon,
1278  std_A, std_B, std_C,
1279  vcl_A, vcl_B, vcl_slice_C) != EXIT_SUCCESS)
1280  {
1281  return EXIT_FAILURE;
1282  }
1283 
1284  std::cout << "Testing A=matrix, B=range, C=matrix ..." << std::endl;
1285  viennacl::copy(std_A, vcl_A);
1286  viennacl::copy(std_B, vcl_range_B);
1287  viennacl::copy(std_C, vcl_C);
1288  if (run_test(epsilon,
1289  std_A, std_B, std_C,
1290  vcl_A, vcl_range_B, vcl_C) != EXIT_SUCCESS)
1291  {
1292  return EXIT_FAILURE;
1293  }
1294 
1295  std::cout << "Testing A=matrix, B=range, C=range ..." << std::endl;
1296  viennacl::copy(std_A, vcl_A);
1297  viennacl::copy(std_B, vcl_range_B);
1298  viennacl::copy(std_C, vcl_range_C);
1299  if (run_test(epsilon,
1300  std_A, std_B, std_C,
1301  vcl_A, vcl_range_B, vcl_range_C) != EXIT_SUCCESS)
1302  {
1303  return EXIT_FAILURE;
1304  }
1305 
1306  std::cout << "Testing A=matrix, B=range, C=slice ..." << std::endl;
1307  viennacl::copy(std_A, vcl_A);
1308  viennacl::copy(std_B, vcl_range_B);
1309  viennacl::copy(std_C, vcl_slice_C);
1310  if (run_test(epsilon,
1311  std_A, std_B, std_C,
1312  vcl_A, vcl_range_B, vcl_slice_C) != EXIT_SUCCESS)
1313  {
1314  return EXIT_FAILURE;
1315  }
1316 
1317 
1318  std::cout << "Testing A=matrix, B=slice, C=matrix ..." << std::endl;
1319  viennacl::copy(std_A, vcl_A);
1320  viennacl::copy(std_B, vcl_slice_B);
1321  viennacl::copy(std_C, vcl_C);
1322  if (run_test(epsilon,
1323  std_A, std_B, std_C,
1324  vcl_A, vcl_slice_B, vcl_C) != EXIT_SUCCESS)
1325  {
1326  return EXIT_FAILURE;
1327  }
1328 
1329  std::cout << "Testing A=matrix, B=slice, C=range ..." << std::endl;
1330  viennacl::copy(std_A, vcl_A);
1331  viennacl::copy(std_B, vcl_slice_B);
1332  viennacl::copy(std_C, vcl_range_C);
1333  if (run_test(epsilon,
1334  std_A, std_B, std_C,
1335  vcl_A, vcl_slice_B, vcl_range_C) != EXIT_SUCCESS)
1336  {
1337  return EXIT_FAILURE;
1338  }
1339 
1340  std::cout << "Testing A=matrix, B=slice, C=slice ..." << std::endl;
1341  viennacl::copy(std_A, vcl_A);
1342  viennacl::copy(std_B, vcl_slice_B);
1343  viennacl::copy(std_C, vcl_slice_C);
1344  if (run_test(epsilon,
1345  std_A, std_B, std_C,
1346  vcl_A, vcl_slice_B, vcl_slice_C) != EXIT_SUCCESS)
1347  {
1348  return EXIT_FAILURE;
1349  }
1350 
1351 
1352 
1354  std::cout << "Testing A=range, B=matrix, C=matrix ..." << std::endl;
1355  viennacl::copy(std_A, vcl_range_A);
1356  viennacl::copy(std_B, vcl_B);
1357  viennacl::copy(std_C, vcl_C);
1358  if (run_test(epsilon,
1359  std_A, std_B, std_C,
1360  vcl_range_A, vcl_B, vcl_C) != EXIT_SUCCESS)
1361  {
1362  return EXIT_FAILURE;
1363  }
1364 
1365  std::cout << "Testing A=range, B=matrix, C=range ..." << std::endl;
1366  viennacl::copy(std_A, vcl_range_A);
1367  viennacl::copy(std_B, vcl_B);
1368  viennacl::copy(std_C, vcl_range_C);
1369  if (run_test(epsilon,
1370  std_A, std_B, std_C,
1371  vcl_range_A, vcl_B, vcl_range_C) != EXIT_SUCCESS)
1372  {
1373  return EXIT_FAILURE;
1374  }
1375 
1376  std::cout << "Testing A=range, B=matrix, C=slice ..." << std::endl;
1377  viennacl::copy(std_A, vcl_range_A);
1378  viennacl::copy(std_B, vcl_B);
1379  viennacl::copy(std_C, vcl_slice_C);
1380  if (run_test(epsilon,
1381  std_A, std_B, std_C,
1382  vcl_range_A, vcl_B, vcl_slice_C) != EXIT_SUCCESS)
1383  {
1384  return EXIT_FAILURE;
1385  }
1386 
1387 
1388 
1389  std::cout << "Testing A=range, B=range, C=matrix ..." << std::endl;
1390  viennacl::copy(std_A, vcl_range_A);
1391  viennacl::copy(std_B, vcl_range_B);
1392  viennacl::copy(std_C, vcl_C);
1393  if (run_test(epsilon,
1394  std_A, std_B, std_C,
1395  vcl_range_A, vcl_range_B, vcl_C) != EXIT_SUCCESS)
1396  {
1397  return EXIT_FAILURE;
1398  }
1399 
1400  std::cout << "Testing A=range, B=range, C=range ..." << std::endl;
1401  viennacl::copy(std_A, vcl_range_A);
1402  viennacl::copy(std_B, vcl_range_B);
1403  viennacl::copy(std_C, vcl_range_C);
1404  if (run_test(epsilon,
1405  std_A, std_B, std_C,
1406  vcl_range_A, vcl_range_B, vcl_range_C) != EXIT_SUCCESS)
1407  {
1408  return EXIT_FAILURE;
1409  }
1410 
1411  std::cout << "Testing A=range, B=range, C=slice ..." << std::endl;
1412  viennacl::copy(std_A, vcl_range_A);
1413  viennacl::copy(std_B, vcl_range_B);
1414  viennacl::copy(std_C, vcl_slice_C);
1415  if (run_test(epsilon,
1416  std_A, std_B, std_C,
1417  vcl_range_A, vcl_range_B, vcl_slice_C) != EXIT_SUCCESS)
1418  {
1419  return EXIT_FAILURE;
1420  }
1421 
1422 
1423 
1424  std::cout << "Testing A=range, B=slice, C=matrix ..." << std::endl;
1425  viennacl::copy(std_A, vcl_range_A);
1426  viennacl::copy(std_B, vcl_slice_B);
1427  viennacl::copy(std_C, vcl_C);
1428  if (run_test(epsilon,
1429  std_A, std_B, std_C,
1430  vcl_range_A, vcl_slice_B, vcl_C) != EXIT_SUCCESS)
1431  {
1432  return EXIT_FAILURE;
1433  }
1434 
1435  std::cout << "Testing A=range, B=slice, C=range ..." << std::endl;
1436  viennacl::copy(std_A, vcl_range_A);
1437  viennacl::copy(std_B, vcl_slice_B);
1438  viennacl::copy(std_C, vcl_range_C);
1439  if (run_test(epsilon,
1440  std_A, std_B, std_C,
1441  vcl_range_A, vcl_slice_B, vcl_range_C) != EXIT_SUCCESS)
1442  {
1443  return EXIT_FAILURE;
1444  }
1445 
1446  std::cout << "Testing A=range, B=slice, C=slice ..." << std::endl;
1447  viennacl::copy(std_A, vcl_range_A);
1448  viennacl::copy(std_B, vcl_slice_B);
1449  viennacl::copy(std_C, vcl_slice_C);
1450  if (run_test(epsilon,
1451  std_A, std_B, std_C,
1452  vcl_range_A, vcl_slice_B, vcl_slice_C) != EXIT_SUCCESS)
1453  {
1454  return EXIT_FAILURE;
1455  }
1456 
1457 
1459  std::cout << "Testing A=slice, B=matrix, C=matrix ..." << std::endl;
1460  viennacl::copy(std_A, vcl_slice_A);
1461  viennacl::copy(std_B, vcl_B);
1462  viennacl::copy(std_C, vcl_C);
1463  if (run_test(epsilon,
1464  std_A, std_B, std_C,
1465  vcl_slice_A, vcl_B, vcl_C) != EXIT_SUCCESS)
1466  {
1467  return EXIT_FAILURE;
1468  }
1469 
1470  std::cout << "Testing A=slice, B=matrix, C=range ..." << std::endl;
1471  viennacl::copy(std_A, vcl_slice_A);
1472  viennacl::copy(std_B, vcl_B);
1473  viennacl::copy(std_C, vcl_range_C);
1474  if (run_test(epsilon,
1475  std_A, std_B, std_C,
1476  vcl_slice_A, vcl_B, vcl_range_C) != EXIT_SUCCESS)
1477  {
1478  return EXIT_FAILURE;
1479  }
1480 
1481  std::cout << "Testing A=slice, B=matrix, C=slice ..." << std::endl;
1482  viennacl::copy(std_A, vcl_slice_A);
1483  viennacl::copy(std_B, vcl_B);
1484  viennacl::copy(std_C, vcl_slice_C);
1485  if (run_test(epsilon,
1486  std_A, std_B, std_C,
1487  vcl_slice_A, vcl_B, vcl_slice_C) != EXIT_SUCCESS)
1488  {
1489  return EXIT_FAILURE;
1490  }
1491 
1492 
1493 
1494  std::cout << "Testing A=slice, B=range, C=matrix ..." << std::endl;
1495  viennacl::copy(std_A, vcl_slice_A);
1496  viennacl::copy(std_B, vcl_range_B);
1497  viennacl::copy(std_C, vcl_C);
1498  if (run_test(epsilon,
1499  std_A, std_B, std_C,
1500  vcl_slice_A, vcl_range_B, vcl_C) != EXIT_SUCCESS)
1501  {
1502  return EXIT_FAILURE;
1503  }
1504 
1505  std::cout << "Testing A=slice, B=range, C=range ..." << std::endl;
1506  viennacl::copy(std_A, vcl_slice_A);
1507  viennacl::copy(std_B, vcl_range_B);
1508  viennacl::copy(std_C, vcl_range_C);
1509  if (run_test(epsilon,
1510  std_A, std_B, std_C,
1511  vcl_slice_A, vcl_range_B, vcl_range_C) != EXIT_SUCCESS)
1512  {
1513  return EXIT_FAILURE;
1514  }
1515 
1516  std::cout << "Testing A=slice, B=range, C=slice ..." << std::endl;
1517  viennacl::copy(std_A, vcl_slice_A);
1518  viennacl::copy(std_B, vcl_range_B);
1519  viennacl::copy(std_C, vcl_slice_C);
1520  if (run_test(epsilon,
1521  std_A, std_B, std_C,
1522  vcl_slice_A, vcl_range_B, vcl_slice_C) != EXIT_SUCCESS)
1523  {
1524  return EXIT_FAILURE;
1525  }
1526 
1527 
1528 
1529  std::cout << "Testing A=slice, B=slice, C=matrix ..." << std::endl;
1530  viennacl::copy(std_A, vcl_slice_A);
1531  viennacl::copy(std_B, vcl_slice_B);
1532  viennacl::copy(std_C, vcl_C);
1533  if (run_test(epsilon,
1534  std_A, std_B, std_C,
1535  vcl_slice_A, vcl_slice_B, vcl_C) != EXIT_SUCCESS)
1536  {
1537  return EXIT_FAILURE;
1538  }
1539 
1540  std::cout << "Testing A=slice, B=slice, C=range ..." << std::endl;
1541  viennacl::copy(std_A, vcl_slice_A);
1542  viennacl::copy(std_B, vcl_slice_B);
1543  viennacl::copy(std_C, vcl_range_C);
1544  if (run_test(epsilon,
1545  std_A, std_B, std_C,
1546  vcl_slice_A, vcl_slice_B, vcl_range_C) != EXIT_SUCCESS)
1547  {
1548  return EXIT_FAILURE;
1549  }
1550 
1551  std::cout << "Testing A=slice, B=slice, C=slice ..." << std::endl;
1552  viennacl::copy(std_A, vcl_slice_A);
1553  viennacl::copy(std_B, vcl_slice_B);
1554  viennacl::copy(std_C, vcl_slice_C);
1555  if (run_test(epsilon,
1556  std_A, std_B, std_C,
1557  vcl_slice_A, vcl_slice_B, vcl_slice_C) != EXIT_SUCCESS)
1558  {
1559  return EXIT_FAILURE;
1560  }
1561 
1562 
1563  return EXIT_SUCCESS;
1564 }
1565 
1566 
viennacl::linalg::element_div
viennacl::vector_expression< const vector_base< T >, const vector_base< T >, op_element_binary< op_div > > element_div(vector_base< T > const &v1, vector_base< T > const &v2)
viennacl::scalar
This class represents a single scalar value on the GPU and behaves mostly like a built-in scalar type...
Definition: forwards.h:227
norm_2.hpp
Generic interface for the l^2-norm. See viennacl/linalg/vector_operations.hpp for implementations...
viennacl::trans
viennacl::enable_if< viennacl::is_any_sparse_matrix< M1 >::value, matrix_expression< const M1, const M1, op_trans > >::type trans(const M1 &mat)
Returns an expression template class representing a transposed matrix.
Definition: sparse_matrix_operations.hpp:376
viennacl::matrix_slice
Class for representing strided submatrices of a bigger matrix A.
Definition: forwards.h:443
prod.hpp
Generic interface for matrix-vector and matrix-matrix products. See viennacl/linalg/vector_operations...
matrix.hpp
Implementation of the dense matrix class.
viennacl::backend::finish
void finish()
Synchronizes the execution. finish() will only return after all compute kernels (CUDA, OpenCL) have completed.
Definition: memory.hpp:54
viennacl::matrix
A dense matrix class.
Definition: forwards.h:375
viennacl::identity_matrix
Represents a vector consisting of 1 at a given index and zeros otherwise. To be used as an initialize...
Definition: matrix_def.hpp:69
viennacl::linalg::detail::max
T max(const T &lhs, const T &rhs)
Maximum.
Definition: util.hpp:59
norm_1.hpp
Generic interface for the l^1-norm. See viennacl/linalg/vector_operations.hpp for implementations...
viennacl::traits::size
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:239
norm_frobenius.hpp
Generic interface for the Frobenius norm.
viennacl::scalar_matrix
Represents a vector consisting of scalars 's' only, i.e. v[i] = s for all i. To be used as an initial...
Definition: matrix_def.hpp:93
GENERATE_UNARY_OP_TEST
#define GENERATE_UNARY_OP_TEST(FUNCNAME)
viennacl::zero_matrix
Represents a vector consisting of zeros only. To be used as an initializer for viennacl::vector, vector_range, or vector_slize only.
Definition: matrix_def.hpp:81
viennacl::result_of::cpu_value_type::type
T::ERROR_CANNOT_DEDUCE_CPU_SCALAR_TYPE_FOR_T type
Definition: result_of.hpp:271
vector_proxy.hpp
Proxy classes for vectors.
matrix_proxy.hpp
Proxy classes for matrices.
check_for_equality
bool check_for_equality(std::vector< std::vector< NumericT > > const &std_A, VCLMatrixType const &vcl_A, double epsilon)
Definition: matrix_float_double.hpp:41
viennacl::linalg::norm_frobenius
scalar_expression< const matrix_base< NumericT >, const matrix_base< NumericT >, op_norm_frobenius > norm_frobenius(const matrix_base< NumericT > &A)
Definition: norm_frobenius.hpp:61
viennacl::copy
void copy(std::vector< NumericT > &cpu_vec, circulant_matrix< NumericT, AlignmentV > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
Definition: circulant_matrix.hpp:150
viennacl::basic_range
A range class that refers to an interval [start, stop), where 'start' is included, and 'stop' is excluded.
Definition: forwards.h:424
ScalarType
float ScalarType
Definition: fft_1d.cpp:42
viennacl::linalg::element_prod
viennacl::vector_expression< const vector_base< T >, const vector_base< T >, op_element_binary< op_prod > > element_prod(vector_base< T > const &v1, vector_base< T > const &v2)
viennacl::matrix_range
Class for representing non-strided submatrices of a bigger matrix A.
Definition: forwards.h:440
viennacl::basic_slice
A slice class that refers to an interval [start, stop), where 'start' is included, and 'stop' is excluded.
Definition: forwards.h:429
scalar.hpp
Implementation of the ViennaCL scalar class.
run_test
int run_test(double epsilon, STLMatrixType &std_A, STLMatrixType &std_B, STLMatrixType &std_C, ViennaCLMatrixType1 &vcl_A, ViennaCLMatrixType2 &vcl_B, ViennaCLMatrixType3 vcl_C)
Definition: matrix_float_double.hpp:72
norm_inf.hpp
Generic interface for the l^infty-norm. See viennacl/linalg/vector_operations.hpp for implementations...