lumped_mass_matrix_test.C

Go to the documentation of this file.

// libmesh includes
#include <libmesh/lumped_mass_matrix.h>
#include <libmesh/parallel.h>
#include <libmesh/id_types.h>
#include <libmesh/int_range.h>
#include <libmesh/numeric_vector.h>
 
#include "libmesh_cppunit.h"
#include "test_comm.h"
 
#include <memory>
 
#define UNUSED 0
 
using namespace libMesh;
 
class LumpedMassMatrixTest : public CppUnit::TestCase
{
public:
  LIBMESH_CPPUNIT_TEST_SUITE(LumpedMassMatrixTest);
 
  CPPUNIT_TEST(testNumerics);
 
  CPPUNIT_TEST_SUITE_END();
 
public:
  void setUp()
  {
    _comm = TestCommWorld;
    _matrix = std::make_unique<LumpedMassMatrix<Number>>(*_comm);
 
    numeric_index_type root_block_size = 2;
    _local_size = root_block_size + static_cast<numeric_index_type>(_comm->rank());
    _global_size = 0;
 
    for (processor_id_type p = 0; p < _comm->size(); ++p)
    {
      numeric_index_type block_size = root_block_size + static_cast<numeric_index_type>(p);
      _global_size += block_size;
    }
 
    _matrix->init(_global_size, UNUSED, _local_size, UNUSED);
  }
 
  void tearDown() {}
 
  void testNumerics()
  {
    LOG_UNIT_TEST;
 
    numeric_index_type beginning_index = _matrix->row_start();
    numeric_index_type end_index = _matrix->row_stop();
 
    CPPUNIT_ASSERT_EQUAL(_local_size, numeric_index_type(_matrix->row_stop() - _matrix->row_start()));
 
    _matrix->zero();
 
    std::vector<Real> gold_values(_local_size, 0);
 
    // Test add
    for (const auto i : make_range(beginning_index, end_index))
      for (const auto j : make_range(beginning_index, end_index))
      {
        const Real sgn = j % 2 ? 1 : -1;
        _matrix->add(i, j, sgn * j);
        gold_values[i - beginning_index] += j;
      }
 
    _matrix->close();
 
    for (const auto i : make_range(beginning_index, end_index))
      for (const auto j : make_range(beginning_index, end_index))
      {
        if (i != j)
          LIBMESH_ASSERT_NUMBERS_EQUAL(0, (*_matrix)(i, j), _tolerance);
        else
          LIBMESH_ASSERT_NUMBERS_EQUAL(gold_values[i - beginning_index], (*_matrix)(i, i), _tolerance);
      }
 
    // Test set
    for (const auto i : make_range(beginning_index, end_index))
    {
      const Real sgn = i % 2 ? 1 : -1;
      _matrix->set(i, i, sgn * gold_values[i - beginning_index]);
    }
 
    _matrix->close();
 
    for (const auto i : make_range(beginning_index, end_index))
      LIBMESH_ASSERT_NUMBERS_EQUAL(gold_values[i - beginning_index], (*_matrix)(i, i), _tolerance);
 
    // Test clone
    auto copy = _matrix->clone();
    // Check that matrices have the same local/global sizes
    CPPUNIT_ASSERT_EQUAL(copy->m(), _matrix->m());
    CPPUNIT_ASSERT_EQUAL(copy->n(), _matrix->n());
    CPPUNIT_ASSERT_EQUAL(copy->local_m(), _matrix->local_m());
    CPPUNIT_ASSERT_EQUAL(copy->row_start(), _matrix->row_start());
    CPPUNIT_ASSERT_EQUAL(copy->row_stop(), _matrix->row_stop());
 
    for (const auto i : make_range(beginning_index, end_index))
      for (const auto j : make_range(beginning_index, end_index))
        LIBMESH_ASSERT_NUMBERS_EQUAL((*_matrix)(i, j), (*copy)(i, j), _tolerance);
 
    // Test zero clone
    auto zero_copy = _matrix->zero_clone();
    // Check that matrices have the same local/global sizes
    CPPUNIT_ASSERT_EQUAL(zero_copy->m(), _matrix->m());
    CPPUNIT_ASSERT_EQUAL(zero_copy->n(), _matrix->n());
    CPPUNIT_ASSERT_EQUAL(zero_copy->local_m(), _matrix->local_m());
    CPPUNIT_ASSERT_EQUAL(zero_copy->row_start(), _matrix->row_start());
    CPPUNIT_ASSERT_EQUAL(zero_copy->row_stop(), _matrix->row_stop());
 
    for (const auto i : make_range(beginning_index, end_index))
      for (const auto j : make_range(beginning_index, end_index))
        LIBMESH_ASSERT_NUMBERS_EQUAL(0, (*zero_copy)(i, j), _tolerance);
  }
 
private:
 
  Parallel::Communicator * _comm;
  std::unique_ptr<LumpedMassMatrix<Number>> _matrix;
  numeric_index_type _local_size, _global_size;
  const Real _tolerance = TOLERANCE * TOLERANCE;
};
 
CPPUNIT_TEST_SUITE_REGISTRATION(LumpedMassMatrixTest);