psymm1

Source: tests/pdense/psymm1.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>
#include <myramath/dense/MatrixRange.h>

// Serial algorithms.
#include <myramath/utility/random.h>
#include <myramath/dense/transpose.h>
#include <myramath/dense/tril.h>
#include <myramath/dense/triu.h>
#include <myramath/dense/symm.h>
#include <myramath/dense/frobenius.h>

// Parallel algorithms.
#include <myramath/pdense/psymm.h>
#include <myramath/pdense/Options.h>

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

using namespace myra;
typedef pdense::Options Options;

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>();
  // Make upper/lower triangles of it.
  Matrix<Number> L = tril(A);
  Matrix<Number> U = triu(A);
  // Initialize options.
  auto options = Options::create().set_nthreads(4).set_blocksize(128);
  // Check psymm('L','L')
    {
    auto B  = Matrix<Number>::random(I,J);
    auto C1 = Matrix<Number>::random(I,J);
    auto C2 = C1;
    symm_inplace('L', 'L', C1, L, B, alpha, beta);
    psymm_inplace('L', 'L', C2, L, B, alpha, beta, options);
    Precision error = frobenius(C1-C2)/frobenius(C1);
    myra::out() << "  |symm('L','L')-psymm('L','L')| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  // Check psymm('L','U')
    {
    auto B  = Matrix<Number>::random(I,J);
    auto C1 = Matrix<Number>::random(I,J);
    auto C2 = C1;
    symm_inplace('L', 'U', C1, U, B, alpha, beta);
    psymm_inplace('L', 'U', C2, U, B, alpha, beta, options);
    Precision error = frobenius(C1-C2)/frobenius(C1);
    myra::out() << "  |symm('L','U')-psymm('L','U')| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  // Check psymm('R','L')
    {
    auto B  = Matrix<Number>::random(J,I);
    auto C1 = Matrix<Number>::random(J,I);
    auto C2 = C1;
    symm_inplace('R', 'L', C1, L, B, alpha, beta);
    psymm_inplace('R', 'L', C2, L, B, alpha, beta, options);
    Precision error = frobenius(C1-C2)/frobenius(C1);
    myra::out() << "  |symm('R','L')-psymm('R','L')| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  // Check psymm('R','U')
    {
    auto B  = Matrix<Number>::random(J,I);
    auto C1 = Matrix<Number>::random(J,I);
    auto C2 = C1;
    symm_inplace('R', 'U', C1, U, B, alpha, beta);
    psymm_inplace('R', 'U', C2, U, B, alpha, beta, options);
    Precision error = frobenius(C1-C2)/frobenius(C1);
    myra::out() << "  |symm('R','U')-psymm('R','U')| = " << error << std::endl;
    REQUIRE(error < tolerance);
    }
  }

} // namespace

ADD_TEST("psymm1","[pdense][parallel]")
  {
  int I = 512;
  int J = 256;
  test<NumberS> (I,J,1.0e-4f);
  test<NumberD> (I,J,1.0e-8);
  test<NumberC> (I,J,1.0e-4f);
  test<NumberZ> (I,J,1.0e-8);
  }

Results: [PASS]

float
  |symm('L','L')-psymm('L','L')| = 2.58582e-07
  |symm('L','U')-psymm('L','U')| = 3.31749e-07
  |symm('R','L')-psymm('R','L')| = 3.79184e-07
  |symm('R','U')-psymm('R','U')| = 3.78143e-07
double
  |symm('L','L')-psymm('L','L')| = 4.82599e-16
  |symm('L','U')-psymm('L','U')| = 6.09224e-16
  |symm('R','L')-psymm('R','L')| = 7.00946e-16
  |symm('R','U')-psymm('R','U')| = 6.97612e-16
std::complex<float>
  |symm('L','L')-psymm('L','L')| = 2.99332e-07
  |symm('L','U')-psymm('L','U')| = 3.5593e-07
  |symm('R','L')-psymm('R','L')| = 4.04806e-07
  |symm('R','U')-psymm('R','U')| = 4.05293e-07
std::complex<double>
  |symm('L','L')-psymm('L','L')| = 5.51966e-16
  |symm('L','U')-psymm('L','U')| = 6.58328e-16
  |symm('R','L')-psymm('R','L')| = 7.21684e-16
  |symm('R','U')-psymm('R','U')| = 7.19744e-16

Go back to Summary of /test programs.