SpECTRE  v2024.09.29
Implementing SpECTRE vectors

Overview of SpECTRE Vectors

In SpECTRE, sets of contiguous or related data are stored in specializations of vector data types. The canonical implementation of this is the DataVector, which is used for storage of a contiguous sequence of doubles which support a wide variety of mathematical operations and represent data on a grid used during an evolution or elliptic solve. However, we support the ability to easily generate similar vector types which can hold data of a different type (e.g. std::complex<double>), or support a different set of mathematical operations. SpECTRE vector classes are derived from the class template VectorImpl. The remainder of this brief guide gives a description of the tools for defining additional vector types.

For reference, all functions described here can also be found in brief in the Doxygen documentation for VectorImpl.hpp, and a simple reference implementation can be found in DataVector.hpp and DataVector.cpp.

The class definition

SpECTRE vector types inherit from vector types implemented in the high-performance arithmetic library Blaze. Using inheritance, SpECTRE vectors gracefully make use of the math functions defined for the Blaze types, but can be customized for the specific needs in SpECTRE computations.

The trio of template parameters for VectorImpl are the type of the stored data (e.g. double for DataVector), the result type for mathematical operations, and the static size. The result type is used by Blaze to ensure that only compatible vector types are used together in mathematical expressions. For example, a vector representing double data on a grid (DataVector) cannot be added to a vector representing spectral coefficients (ModalVector). This avoids subtle bugs that arise when vector types are unintentionally mixed. In nearly all cases the result type will be the vector type that is being defined, so, for instance, DataVector is a derived class of VectorImpl<double, DataVector, 5>. This template pattern is known as the "Curiously Recurring Template Pattern" (CRTP). The static size is used as an optimization for small vector sizes. If your vector is small, rather than doing heap allocations, it will use stack allocations in order to be efficient. The default static size is set by a global constexpr bool default_vector_impl_static_size.

For the Blaze system to use the CRTP inheritance appropriately, it requires the specification of separate type traits in the blaze namespace. These traits can usually be declared in a standard form, so are abstracted in a macro. For any new vector MyNewVector, the pattern that must appear at the beginning of the file (i.e. before the class definition) is:

/// \cond
class MyNewVector;
/// \endcond
namespace blaze {
DECLARE_GENERAL_VECTOR_BLAZE_TRAITS(MyNewVector);
}

The class template VectorImpl defines various constructors, assignment operators, and iterator generation members. Most of these are inherited from Blaze types, but in addition, the methods set_data_ref, and pup are defined for use in SpECTRE. All except for the assignment operators and constructors will be implicitly inherited from VectorImpl. The assignment and constructors may be inherited calling the following alias code in the vector class definition:

using VectorImpl<T,VectorType,StaticSize>::operator=;
using VectorImpl<T,VectorType,StaticSize>::VectorImpl;

Only the mathematical operations supported on the base Blaze types are supported by default. Those operations are determined by the storage type T and by the Blaze library. See blaze-wiki/Vector_Operations.

Allowed operator specification

Blaze keeps track of the return type of unary and binary operations using "type trait" structs. These specializations for vector types should be placed in the header file associated with the VectorImpl specialization. For DataVector, the specializations are defined in DataStructures/DataVector.hpp. The presence or absence of template specializations of these structs also determines the set of allowed operations between the vector type and other types. For example, if adding a double to a DataVector should be allowed and the result should be treated as a DataVector for subsequent operations, then the struct blaze::AddTrait<DataVector, double> needs to be defined as follows:

namespace blaze {
// the `template <>` head tells the compiler that
// `AddTrait<DataVector, double>` is a class template specialization
template <>
struct AddTrait<DataVector, double> {
// the `Type` alias tells blaze that the result should be treated like a
// `DataVector` for any further operations
using Type = DataVector;
};
} // namespace blaze

Note that this only adds support for DataVector + double, not double + DataVector. To get the latter the following AddTrait specialization must be defined

namespace blaze {
// the `template <>` head tells the compiler that
// `AddTrait<double, DataVector>` is a class template specialization
template <>
struct AddTrait<double, DataVector> {
// the `Type` alias tells blaze that the result should be treated like a
// `DataVector` for any further operations
using Type = DataVector;
};
} // namespace blaze

Four helper macros are defined to assist with generating the many specializations that binary operations may require. Both of these macros must be put inside the blaze namespace for them to work correctly.

The first of these helper macros is BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, BLAZE_MATH_TRAIT), which will define all of the pairwise operations (BLAZE_MATH_TRAIT) for the vector type (VECTOR_TYPE) with itself and for the vector type with its value_type. This reduces the three specializations similar to the above code blocks to a single line call,

