zsymmU

Source: tests/dense/symm2.cpp

// ========================================================================= //
// This file is part of MyraMath, copyright (c) 2014-2019 by Ryan A Chilton  //
// and distributed by MyraCore, LLC. See LICENSE.txt for license terms.      //
// ========================================================================= //

// Containers.
#include <myramath/utility/Number.h>
#include <myramath/dense/Matrix.h>

// Algorithms.
#include <myramath/utility/random.h>
#include <myramath/dense/gemm.h>
#include <myramath/dense/symm.h>
#include <myramath/dense/transpose.h>
#include <myramath/dense/frobenius.h>
#include <myramath/dense/triu.h>

// Reporting.
#include <tests/myratest.h>

using namespace myra;

namespace {

template<class Number> void test(int I, int J, typename ReflectPrecision<Number>::type tolerance)
  {  
  typedef typename ReflectPrecision<Number>::type Precision;
  myra::out() << typestring<Number>() << std::endl;
  // Make random symmetric matrix A.
  auto A = Matrix<Number>::random(I,I);
  A = A + transpose(A);
  Number alpha = random<Number>();
  Number beta  = random<Number>();
  // Extract upper triangle, fill lower triangle with random data.
  Matrix<Number> U = triu(A);
  for (int j = 0; j < I; ++j)
    for (int i = j+1; i < I; ++i)
      U(i,j) = random<Number>();
  // Check symm('L','U')
    {
    auto B = Matrix<Number>::random(I,J);
    auto C1 = Matrix<Number>::random(I,J);
    Matrix<Number> C2 = C1;
    gemm_inplace(C1, A, 'N', B, 'N', alpha, beta);
    symm_inplace('L', 'U', C2, U, B, alpha, beta);
    Precision error = frobenius(C1-C2);
    myra::out() << "  |symm('L','U')-gemm()| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  // Check symm('R','U')
    {
    auto B = Matrix<Number>::random(J,I);
    auto C1 = Matrix<Number>::random(J,I);
    Matrix<Number> C2 = C1;
    gemm_inplace(C1, B, 'N', A, 'N', alpha, beta);
    symm_inplace('R', 'U', C2, U, B, alpha, beta);
    Precision error = frobenius(C1-C2);
    myra::out() << "  |symm('R','U')-gemm()| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  }

} // namespace

ADD_TEST("ssymmU","[dense][blas]")
  { test<NumberS>(51,32, 1.0e-4f); }

ADD_TEST("dsymmU","[dense][blas]")
  { test<NumberD>(51,32, 1.0e-8); }

ADD_TEST("csymmU","[dense][blas]")
  { test<NumberC>(51,32, 1.0e-4f); }

ADD_TEST("zsymmU","[dense][blas]")
  { test<NumberZ>(51,32, 1.0e-8); }

Results: [PASS]

std::complex<double>
  |symm('L','U')-gemm()| = 1.10954e-13
  |symm('R','U')-gemm()| = 1.3004e-13

Go back to Summary of /test programs.