26 #include <tests/myratest.h> 33 template<
class Precision>
void test(
int I,
int J, Precision tolerance)
35 typedef std::complex<Precision> Number;
36 myra::out() << typestring<Number>() << std::endl;
40 Number alpha = random<Number>();
41 Precision beta = random<Precision>();
43 auto options = Options::create().set_nthreads(4).set_blocksize(128);
47 for (
int i = 0; i < I; ++i)
48 C_serial(i,i) = std::real(C_serial(i,i));
49 auto C_parallel = C_serial;
50 her2k_inplace(C_serial,
'U',A,B,
'N',alpha,beta);
51 pher2k_inplace(C_parallel,
'U',A,B,
'N',alpha,beta,options);
52 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
53 myra::out() <<
" |her2k('U','N')-pher2k('U','N')| = " << error << std::endl;
54 REQUIRE(error < tolerance);
59 for (
int j = 0; j < J; ++j)
60 C_serial(j,j) = std::real(C_serial(j,j));
61 auto C_parallel = C_serial;
62 her2k_inplace(C_serial,
'U',A,B,
'T',alpha,beta);
63 pher2k_inplace(C_parallel,
'U',A,B,
'T',alpha,beta,options);
64 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
65 myra::out() <<
" |her2k('U','T')-pher2k('U','T')| = " << error << std::endl;
66 REQUIRE(error < tolerance);
71 for (
int j = 0; j < J; ++j)
72 C_serial(j,j) = std::real(C_serial(j,j));
73 auto C_parallel = C_serial;
74 her2k_inplace(C_serial,
'U',A,B,
'H',alpha,beta);
75 pher2k_inplace(C_parallel,
'U',A,B,
'H',alpha,beta,options);
76 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
77 myra::out() <<
" |her2k('U','H')-pher2k('U','H')| = " << error << std::endl;
78 REQUIRE(error < tolerance);
83 for (
int i = 0; i < I; ++i)
84 C_serial(i,i) = std::real(C_serial(i,i));
85 auto C_parallel = C_serial;
86 her2k_inplace(C_serial,
'U',A,B,
'C',alpha,beta);
87 pher2k_inplace(C_parallel,
'U',A,B,
'C',alpha,beta,options);
88 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
89 myra::out() <<
" |her2k('U','C')-pher2k('U','C')| = " << error << std::endl;
90 REQUIRE(error < tolerance);
95 for (
int i = 0; i < I; ++i)
96 C_serial(i,i) = std::real(C_serial(i,i));
97 auto C_parallel = C_serial;
98 her2k_inplace(C_serial,
'L',A,B,
'N',alpha,beta);
99 pher2k_inplace(C_parallel,
'L',A,B,
'N',alpha,beta);
100 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
101 myra::out() <<
" |her2k('L','N')-pher2k('L','N')| = " << error << std::endl;
102 REQUIRE(error < tolerance);
107 for (
int j = 0; j < J; ++j)
108 C_serial(j,j) = std::real(C_serial(j,j));
109 auto C_parallel = C_serial;
110 her2k_inplace(C_serial,
'L',A,B,
'T',alpha,beta);
111 pher2k_inplace(C_parallel,
'L',A,B,
'T',alpha,beta);
112 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
113 myra::out() <<
" |her2k('L','T')-pher2k('L','T')| = " << error << std::endl;
114 REQUIRE(error < tolerance);
119 for (
int j = 0; j < J; ++j)
120 C_serial(j,j) = std::real(C_serial(j,j));
121 auto C_parallel = C_serial;
122 her2k_inplace(C_serial,
'L',A,B,
'H',alpha,beta);
123 pher2k_inplace(C_parallel,
'L',A,B,
'H',alpha,beta);
124 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
125 myra::out() <<
" |her2k('L','H')-pher2k('L','H')| = " << error << std::endl;
126 REQUIRE(error < tolerance);
131 for (
int i = 0; i < I; ++i)
132 C_serial(i,i) = std::real(C_serial(i,i));
133 auto C_parallel = C_serial;
134 her2k_inplace(C_serial,
'L',A,B,
'C',alpha,beta);
135 pher2k_inplace(C_parallel,
'L',A,B,
'C',alpha,beta);
136 Precision error = frobenius(C_serial-C_parallel) / frobenius(C_serial);
137 myra::out() <<
" |her2k('L','C')-pher2k('L','C')| = " << error << std::endl;
138 REQUIRE(error < tolerance);
144 ADD_TEST(
"pher2k1",
"[pdense][parallel]")
148 test<float > (I,J,1.0e-4f);
149 test<double> (I,J,1.0e-8);
Thread parallel version of dense/her2k.h, hermitian rank-2k updates.
Interface class for representing subranges of dense Matrix's.
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
Options pack for routines in /pdense.
Definition: Options.h:24
Routines for hermitian rank-2k updates, a specialized form of Matrix*Matrix multiplication.
Various utility functions/classes related to scalar Number types.
General purpose dense matrix container, O(i*j) storage.
Options pack for routines in /pdense.
Simplistic random number functions.
std::complex<float>
|her2k('U','N')-pher2k('U','N')| = 0
|her2k('U','T')-pher2k('U','T')| = 0
|her2k('U','H')-pher2k('U','H')| = 0
|her2k('U','C')-pher2k('U','C')| = 0
|her2k('L','N')-pher2k('L','N')| = 3.78802e-08
|her2k('L','T')-pher2k('L','T')| = 0
|her2k('L','H')-pher2k('L','H')| = 0
|her2k('L','C')-pher2k('L','C')| = 3.76193e-08
std::complex<double>
|her2k('U','N')-pher2k('U','N')| = 1.11281e-16
|her2k('U','T')-pher2k('U','T')| = 1.16702e-16
|her2k('U','H')-pher2k('U','H')| = 1.17732e-16
|her2k('U','C')-pher2k('U','C')| = 1.11603e-16
|her2k('L','N')-pher2k('L','N')| = 1.08926e-16
|her2k('L','T')-pher2k('L','T')| = 0
|her2k('L','H')-pher2k('L','H')| = 0
|her2k('L','C')-pher2k('L','C')| = 1.08844e-16