namespace blaze {
BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(DataVector, AddTrait)
} // namespace blaze

The second helper macro is provided to easily define all of the arithmetic operations that will typically be supported for a vector type with its value type. The macro is VECTOR_BLAZE_TRAIT_SPECIALIZE_ARITHMETIC_TRAITS(VECTOR_TYPE), and defines all of:

  • IsVector<VECTOR_TYPE> to std::true_type
  • TransposeFlag<VECTOR_TYPE>, which informs Blaze of the interpretation of the data as a "column" or "row" vector
  • AddTrait for the VECTOR_TYPE and its value type (3 Blaze struct specializations)
  • SubTrait for the VECTOR_TYPE and its value type (3 Blaze struct specializations)
  • MultTrait for the VECTOR_TYPE and its value type (3 Blaze struct specializations)
  • DivTrait for the VECTOR_TYPE and its value type (3 Blaze struct specializations)

This macro is similarly intended to be used in the blaze namespace and can substantially simplify these specializations for new vector types. For instance, the call for DataVector is:

namespace blaze {
VECTOR_BLAZE_TRAIT_SPECIALIZE_ARITHMETIC_TRAITS(DataVector)
} // namespace blaze

The third helper macro is provided to define a combination of Blaze traits for symmetric operations of a vector type with a second type (which may or may not be a vector type). The macro is BLAZE_TRAIT_SPECIALIZE_COMPATIBLE_BINARY_TRAIT(VECTOR, COMPATIBLE, TRAIT), and defines the appropriate trait for the two combinations <VECTOR, COMPATIBLE> and <COMPATIBLE, VECTOR>, and defines the result type to be VECTOR. For instance, to support the multiplication of a ComplexDataVector with a DataVector and have the result be a ComplexDataVector, the following macro call should be included in the blaze namespace:

namespace blaze {
BLAZE_TRAIT_SPECIALIZE_COMPATIBLE_BINARY_TRAIT(ComplexDataVector, DataVector,
MultTrait);
} // namespace blaze

Finally, the fourth helper macro is provided to define all of the blaze traits which are considered either unary or binary maps. This comprises most named unary functions (like sin() or sqrt()) and named binary functions (like hypot() and atan2()). The macro VECTOR_BLAZE_TRAIT_SPECIALIZE_ALL_MAP_TRAITS(VECTOR_TYPE) broadly specializes all blaze-defined maps in which the given VECTOR_TYPE as the sole argument (for unary maps) or both arguments (for binary maps). This macro is also intended to be used in the blaze namespace. The call for DataVector is:

namespace blaze {
VECTOR_BLAZE_TRAIT_SPECIALIZE_ALL_MAP_TRAITS(DataVector)
} // namespace blaze

Supporting operations for <tt>std::array</tt>s of vectors

In addition to operations between SpECTRE vectors, it is useful to gracefully handle operations between std::arrays of vectors element-wise. There are general macros defined for handling operations between array specializations: DEFINE_STD_ARRAY_BINOP and DEFINE_STD_ARRAY_INPLACE_BINOP from Utilities/StdArrayHelpers.hpp.

In addition, there is a macro for rapidly generating addition and subtraction between arrays of vectors and arrays of their data types. The macro MAKE_STD_ARRAY_VECTOR_BINOPS(VECTOR_TYPE) will define:

Equivalence operators

Equivalence operators are supported by the Blaze type inheritance. The equivalence operator == evaluates to true on a pair of vectors if they are the same size and contain the same values, regardless of ownership.

MakeWithValueImpl

SpECTRE offers the convenience function make_with_value for various types. The typical behavior for a SpECTRE vector type is to create a new vector type of the same type and length initialized with the value provided as the second argument in all entries. This behavior may be created by placing the macro MAKE_WITH_VALUE_IMPL_DEFINITION_FOR(VECTOR_TYPE) in the .hpp file. Any other specializations of MakeWithValueImpl will need to be written manually.

Interoperability with other data types

When additional vector types are added, small changes are necessary if they are to be used as the base container type either for Tensors or for Variables (a Variables contains Tensors), which contain some vector type.

In Tensor.hpp, there is a static_assert which white-lists the possible types that can be used as the storage type in Tensors. Any new vectors must be added to that white-list if they are to be used within Tensors.

Variables is templated on the storage type of the stored Tensors. However, any new data type should be appropriately tested. New vector types should be tested by invoking new versions of existing testing functions templated on the new vector type, rather than DataVector.

Writing tests

