rcholesky2_constructors

Source: tests/multifrontal/rcholesky/rcholesky2_constructors.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/Vector.h>
#include <myramath/dense/VectorRange.h>
#include <myramath/sparse/SparseMatrix.h>
#include <myramath/sparse/SparseMatrixRange.h>
#include <myramath/sparse/Permutation.h>

// Algorithms.
#include <myramath/dense/euclidean.h>
#include <myramath/sparse/gemv.h>
#include <myramath/sparse/laplacian2.h>

// Solver under test.
#include <myramath/multifrontal/rcholesky/SparseRCholeskySolver.h>

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

using namespace myra;

ADD_TEST("rcholesky2_constructors","[multifrontal][parallel]")
  {
  typedef double Precision;
  // Make laplacian matrix A.  
  int I = 64;
  int J = 64;
  int N = I*J;
  SparseMatrix<Precision> A = laplacian2<Precision>(I,J);
  // Define solver options.
  typedef SparseRCholeskySolver<Precision> Solver;
  typedef ::multifrontal::Options Options;
  Options options = Options::create().set_blocksize(16).set_globsize(4).set_nthreads(1); 
  // Three construction options ..
  Permutation perm = bisect2(I,J);
  Solver solver1(A,options);                // .. factor A, using internal reorder()'ing
  Solver solver2(A,perm,options);           // .. factor A, using caller-provided ordering
  Solver solver3(A,solver1.tree(),options); // .. factor A, reusing symbolic factorization
  // Make random b and three solution vectors.
  auto b = Vector<Precision>::random(N);
  auto x1 = b;
  auto x2 = b;
  auto x3 = b;
  // Solve system using each solver.
  solver1.solve(x1.column());
  solver2.solve(x2.column());
  solver3.solve(x3.column());
  // Check residuals.
  Precision error1 = euclidean(A*x1-b);
  Precision error2 = euclidean(A*x2-b);
  Precision error3 = euclidean(A*x3-b);
  myra::out() << "|A*x1-b| = " << error1 << std::endl;
  myra::out() << "|A*x2-b| = " << error2 << std::endl;
  myra::out() << "|A*x3-b| = " << error3 << std::endl;
  }
  

Results: [PASS]

|A*x1-b| = 1.02357e-14
|A*x2-b| = 1.01954e-14
|A*x3-b| = 1.02611e-14

Go back to Summary of /test programs.