MyraMath
cpinv


Source: tests/dense/pinv.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.
12 #include <myramath/dense/Matrix.h>
13 
14 // Algorithms.
15 #include <myramath/dense/pinv.h>
16 #include <myramath/dense/gemm.h>
19 
20 // Reporting.
21 #include <tests/myratest.h>
22 
23 using namespace myra;
24 
25 namespace {
26 
27 template<class Number> void test(int I, int J, int K, typename ReflectPrecision<Number>::type tolerance)
28  {
29  myra::out() << typestring<Number>() << std::endl;
30  typedef typename ReflectPrecision<Number>::type Precision;
31  // Construct random A with size IxJ and rank K
32  auto U = Matrix<Number>::random(I,K);
33  auto V = Matrix<Number>::random(K,J);
34  auto A = U*V;
35  // Compute B = pinv(A)
36  auto B = pinv(A);
37  // Examine the psuedoinverse properties:
38  Precision error1 = frobenius(A*B*A - A);
39  Precision error2 = frobenius(B*A*B - B);
40  Precision error3 = frobenius(A*B - hermitian(B)*hermitian(A));
41  Precision error4 = frobenius(B*A - hermitian(A)*hermitian(B));
42  myra::out() << " |A*B*A - A| = " << error1 << std::endl;
43  myra::out() << " |B*A*B - B| = " << error2 << std::endl;
44  myra::out() << " |A*B - B'A'| = " << error3 << std::endl;
45  myra::out() << " |B*A - A'B'| = " << error4 << std::endl;
46  // Assert everything.
47  REQUIRE(error1 < tolerance);
48  REQUIRE(error2 < tolerance);
49  REQUIRE(error3 < tolerance);
50  REQUIRE(error4 < tolerance);
51  }
52 
53 } // namespace
54 
55 ADD_TEST("spinv","[dense]")
56  { test<NumberS>(10,12,4,1.0e-4f); }
57 
58 ADD_TEST("dpinv","[dense]")
59  { test<NumberD>(10,12,4,1.0e-10); }
60 
61 ADD_TEST("cpinv","[dense]")
62  { test<NumberC>(10,12,4,1.0e-4f); }
63 
64 ADD_TEST("zpinv","[dense]")
65  { test<NumberZ>(10,12,4,1.0e-10); }
66 
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 Moore-Penrose pseudoinverse of A.
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.
Variety of routines all for dense Matrix*Matrix multiplies. Delegates to the BLAS.


Results: [PASS]

std::complex<float>
|A*B*A - A| = 3.16424e-06
|B*A*B - B| = 1.05047e-07
|A*B - B'A'| = 9.02602e-07
|B*A - A'B'| = 2.19357e-06


Go back to Summary of /test programs.