In addition to the utilities for generating new vector types, there are a number of convenience functions and utilities for easily generating the tests necessary to verify that the vectors function appropriately. These utilities are in VectorImplTestHelper.hpp, and documented individually in the TestingFrameworkGroup. Presented here are the salient details for rapidly assembling basic tests for vectors.

Utility check functions

Each of these functions is intended to encapsulate a single frequently used unit test and is templated (in order) on the vector type and the value type to be generated. The default behavior is to uniformly sample values between -100 and 100, but alternative bounds may be passed in via the function arguments.

<tt>TestHelpers::VectorImpl::vector_test_construct_and_assign()</tt>

This function tests a battery of construction and assignment operators for the vector type.

<tt>TestHelpers::VectorImpl::vector_test_serialize()</tt>

This function tests that vector types can be serialized and deserialized, retaining their data.

<tt>TestHelpers::VectorImpl::vector_test_ref()</tt>

This function tests the set_data_ref method of sharing data between vectors, and that the appropriate owning flags and move operations are handled correctly.

<tt>TestHelpers::VectorImpl::vector_test_math_after_move()</tt>

Tests several combinations of math operations and ownership before and after use of std::move.

<tt>TestHelpers::VectorImpl::vector_ref_test_size_error()</tt>

This function intentionally generates an error when assigning values from one vector to a differently sized, non-owning vector (made non-owning by use of set_data_ref). The assertion test which calls this function should search for the string "Must copy/move/assign into same size". Three forms of the test are provided, which are switched between using a value from the enum RefSizeErrorTestKind in the first function argument:

  • RefSizeErrorTestKind::Copy: tests that the size error is appropriately generated when copying to a non-owning vector of the wrong size. This has "copy" in the message.
  • RefSizeErrorTestKind::ExpressionAssign: tests that the size error is appropriately generated when assigning the result of a mathematical expression to a non-owning vector of the wrong size. This has "assign" in the message.
  • RefSizeErrorTestKind::Move: tests that the size error is appropriately generated when a vector is std::moved into a non-owning vector of the wrong size. This has "move" in the message.

<tt>TestHelpers::VectorImpl::test_functions_with_vector_arguments()</tt>

This is a general function for testing the mathematical operation of vector types with other vector types and/or their base types, with or without various reference wrappers. This may be used to efficiently test the full set of permitted math operations on a vector. See the documentation of test_functions_with_vector_arguments() for full usage details.

An example simple use case for the math test utility:

const TestHelpers::VectorImpl::Bound generic{{-100.0, 100.0}};
const TestHelpers::VectorImpl::Bound mone_one{{-1.0, 1.0}};
const TestHelpers::VectorImpl::Bound gt_one{{1.0, 100.0}};
const TestHelpers::VectorImpl::Bound positive{{0.01, 100.0}};
const auto unary_ops = std::make_tuple(
std::make_tuple(funcl::Abs<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Acos<>{}, std::make_tuple(mone_one)),
std::make_tuple(funcl::Acosh<>{}, std::make_tuple(gt_one)),
std::make_tuple(funcl::Asin<>{}, std::make_tuple(mone_one)),
std::make_tuple(funcl::Asinh<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Atan<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Atanh<>{}, std::make_tuple(mone_one)),
std::make_tuple(funcl::Cbrt<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Cos<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Cosh<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Erf<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Exp<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Exp2<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Fabs<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::InvCbrt<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::InvSqrt<>{}, std::make_tuple(positive)),
std::make_tuple(funcl::Log<>{}, std::make_tuple(positive)),
std::make_tuple(funcl::Log10<>{}, std::make_tuple(positive)),
std::make_tuple(funcl::Log2<>{}, std::make_tuple(positive)),
std::make_tuple(funcl::Sin<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Sinh<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::StepFunction<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Square<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Sqrt<>{}, std::make_tuple(positive)),
std::make_tuple(funcl::Tan<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::Tanh<>{}, std::make_tuple(generic)),
std::make_tuple(funcl::UnaryPow<1>{}, std::make_tuple(generic)),
std::make_tuple(funcl::UnaryPow<-2>{}, std::make_tuple(generic)),
std::make_tuple(funcl::UnaryPow<3>{}, std::make_tuple(generic)));
TestHelpers::VectorImpl::TestKind::Normal, DataVector>(unary_ops);
Stores a collection of function values.
Definition: DataVector.hpp:48
void test_functions_with_vector_arguments(const std::tuple< FunctionsAndArgumentBounds... > &tuple_of_functions_and_argument_bounds)
General entry function for testing arbitrary math functions on vector types.
Definition: VectorImplTestHelper.hpp:1280
Functional for computing abs on an object.
Definition: Functional.hpp:249
Functional for computing acos on an object.
Definition: Functional.hpp:250
Functional for computing acosh on an object.
Definition: Functional.hpp:251
Functional for computing asin on an object.
Definition: Functional.hpp:252
Functional for computing asinh on an object.
Definition: Functional.hpp:253
Functional for computing atan on an object.
Definition: Functional.hpp:254
Functional for computing atanh on an object.
Definition: Functional.hpp:255
Functional for computing cbrt on an object.
Definition: Functional.hpp:256
Functional for computing cos on an object.
Definition: Functional.hpp:258
Functional for computing cosh on an object.
Definition: Functional.hpp:259
Functional for computing erf on an object.
Definition: Functional.hpp:260
Functional for computing exp2 on an object.
Definition: Functional.hpp:262
Functional for computing exp on an object.
Definition: Functional.hpp:261
Functional for computing fabs on an object.
Definition: Functional.hpp:263
Functional for computing invcbrt on an object.
Definition: Functional.hpp:265
Functional for computing invsqrt on an object.
Definition: Functional.hpp:266
Functional for computing log10 on an object.
Definition: Functional.hpp:268
Functional for computing log2 on an object.
Definition: Functional.hpp:269
Functional for computing log on an object.
Definition: Functional.hpp:267
Functional for computing sin on an object.
Definition: Functional.hpp:271
Functional for computing sinh on an object.
Definition: Functional.hpp:272
Functional for computing sqrt on an object.
Definition: Functional.hpp:273
Function for squaring a quantity.
Definition: Functional.hpp:303
Functional for computing step_function on an object.
Definition: Functional.hpp:274
Functional for computing tan on an object.
Definition: Functional.hpp:275
Functional for computing tanh on an object.
Definition: Functional.hpp:276
Function for computing an integer power, forwards to template pow<N>()
Definition: Functional.hpp:281

