MyraMath
csyrkL


Source: tests/dense/syrk1.cpp

1 // ========================================================================= //
2 // This file is part of MyraMath, copyright (c) 2014-2019 by Ryan A Chilton //
3 // and distributed by MyraCore, LLC. See LICENSE.txt for license terms. //
4 // ========================================================================= //
5 
11 // Containers.
13 #include <myramath/dense/Matrix.h>
14 
15 // Algorithms.
17 #include <myramath/dense/gemm.h>
18 #include <myramath/dense/syrk.h>
23 #include <myramath/dense/tril.h>
24 
25 // Reporting.
26 #include <tests/myratest.h>
27 
28 using namespace myra;
29 
30 namespace {
31 
32 template<class Number> void test(int I, int J, typename ReflectPrecision<Number>::type tolerance)
33  {
34  // Useful typedef.
35  typedef typename ReflectPrecision<Number>::type Precision;
36  myra::out() << typestring<Number>() << std::endl;
37  // Make random nonsquare matrix A.
38  auto A = Matrix<Number>::random(I,J);
39  Number alpha = random<Number>();
40  Number beta = random<Number>();
41  // Compare syrk_inplace('N') to gemm_inplace(A,'N',A,'T')
42  {
43  auto C1 = Matrix<Number>::random(I,I);
44  Matrix<Number> C2 = C1;
45  gemm_inplace(C1,A,'N',A,'T',alpha,beta);
46  syrk_inplace(C2,'L',A,'N',alpha,beta);
47  Precision error = frobenius(tril(C1-C2));
48  myra::out() << " |syrk('N')-gemm('N','T')| = " << error << std::endl;
49  REQUIRE(error < tolerance);
50  }
51  // Compare syrk_inplace('T') to gemm_inplace(A,'T',A,'N')
52  {
53  auto C1 = Matrix<Number>::random(J,J);
54  Matrix<Number> C2 = C1;
55  gemm_inplace(C1,A,'T',A,'N',alpha,beta);
56  syrk_inplace(C2,'L',A,'T',alpha,beta);
57  Precision error = frobenius(tril(C1-C2));
58  myra::out() << " |syrk('T')-gemm('T','N')| = " << error << std::endl;
59  REQUIRE(error < tolerance);
60  }
61  // Compare syrk_inplace('C') to gemm_inplace(A,'C',A,'H')
62  {
63  auto C1 = Matrix<Number>::random(I,I);
64  Matrix<Number> C2 = C1;
65  gemm_inplace(C1,A,'C',A,'H',alpha,beta);
66  syrk_inplace(C2,'L',A,'C',alpha,beta);
67  Precision error = frobenius(tril(C1-C2));
68  myra::out() << " |syrk('C')-gemm('C','H')| = " << error << std::endl;
69  REQUIRE(error < tolerance);
70  }
71  // Compare syrk_inplace('H') to gemm(A,'H',A,'C')
72  {
73  auto C1 = Matrix<Number>::random(J,J);
74  Matrix<Number> C2 = C1;
75  gemm_inplace(C1,A,'H',A,'C',alpha,beta);
76  syrk_inplace(C2,'L',A,'H',alpha,beta);
77  Precision error = frobenius(tril(C1-C2));
78  myra::out() << " |syrk('H')-gemm('H','C')| = " << error << std::endl;
79  REQUIRE(error < tolerance);
80  }
81  }
82 
83 } // namespace
84 
85 ADD_TEST("ssyrkL","[dense][blas]")
86  { test<NumberS>(57,24,1.0e-4f); }
87 
88 ADD_TEST("dsyrkL","[dense][blas]")
89  { test<NumberD>(57,24,1.0e-8); }
90 
91 ADD_TEST("csyrkL","[dense][blas]")
92  { test<NumberC>(57,24,1.0e-4f); }
93 
94 ADD_TEST("zhsykL","[dense][blas]")
95  { test<NumberZ>(57,24,1.0e-8); }
96 
Returns a conjugated copy of a Matrix or Vector. Or, conjugate one inplace.
Tabulates an IxJ matrix. Allows random access, has column major layout to be compatible with BLAS/LAP...
Definition: bdsqr.h:20
Routines for computing Frobenius norms of various algebraic containers.
static Matrix< Number > random(int I, int J)
Generates a random Matrix of specified size.
Definition: Matrix.cpp:353
Definition: syntax.dox:1
Returns a transposed copy of a Matrix. The inplace version only works on a square operand...
Returns the lower triangle of a dense Matrix.
Various utility functions/classes related to scalar Number types.
General purpose dense matrix container, O(i*j) storage.
Reflects Precision trait for a Number, scalar Number types should specialize it.
Definition: Number.h:33
Returns a hermitian copy of a Matrix. The inplace version only works on a square operand.
Simplistic random number functions.
Routines for symmetric rank-k updates, a specialized form of Matrix*Matrix multiplication.
Variety of routines all for dense Matrix*Matrix multiplies. Delegates to the BLAS.


Results: [PASS]

std::complex<float>
|syrk('N')-gemm('N','T')| = 0
|syrk('T')-gemm('T','N')| = 0
|syrk('C')-gemm('C','H')| = 0
|syrk('H')-gemm('H','C')| = 0


Go back to Summary of /test programs.