petsc_vector_test.C

Go to the documentation of this file.

#include <libmesh/petsc_vector.h>
 
#ifdef LIBMESH_HAVE_PETSC
 
#include "numeric_vector_test.h"
 
 
using namespace libMesh;
 
class PetscVectorTest : public NumericVectorTest<PetscVector<Number>> {
public:
  PetscVectorTest() :
    NumericVectorTest<PetscVector<Number>>() {
    if (unitlog->summarized_logs_enabled())
      this->libmesh_suite_name = "NumericVectorTest";
    else
      this->libmesh_suite_name = "PetscVectorTest";
  }
 
  CPPUNIT_TEST_SUITE( PetscVectorTest );
 
  NUMERICVECTORTEST
 
  CPPUNIT_TEST( testGetArray );
  CPPUNIT_TEST( testPetscOperations );
  CPPUNIT_TEST( testSubvectors );
  CPPUNIT_TEST( testSubvectorsBase );
 
  CPPUNIT_TEST_SUITE_END();
 
  void testGetArray()
  {
    LOG_UNIT_TEST;
 
    unsigned int min_block_size  = 2;
 
    // a different size on each processor.
    unsigned int my_p = my_comm->rank();
    unsigned int local_size  = (min_block_size + my_p);
    unsigned int global_size = 0;
    unsigned int my_offset = 0;
 
    for (libMesh::processor_id_type p=0; p<my_comm->size(); p++)
      {
        const unsigned int p_size =
          (min_block_size + static_cast<unsigned int>(p));
        global_size += p_size;
        if (p < my_p)
          my_offset += p_size;
      }
 
    PetscVector<Number> v(*my_comm, global_size, local_size);
 
    PetscScalar * values = v.get_array();
 
    for (unsigned int i=0; i<local_size; i++)
      values[i] = i;
 
    v.restore_array();
 
    v.close();
 
    // Check the values through the interface
    for (unsigned int i=0; i<local_size; i++)
      LIBMESH_ASSERT_FP_EQUAL(i, std::abs(v(my_offset + i)), TOLERANCE*TOLERANCE);
 
    // Check that we can see the same thing with get_array_read
    const PetscScalar * read_only_values = v.get_array_read();
 
    for (unsigned int i=0; i<local_size; i++)
      LIBMESH_ASSERT_FP_EQUAL(i, std::abs(read_only_values[i]), TOLERANCE*TOLERANCE);
 
    v.restore_array();
 
    // Test getting a read only array after getting a writable array
    values = v.get_array();
    read_only_values = v.get_array_read();
    CPPUNIT_ASSERT_EQUAL(read_only_values, const_cast<const PetscScalar *>(values));
 
    v.restore_array();
 
    // Test to make sure we can get arrays after other operators
    for (unsigned int i = 0; i < local_size; i++)
      v.set(my_offset + i, i * 2.0);
    v.close();
    for (unsigned int i = 0; i < local_size; i++)
      LIBMESH_ASSERT_FP_EQUAL(i * 2.0, std::abs(v(my_offset + i)), TOLERANCE * TOLERANCE);
    values = v.get_array();
    read_only_values = v.get_array_read();
    for (unsigned int i = 0; i < local_size; i++)
    {
      LIBMESH_ASSERT_FP_EQUAL(i * 2.0, std::abs(values[i]), TOLERANCE * TOLERANCE);
      LIBMESH_ASSERT_FP_EQUAL(i * 2.0, std::abs(read_only_values[i]), TOLERANCE * TOLERANCE);
    }
    v.restore_array();
  }
 
  void testPetscOperations()
  {
    PetscVector<Number> v1(*my_comm, global_size, local_size);
    auto v2 = v1.clone();
 
    const libMesh::dof_id_type
      first = v1.first_local_index(),
      last  = v1.last_local_index();
 
    for (libMesh::dof_id_type n=first; n != last; n++)
    {
      v1.set (n, static_cast<libMesh::Number>(n+1));
      v2->set (n, static_cast<libMesh::Number>(2*(n+1)));
    }
    v1.close();
    v2->close();
 
    auto v_working_ptr = v1.clone();
    auto & v_working = *v_working_ptr;
 
    v_working.pointwise_mult(v1, *v2);
 
    for (libMesh::dof_id_type n=first; n != last; n++)
      LIBMESH_ASSERT_NUMBERS_EQUAL
        (v_working(n), libMesh::Real((n+1)*2*(n+1)),
         libMesh::TOLERANCE*libMesh::TOLERANCE);
 
    v_working.pointwise_divide(v1, *v2);
 
    for (libMesh::dof_id_type n=first; n != last; n++)
      LIBMESH_ASSERT_NUMBERS_EQUAL
        (v_working(n), libMesh::Real(0.5),
         libMesh::TOLERANCE*libMesh::TOLERANCE);
  }
 
};
 
CPPUNIT_TEST_SUITE_REGISTRATION( PetscVectorTest );
 
#endif // #ifdef LIBMESH_HAVE_PETSC