aboutsummaryrefslogtreecommitdiffhomepage
diff options
context:
space:
mode:
-rw-r--r--Eigen/src/Core/Matrix.h5
-rw-r--r--Eigen/src/Core/util/Memory.h2
-rw-r--r--Eigen/src/Geometry/AlignedBox.h2
-rw-r--r--Eigen/src/Geometry/Hyperplane.h2
-rw-r--r--Eigen/src/Geometry/ParametrizedLine.h2
-rw-r--r--Eigen/src/Geometry/Quaternion.h21
-rw-r--r--Eigen/src/Geometry/Scaling.h2
-rw-r--r--Eigen/src/Geometry/Transform.h5
-rw-r--r--Eigen/src/Geometry/Translation.h2
-rw-r--r--Eigen/src/StdVector/StdVector.h46
-rw-r--r--Eigen/src/StdVector/UnalignedType.h34
-rw-r--r--test/geometry.cpp8
-rw-r--r--test/stdvector.cpp105
13 files changed, 192 insertions, 44 deletions
diff --git a/Eigen/src/Core/Matrix.h b/Eigen/src/Core/Matrix.h
index 5349029bf..75c79d2fd 100644
--- a/Eigen/src/Core/Matrix.h
+++ b/Eigen/src/Core/Matrix.h
@@ -132,8 +132,6 @@ class Matrix
protected:
ei_matrix_storage<Scalar, MaxSizeAtCompileTime, RowsAtCompileTime, ColsAtCompileTime, Options> m_storage;
- Matrix(ei_constructor_without_unaligned_array_assert)
- : m_storage(ei_constructor_without_unaligned_array_assert()) {}
public:
enum { NeedsToAlign = (Options&AutoAlign) == AutoAlign
@@ -339,6 +337,9 @@ class Matrix
ei_assert(RowsAtCompileTime > 0 && ColsAtCompileTime > 0);
}
+ Matrix(ei_constructor_without_unaligned_array_assert)
+ : m_storage(ei_constructor_without_unaligned_array_assert()) {}
+
/** Constructs a vector or row-vector with given dimension. \only_for_vectors
*
* Note that this is only useful for dynamic-size vectors. For fixed-size vectors,
diff --git a/Eigen/src/Core/util/Memory.h b/Eigen/src/Core/util/Memory.h
index 21d9491d3..67a023bb0 100644
--- a/Eigen/src/Core/util/Memory.h
+++ b/Eigen/src/Core/util/Memory.h
@@ -303,7 +303,7 @@ inline static int ei_alignmentOffset(const Scalar* ptr, int maxOffset)
EIGEN_WORKAROUND_FOR_QT_BUG_CALLING_WRONG_OPERATOR_NEW
#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF(true)
-#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(Scalar,Size) \
+#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(Scalar,Size) \
EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF(((Size)!=Eigen::Dynamic) && ((sizeof(Scalar)*(Size))%16==0))
/** Deprecated, use the EIGEN_MAKE_ALIGNED_OPERATOR_NEW macro instead in your own class */
diff --git a/Eigen/src/Geometry/AlignedBox.h b/Eigen/src/Geometry/AlignedBox.h
index 1c2635cd5..652e72768 100644
--- a/Eigen/src/Geometry/AlignedBox.h
+++ b/Eigen/src/Geometry/AlignedBox.h
@@ -41,7 +41,7 @@ template <typename _Scalar, int _AmbientDim>
class AlignedBox
{
public:
-EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_AmbientDim==Dynamic ? Dynamic : _AmbientDim+1)
+EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_AmbientDim==Dynamic ? Dynamic : _AmbientDim+1)
enum { AmbientDimAtCompileTime = _AmbientDim };
typedef _Scalar Scalar;
typedef typename NumTraits<Scalar>::Real RealScalar;
diff --git a/Eigen/src/Geometry/Hyperplane.h b/Eigen/src/Geometry/Hyperplane.h
index a3425f6cb..44ffad36d 100644
--- a/Eigen/src/Geometry/Hyperplane.h
+++ b/Eigen/src/Geometry/Hyperplane.h
@@ -47,7 +47,7 @@ template <typename _Scalar, int _AmbientDim>
class Hyperplane
{
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_AmbientDim==Dynamic ? Dynamic : _AmbientDim+1)
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_AmbientDim==Dynamic ? Dynamic : _AmbientDim+1)
enum { AmbientDimAtCompileTime = _AmbientDim };
typedef _Scalar Scalar;
typedef typename NumTraits<Scalar>::Real RealScalar;
diff --git a/Eigen/src/Geometry/ParametrizedLine.h b/Eigen/src/Geometry/ParametrizedLine.h
index da30c8e82..7897382df 100644
--- a/Eigen/src/Geometry/ParametrizedLine.h
+++ b/Eigen/src/Geometry/ParametrizedLine.h
@@ -43,7 +43,7 @@ template <typename _Scalar, int _AmbientDim>
class ParametrizedLine
{
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_AmbientDim)
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_AmbientDim)
enum { AmbientDimAtCompileTime = _AmbientDim };
typedef _Scalar Scalar;
typedef typename NumTraits<Scalar>::Real RealScalar;
diff --git a/Eigen/src/Geometry/Quaternion.h b/Eigen/src/Geometry/Quaternion.h
index f595adcb9..d92ed48c9 100644
--- a/Eigen/src/Geometry/Quaternion.h
+++ b/Eigen/src/Geometry/Quaternion.h
@@ -61,17 +61,17 @@ template<typename _Scalar>
class Quaternion : public RotationBase<Quaternion<_Scalar>,3>
{
typedef RotationBase<Quaternion<_Scalar>,3> Base;
- typedef Matrix<_Scalar, 4, 1> Coefficients;
- Coefficients m_coeffs;
-
+
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,4)
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,4)
using Base::operator*;
-
+
/** the scalar type of the coefficients */
typedef _Scalar Scalar;
+ /** the type of the Coefficients 4-vector */
+ typedef Matrix<Scalar, 4, 1> Coefficients;
/** the type of a 3D vector */
typedef Matrix<Scalar,3,1> Vector3;
/** the equivalent rotation matrix type */
@@ -110,8 +110,11 @@ public:
inline Coefficients& coeffs() { return m_coeffs; }
/** Default constructor and initializing an identity quaternion. */
- inline Quaternion()
- { m_coeffs << 0, 0, 0, 1; }
+ inline Quaternion() {}
+
+ inline Quaternion(ei_constructor_without_unaligned_array_assert)
+ : m_coeffs(ei_constructor_without_unaligned_array_assert()) {}
+
/** Constructs and initializes the quaternion \f$ w+xi+yj+zk \f$ from
* its four coefficients \a w, \a x, \a y and \a z.
@@ -149,7 +152,7 @@ public:
/** \sa Quaternion::Identity(), MatrixBase::setIdentity()
*/
- inline Quaternion& setIdentity() { m_coeffs << 1, 0, 0, 0; return *this; }
+ inline Quaternion& setIdentity() { m_coeffs << 0, 0, 0, 1; return *this; }
/** \returns the squared norm of the quaternion's coefficients
* \sa Quaternion::norm(), MatrixBase::squaredNorm()
@@ -214,6 +217,8 @@ public:
bool isApprox(const Quaternion& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const
{ return m_coeffs.isApprox(other.m_coeffs, prec); }
+protected:
+ Coefficients m_coeffs;
};
/** \ingroup GeometryModule
diff --git a/Eigen/src/Geometry/Scaling.h b/Eigen/src/Geometry/Scaling.h
index d46296707..186cf5b33 100644
--- a/Eigen/src/Geometry/Scaling.h
+++ b/Eigen/src/Geometry/Scaling.h
@@ -43,7 +43,7 @@ template<typename _Scalar, int _Dim>
class Scaling
{
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_Dim)
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_Dim)
/** dimension of the space */
enum { Dim = _Dim };
/** the scalar type of the coefficients */
diff --git a/Eigen/src/Geometry/Transform.h b/Eigen/src/Geometry/Transform.h
index cb092985f..83bad3afd 100644
--- a/Eigen/src/Geometry/Transform.h
+++ b/Eigen/src/Geometry/Transform.h
@@ -63,7 +63,7 @@ template<typename _Scalar, int _Dim>
class Transform
{
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_Dim==Dynamic ? Dynamic : (_Dim+1)*(_Dim+1))
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_Dim==Dynamic ? Dynamic : (_Dim+1)*(_Dim+1))
enum {
Dim = _Dim, ///< space dimension in which the transformation holds
HDim = _Dim+1 ///< size of a respective homogeneous vector
@@ -94,6 +94,9 @@ public:
/** Default constructor without initialization of the coefficients. */
inline Transform() { }
+ inline Transform(ei_constructor_without_unaligned_array_assert)
+ : m_matrix(ei_constructor_without_unaligned_array_assert()) {}
+
inline Transform(const Transform& other)
{
m_matrix = other.m_matrix;
diff --git a/Eigen/src/Geometry/Translation.h b/Eigen/src/Geometry/Translation.h
index ab8ce3899..ba8f728c3 100644
--- a/Eigen/src/Geometry/Translation.h
+++ b/Eigen/src/Geometry/Translation.h
@@ -43,7 +43,7 @@ template<typename _Scalar, int _Dim>
class Translation
{
public:
- EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE(_Scalar,_Dim)
+ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_Dim)
/** dimension of the space */
enum { Dim = _Dim };
/** the scalar type of the coefficients */
diff --git a/Eigen/src/StdVector/StdVector.h b/Eigen/src/StdVector/StdVector.h
index 432582c52..a45036736 100644
--- a/Eigen/src/StdVector/StdVector.h
+++ b/Eigen/src/StdVector/StdVector.h
@@ -26,6 +26,20 @@
#ifndef EIGEN_STDVECTOR_H
#define EIGEN_STDVECTOR_H
+#define EIGEN_STD_VECTOR_SPECIALIZATION_BODY \
+ typedef Eigen::aligned_allocator<value_type> allocator_type; \
+ typedef vector<value_type, allocator_type > unaligned_base; \
+ typedef typename unaligned_base::size_type size_type; \
+ typedef typename unaligned_base::iterator iterator; \
+ explicit vector(const allocator_type& __a = allocator_type()) : unaligned_base(__a) {} \
+ vector(const vector& c) : unaligned_base(c) {} \
+ vector(size_type num, const value_type& val = value_type()) : unaligned_base(num, val) {}\
+ vector(iterator start, iterator end) : unaligned_base(start, end) {} \
+ vector& operator=(const vector& __x) { \
+ unaligned_base::operator=(__x); \
+ return *this; \
+ }
+
template <typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols, typename _Alloc>
class vector<Eigen::Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols>, _Alloc>
: public vector<Eigen::ei_unaligned_type<Eigen::Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols> >,
@@ -33,19 +47,27 @@ class vector<Eigen::Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols>, _All
{
public:
typedef Eigen::ei_unaligned_type<Eigen::Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols> > value_type;
- typedef Eigen::aligned_allocator<value_type> allocator_type;
- typedef vector<value_type, allocator_type > unaligned_base;
- typedef typename unaligned_base::size_type size_type;
- typedef typename unaligned_base::iterator iterator;
+ EIGEN_STD_VECTOR_SPECIALIZATION_BODY
+};
- explicit vector(const allocator_type& __a = allocator_type()) : unaligned_base(__a) {}
- vector(const vector& c) : unaligned_base(c) {}
- vector(size_type num, const value_type& val = value_type()) : unaligned_base(num, val) {}
- vector(iterator start, iterator end) : unaligned_base(start, end) {}
- vector& operator=(const vector& __x) {
- unaligned_base::operator=(__x);
- return *this;
- }
+template <typename _Scalar, int _Dim, typename _Alloc>
+class vector<Eigen::Transform<_Scalar,_Dim>, _Alloc>
+ : public vector<Eigen::ei_unaligned_type<Eigen::Transform<_Scalar,_Dim> >,
+ Eigen::aligned_allocator<Eigen::ei_unaligned_type<Eigen::Transform<_Scalar,_Dim> > > >
+{
+public:
+ typedef Eigen::ei_unaligned_type<Eigen::Transform<_Scalar,_Dim> > value_type;
+ EIGEN_STD_VECTOR_SPECIALIZATION_BODY
+};
+
+template <typename _Scalar, typename _Alloc>
+class vector<Eigen::Quaternion<_Scalar>, _Alloc>
+ : public vector<Eigen::ei_unaligned_type<Eigen::Quaternion<_Scalar> >,
+ Eigen::aligned_allocator<Eigen::ei_unaligned_type<Eigen::Quaternion<_Scalar> > > >
+{
+public:
+ typedef Eigen::ei_unaligned_type<Eigen::Quaternion<_Scalar> > value_type;
+ EIGEN_STD_VECTOR_SPECIALIZATION_BODY
};
#endif // EIGEN_STDVECTOR_H
diff --git a/Eigen/src/StdVector/UnalignedType.h b/Eigen/src/StdVector/UnalignedType.h
index 73cb5e9a0..728235ac1 100644
--- a/Eigen/src/StdVector/UnalignedType.h
+++ b/Eigen/src/StdVector/UnalignedType.h
@@ -35,12 +35,42 @@ class ei_unaligned_type<Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols> >
public:
typedef Matrix<_Scalar,_Rows,_Cols,_Options,_MaxRows,_MaxCols> aligned_base;
ei_unaligned_type() : aligned_base(ei_constructor_without_unaligned_array_assert()) {}
- ei_unaligned_type(const aligned_base& other)
- : aligned_base(ei_constructor_without_unaligned_array_assert())
+ ei_unaligned_type(const aligned_base& other) : aligned_base(ei_constructor_without_unaligned_array_assert())
{
resize(other.rows(), other.cols());
ei_assign_impl<ei_unaligned_type,aligned_base,NoVectorization>::run(*this, other);
}
};
+template<typename _Scalar, int _Dim>
+class ei_unaligned_type<Transform<_Scalar,_Dim> >
+ : public Transform<_Scalar,_Dim>
+{
+ public:
+ typedef Transform<_Scalar,_Dim> aligned_base;
+ typedef typename aligned_base::MatrixType MatrixType;
+ ei_unaligned_type() : aligned_base(ei_constructor_without_unaligned_array_assert()) {}
+ ei_unaligned_type(const aligned_base& other) : aligned_base(ei_constructor_without_unaligned_array_assert())
+ {
+ // no resizing here, it's fixed-size anyway
+ ei_assign_impl<MatrixType,MatrixType,NoVectorization>::run(this->matrix(), other.matrix());
+ }
+};
+
+template<typename _Scalar>
+class ei_unaligned_type<Quaternion<_Scalar> >
+ : public Quaternion<_Scalar>
+{
+ public:
+ typedef Quaternion<_Scalar> aligned_base;
+ typedef typename aligned_base::Coefficients Coefficients;
+ ei_unaligned_type() : aligned_base(ei_constructor_without_unaligned_array_assert()) {}
+ ei_unaligned_type(const aligned_base& other) : aligned_base(ei_constructor_without_unaligned_array_assert())
+ {
+ // no resizing here, it's fixed-size anyway
+ ei_assign_impl<Coefficients,Coefficients,NoVectorization>::run(this->coeffs(), other.coeffs());
+ }
+};
+
+
#endif // EIGEN_UNALIGNEDTYPE_H
diff --git a/test/geometry.cpp b/test/geometry.cpp
index f092d472a..cd2dbc7ab 100644
--- a/test/geometry.cpp
+++ b/test/geometry.cpp
@@ -52,7 +52,6 @@ template<typename Scalar> void geometry(void)
if (ei_is_same_type<Scalar,float>::ret)
largeEps = 1e-3f;
- Quaternionx q1, q2;
Vector3 v0 = Vector3::Random(),
v1 = Vector3::Random(),
v2 = Vector3::Random();
@@ -69,6 +68,13 @@ template<typename Scalar> void geometry(void)
(v0.cross(v1).cross(v0)).normalized();
VERIFY(m.isUnitary());
+ // Quaternion: Identity(), setIdentity();
+ Quaternionx q1, q2;
+ q2.setIdentity();
+ VERIFY_IS_APPROX(Quaternionx(Quaternionx::Identity()).coeffs(), q2.coeffs());
+ q1.coeffs().setRandom();
+ VERIFY_IS_APPROX(q1.coeffs(), (q1*q2).coeffs());
+
// unitOrthogonal
VERIFY_IS_MUCH_SMALLER_THAN(u0.unitOrthogonal().dot(u0), Scalar(1));
VERIFY_IS_MUCH_SMALLER_THAN(v0.unitOrthogonal().dot(v0), Scalar(1));
diff --git a/test/stdvector.cpp b/test/stdvector.cpp
index 024f97287..ac7a2967e 100644
--- a/test/stdvector.cpp
+++ b/test/stdvector.cpp
@@ -24,9 +24,10 @@
#include "main.h"
#include <Eigen/StdVector>
+#include <Eigen/Geometry>
template<typename MatrixType>
-void check_stdvector(const MatrixType& m)
+void check_stdvector_matrix(const MatrixType& m)
{
int rows = m.rows();
int cols = m.cols();
@@ -61,22 +62,102 @@ void check_stdvector(const MatrixType& m)
}
}
+template<typename TransformType>
+void check_stdvector_transform(const TransformType&)
+{
+ typedef typename TransformType::MatrixType MatrixType;
+ TransformType x(MatrixType::Random()), y(MatrixType::Random());
+ std::vector<TransformType> v(10), w(20, y);
+ v[5] = x;
+ w[6] = v[5];
+ VERIFY_IS_APPROX(w[6], v[5]);
+ v = w;
+ for(int i = 0; i < 20; i++)
+ {
+ VERIFY_IS_APPROX(w[i], v[i]);
+ }
+
+ v.resize(21);
+ v[20] = x;
+ VERIFY_IS_APPROX(v[20], x);
+ v.resize(22,y);
+ VERIFY_IS_APPROX(v[21], y);
+ v.push_back(x);
+ VERIFY_IS_APPROX(v[22], x);
+ VERIFY((size_t)&(v[22]) == (size_t)&(v[21]) + sizeof(TransformType));
+
+ // do a lot of push_back such that the vector gets internally resized
+ // (with memory reallocation)
+ TransformType* ref = &w[0];
+ for(int i=0; i<30 || ((ref==&w[0]) && i<300); ++i)
+ v.push_back(w[i%w.size()]);
+ for(unsigned int i=23; i<v.size(); ++i)
+ {
+ VERIFY(v[i].matrix()==w[(i-23)%w.size()].matrix());
+ }
+}
+
+template<typename QuaternionType>
+void check_stdvector_quaternion(const QuaternionType&)
+{
+ typedef typename QuaternionType::Coefficients Coefficients;
+ QuaternionType x(Coefficients::Random()), y(Coefficients::Random());
+ std::vector<QuaternionType> v(10), w(20, y);
+ v[5] = x;
+ w[6] = v[5];
+ VERIFY_IS_APPROX(w[6], v[5]);
+ v = w;
+ for(int i = 0; i < 20; i++)
+ {
+ VERIFY_IS_APPROX(w[i], v[i]);
+ }
+
+ v.resize(21);
+ v[20] = x;
+ VERIFY_IS_APPROX(v[20], x);
+ v.resize(22,y);
+ VERIFY_IS_APPROX(v[21], y);
+ v.push_back(x);
+ VERIFY_IS_APPROX(v[22], x);
+ VERIFY((size_t)&(v[22]) == (size_t)&(v[21]) + sizeof(QuaternionType));
+
+ // do a lot of push_back such that the vector gets internally resized
+ // (with memory reallocation)
+ QuaternionType* ref = &w[0];
+ for(int i=0; i<30 || ((ref==&w[0]) && i<300); ++i)
+ v.push_back(w[i%w.size()]);
+ for(unsigned int i=23; i<v.size(); ++i)
+ {
+ VERIFY(v[i].coeffs()==w[(i-23)%w.size()].coeffs());
+ }
+}
+
void test_stdvector()
{
// some non vectorizable fixed sizes
- CALL_SUBTEST(check_stdvector(Vector2f()));
- CALL_SUBTEST(check_stdvector(Matrix3f()));
- CALL_SUBTEST(check_stdvector(Matrix3d()));
+ CALL_SUBTEST(check_stdvector_matrix(Vector2f()));
+ CALL_SUBTEST(check_stdvector_matrix(Matrix3f()));
+ CALL_SUBTEST(check_stdvector_matrix(Matrix3d()));
// some vectorizable fixed sizes
- CALL_SUBTEST(check_stdvector(Matrix2f()));
- CALL_SUBTEST(check_stdvector(Vector4f()));
- CALL_SUBTEST(check_stdvector(Matrix4f()));
- CALL_SUBTEST(check_stdvector(Matrix4d()));
+ CALL_SUBTEST(check_stdvector_matrix(Matrix2f()));
+ CALL_SUBTEST(check_stdvector_matrix(Vector4f()));
+ CALL_SUBTEST(check_stdvector_matrix(Matrix4f()));
+ CALL_SUBTEST(check_stdvector_matrix(Matrix4d()));
// some dynamic sizes
- CALL_SUBTEST(check_stdvector(MatrixXd(1,1)));
- CALL_SUBTEST(check_stdvector(VectorXd(20)));
- CALL_SUBTEST(check_stdvector(RowVectorXf(20)));
- CALL_SUBTEST(check_stdvector(MatrixXcf(10,10)));
+ CALL_SUBTEST(check_stdvector_matrix(MatrixXd(1,1)));
+ CALL_SUBTEST(check_stdvector_matrix(VectorXd(20)));
+ CALL_SUBTEST(check_stdvector_matrix(RowVectorXf(20)));
+ CALL_SUBTEST(check_stdvector_matrix(MatrixXcf(10,10)));
+
+ // some Transform
+ CALL_SUBTEST(check_stdvector_transform(Transform2f()));
+ CALL_SUBTEST(check_stdvector_transform(Transform3f()));
+ CALL_SUBTEST(check_stdvector_transform(Transform3d()));
+ //CALL_SUBTEST(check_stdvector_transform(Transform4d()));
+
+ // some Quaternion
+ CALL_SUBTEST(check_stdvector_quaternion(Quaternionf()));
+ CALL_SUBTEST(check_stdvector_quaternion(Quaternionf()));
}