More use cases of this functionality can be found in Test_DataVector.cpp.

Vector storage nuts and bolts

Internally, all vector classes inherit from the templated VectorImpl, which inherits from a blaze::CustomVector. Most of the mathematical operations are supported through the Blaze inheritance, which ensures that the math operations execute the optimized forms in Blaze.

Blaze also offers the possibility of restricting operations via groups in the blaze::CustomVector template arguments. Currently, we do not use the blaze::GroupTag functionality to determine available operations for vectors, but in principle this feature could allow us to further simplify our operator choice logic in the SpECTRE vector code.

SpECTRE vectors can be either "owning" or "non-owning". If a vector is owning, it allocates and controls the data it has access to, and is responsible for eventually freeing that data when the vector goes out of scope. If the vector is non-owning, it acts as a (possibly complete) "view" of otherwise allocated memory. Non-owning vectors do not manage memory, nor can they change size. The two cases of data ownership cause the underlying data to be handled fairly differently, so we will discuss each in turn.

When a SpECTRE vector is constructed as owning, or becomes owning, its memory is allocated in one of two ways.

  1. The size of the vector is larger than the StaticSize template parameter to VectorImpl. In that case, it allocates its own block of memory of appropriate size, and stores a pointer to that memory in a std::unique_ptr named owned_data_. The std::unique_ptr ensures that the SpECTRE vector needs to perform no further direct memory management, and that the memory will be appropriately managed whenever the std::unique_ptr owned_data_ member is deleted or moved.
  2. The size of the vector is less than or equal to the StaticSize template. In this case, the data is stored on the stack in a std::array<T, StaticSize> member variable called static_owned_data_. Since it is on the stack, this doesn't require any memory management by the user.

In either case, the base blaze::CustomVector must also be told about the pointer, which is always accomplished by calling the protected function VectorImpl.reset_pointer_vector(const size_t set_size), which sets the blaze::CustomVector internal pointer to either the pointer obtained by std::unique_pointer.get() or the pointer obtained by std::array.data() depending on the size of the vector.

When a SpECTRE vector is constructed as non-owning by the VectorImpl(ValueType* start, size_t set_size) constructor, or becomes non-owning by the set_data_ref function, neither the internal std::unique_ptr named owned_data_ nor the internal std::array named static_owned_data_ points to the data represented by the vector and both can be thought of as "inactive" for the purposes of computation and memory management. This behavior is desirable, because otherwise the std::unique_ptr would attempt to free memory that is presumed to be also used elsewhere, causing difficult to diagnose memory errors. And we needn't worry about the std::array because it's allocated on the stack. The non-owning SpECTRE vector updates the base blaze::CustomVector pointer directly by calling blaze::CustomVector.reset from the derived class (on itself).