Initial proptest implementation for nalgebra

This introduces functionality for creating strategies for matrices and
vectors, as well as an implementation of Arbitrary.

Strategies for the geometric types (Point3, Quaternion etc.) are not
currently part of this contribution.

The current strategy implementation for matrices has some limitations
that lead to suboptimal shrinking behavior. This is documented in the
module-level docs, with some additional comments in the code.

.circleci/config.yml

@@ -57,7 +57,7 @@ jobs:
       - checkout
       - run:
           name: test
-          command: cargo test --features arbitrary --features serde-serialize --features abomonation-serialize --features sparse --features debug --features io --features compare --features libm
+          command: cargo test --features arbitrary --features serde-serialize --features abomonation-serialize --features sparse --features debug --features io --features compare --features libm --features proptest
       - run:
           name: test nalgebra-glm
           command: cargo test -p nalgebra-glm  --features arbitrary --features serde-serialize --features abomonation-serialize --features sparse --features debug --features io --features compare --features libm

Cargo.toml

@@ -54,6 +54,7 @@ quickcheck     = { version = "0.9", optional = true }
 pest           = { version = "2", optional = true }
 pest_derive    = { version = "2", optional = true }
 matrixcompare-core = { version = "0.1", optional = true }
+proptest = { version = "0.10", optional = true, default-features = false, features = ["std"] }
 
 [dev-dependencies]
 serde_json = "1.0"
@@ -68,6 +69,11 @@ rand_isaac = "0.2"
 # For matrix comparison macro
 matrixcompare = "0.1.3"
 
+# Make sure that we use a specific version of proptest for tests. The reason is that we use a deterministic
+# RNG for certain tests. However, different versions of proptest may give different sequences of numbers,
+# which may cause more brittle tests (although ideally they should take enough samples for it not to matter).
+proptest = { version = "=0.10.1" }
+
 [workspace]
 members = [ "nalgebra-lapack", "nalgebra-glm" ]
 
@@ -78,3 +84,7 @@ path = "benches/lib.rs"
 
 [profile.bench]
 lto = true
+
+[package.metadata.docs.rs]
+# Enable certain features when building docs for docs.rs
+features = [ "proptest" ]

src/lib.rs

@@ -127,6 +127,8 @@ pub mod geometry;
 #[cfg(feature = "io")]
 pub mod io;
 pub mod linalg;
+#[cfg(feature = "proptest")]
+pub mod proptest;
 #[cfg(feature = "sparse")]
 pub mod sparse;
 

src/proptest/mod.rs

@@ -0,0 +1,469 @@
+//! `proptest`-related features for `nalgebra` data structures.
+//!
+//! **This module is only available when the `proptest` feature is enabled in `nalgebra`**.
+//!
+//! `proptest` is a library for *property-based testing*. While similar to QuickCheck,
+//! which may be more familiar to some users, it has a more sophisticated design that
+//! provides users with automatic invariant-preserving shrinking. This means that when using
+//! `proptest`, you rarely need to write your own shrinkers - which is usually very difficult -
+//! and can instead get this "for free". Moreover, `proptest` does not rely on a canonical
+//! `Arbitrary` trait implementation like QuickCheck, though it does also provide this. For
+//! more information, check out the [proptest docs](https://docs.rs/proptest/0.10.1/proptest/)
+//! and the [proptest book](https://altsysrq.github.io/proptest-book/intro.html).
+//!
+//! This module provides users of `nalgebra` with tools to work with `nalgebra` types in
+//! `proptest` tests. At present, this integration is at an early stage, and only
+//! provides tools for generating matrices and vectors, and not any of the geometry types.
+//! There are essentially two ways of using this functionality:
+//!
+//! - Using the [matrix](fn.matrix.html) function to generate matrices with constraints
+//!   on dimensions and elements.
+//! - Relying on the `Arbitrary` implementation of `MatrixMN`.
+//!
+//! The first variant is almost always preferred in practice. Read on to discover why.
+//!
+//! ### Using free function strategies
+//!
+//! In `proptest`, it is usually preferable to have free functions that generate *strategies*.
+//! Currently, the [matrix](fn.matrix.html) function fills this role. The analogous function for
+//! column vectors is [vector](fn.vector.html). Let's take a quick look at how it may be used:
+//! ```rust
+//! use nalgebra::proptest::matrix;
+//! use proptest::prelude::*;
+//!
+//! proptest! {
+//!     # /*
+//!     #[test]
+//!     # */
+//!     fn my_test(a in matrix(-5 ..= 5, 2 ..= 4, 1..=4)) {
+//!         // Generates matrices with elements in the range -5 ..= 5, rows in 2..=4 and
+//!         // columns in 1..=4.
+//!     }
+//! }
+//!
+//! # fn main() { my_test(); }
+//! ```
+//!
+//! In the above example, we generate matrices with constraints on the elements, as well as the
+//! on the allowed dimensions. When a failing example is found, the resulting shrinking process
+//! will preserve these invariants. We can use this to compose more advanced strategies.
+//! For example, let's consider a toy example where we need to generate pairs of matrices
+//! with exactly 3 rows fixed at compile-time and the same number of columns, but we want the
+//! number of columns to vary. One way to do this is to use `proptest` combinators in combination
+//! with [matrix](fn.matrix.html) as follows:
+//!
+//! ```rust
+//! use nalgebra::{Dynamic, MatrixMN, U3};
+//! use nalgebra::proptest::matrix;
+//! use proptest::prelude::*;
+//!
+//! type MyMatrix = MatrixMN<i32, U3, Dynamic>;
+//!
+//! /// Returns a strategy for pairs of matrices with `U3` rows and the same number of
+//! /// columns.
+//! fn matrix_pairs() -> impl Strategy<Value=(MyMatrix, MyMatrix)> {
+//!     matrix(-5 ..= 5, U3, 0 ..= 10)
+//!         // We first generate the initial matrix `a`, and then depending on the concrete
+//!         // instances of `a`, we pick a second matrix with the same number of columns
+//!         .prop_flat_map(|a| {
+//!             let b = matrix(-5 .. 5, U3, a.ncols());
+//!             // This returns a new tuple strategy where we keep `a` fixed while
+//!             // the second item is a strategy that generates instances with the same
+//!             // dimensions as `a`
+//!             (Just(a), b)
+//!         })
+//! }
+//!
+//! proptest! {
+//!     # /*
+//!     #[test]
+//!     # */
+//!     fn my_test((a, b) in matrix_pairs()) {
+//!         // Let's double-check that the two matrices do indeed have the same number of
+//!         // columns
+//!         prop_assert_eq!(a.ncols(), b.ncols());
+//!     }
+//! }
+//!
+//! # fn main() { my_test(); }
+//! ```
+//!
+//! ### The `Arbitrary` implementation
+//!
+//! If you don't care about the dimensions of matrices, you can write tests like these:
+//!
+//! ```rust
+//! use nalgebra::{DMatrix, DVector, Dynamic, Matrix3, MatrixMN, Vector3, U3};
+//! use proptest::prelude::*;
+//!
+//! proptest! {
+//!     # /*
+//!     #[test]
+//!     # */
+//!     fn test_dynamic(matrix: DMatrix<i32>) {
+//!         // This will generate arbitrary instances of `DMatrix` and also attempt
+//!         // to shrink/simplify them when test failures are encountered.
+//!     }
+//!
+//!     # /*
+//!     #[test]
+//!     # */
+//!     fn test_static_and_mixed(matrix: Matrix3<i32>, matrix2: MatrixMN<i32, U3, Dynamic>) {
+//!         // Test some property involving these matrices
+//!     }
+//!
+//!     # /*
+//!     #[test]
+//!     # */
+//!     fn test_vectors(fixed_size_vector: Vector3<i32>, dyn_vector: DVector<i32>) {
+//!         // Test some property involving these vectors
+//!     }
+//! }
+//!
+//! # fn main() { test_dynamic(); test_static_and_mixed(); test_vectors(); }
+//! ```
+//!
+//! While this may be convenient, the default strategies for built-in types in `proptest` can
+//! generate *any* number, including integers large enough to easily lead to overflow when used in
+//! matrix operations, or even infinity or NaN values for floating-point types. Therefore
+//! `Arbitrary` is rarely the method of choice for writing property-based tests.
+//!
+//! ### Notes on shrinking
+//!
+//! Due to some limitations of the current implementation, shrinking takes place by first
+//! shrinking the matrix elements before trying to shrink the dimensions of the matrix.
+//! This unfortunately often leads to the fact that a large number of shrinking iterations
+//! are necessary to find a (nearly) minimal failing test case. As a workaround for this,
+//! you can increase the maximum number of shrinking iterations when debugging. To do this,
+//! simply set the `PROPTEST_MAX_SHRINK_ITERS` variable to a high number. For example:
+//!
+//! ```text
+//! PROPTEST_MAX_SHRINK_ITERS=100000 cargo test my_failing_test
+//! ```
+use crate::allocator::Allocator;
+use crate::{DefaultAllocator, Dim, DimName, Dynamic, MatrixMN, Scalar, U1};
+use proptest::arbitrary::Arbitrary;
+use proptest::collection::vec;
+use proptest::strategy::{BoxedStrategy, Just, NewTree, Strategy, ValueTree};
+use proptest::test_runner::TestRunner;
+
+use std::ops::RangeInclusive;
+
+/// Parameters for arbitrary matrix generation.
+#[derive(Debug, Clone)]
+#[non_exhaustive]
+pub struct MatrixParameters<NParameters, R, C> {
+    /// The range of rows that may be generated.
+    pub rows: DimRange<R>,
+    /// The range of columns that may be generated.
+    pub cols: DimRange<C>,
+    /// Parameters for the `Arbitrary` implementation of the scalar values.
+    pub value_parameters: NParameters,
+}
+
+/// A range of allowed dimensions for use in generation of matrices.
+///
+/// The `DimRange` type is used to encode the range of dimensions that can be used for generation
+/// of matrices with `proptest`. In most cases, you do not need to concern yourself with
+/// `DimRange` directly, as it supports conversion from other types such as `U3` or inclusive
+/// ranges such as `5 ..= 6`. The latter example corresponds to dimensions from (inclusive)
+/// `Dynamic::new(5)` to `Dynamic::new(6)` (inclusive).
+#[derive(Debug, Clone, PartialEq, Eq)]
+pub struct DimRange<D>(RangeInclusive<D>);
+
+impl<D: Dim> DimRange<D> {
+    /// The lower bound for dimensions generated.
+    pub fn lower_bound(&self) -> D {
+        *self.0.start()
+    }
+
+    /// The upper bound for dimensions generated.
+    pub fn upper_bound(&self) -> D {
+        *self.0.end()
+    }
+}
+
+impl<D: Dim> From<D> for DimRange<D> {
+    fn from(dim: D) -> Self {
+        DimRange(dim..=dim)
+    }
+}
+
+impl<D: Dim> From<RangeInclusive<D>> for DimRange<D> {
+    fn from(range: RangeInclusive<D>) -> Self {
+        DimRange(range)
+    }
+}
+
+impl From<RangeInclusive<usize>> for DimRange<Dynamic> {
+    fn from(range: RangeInclusive<usize>) -> Self {
+        DimRange::from(Dynamic::new(*range.start())..=Dynamic::new(*range.end()))
+    }
+}
+
+impl From<usize> for DimRange<Dynamic> {
+    fn from(dim: usize) -> Self {
+        DimRange::from(Dynamic::new(dim))
+    }
+}
+
+/// The default range used for Dynamic dimensions when generating arbitrary matrices.
+fn dynamic_dim_range() -> DimRange<Dynamic> {
+    DimRange::from(0..=6)
+}
+
+/// Create a strategy to generate matrices containing values drawn from the given strategy,
+/// with rows and columns in the provided ranges.
+///
+/// ## Examples
+/// ```
+/// use nalgebra::proptest::matrix;
+/// use nalgebra::{MatrixMN, U3, Dynamic};
+/// use proptest::prelude::*;
+///
+/// proptest! {
+///     # /*
+///     #[test]
+///     # */
+///     fn my_test(a in matrix(0 .. 5i32, U3, 0 ..= 5)) {
+///         // Let's make sure we've got the correct type first
+///         let a: MatrixMN<_, U3, Dynamic> = a;
+///         prop_assert!(a.nrows() == 3);
+///         prop_assert!(a.ncols() <= 5);
+///         prop_assert!(a.iter().all(|x_ij| *x_ij >= 0 && *x_ij < 5));
+///     }
+/// }
+///
+/// # fn main() { my_test(); }
+/// ```
+///
+/// ## Limitations
+/// The current implementation has some limitations that lead to suboptimal shrinking behavior.
+/// See the [module-level documentation](index.html) for more.
+pub fn matrix<R, C, ScalarStrategy>(
+    value_strategy: ScalarStrategy,
+    rows: impl Into<DimRange<R>>,
+    cols: impl Into<DimRange<C>>,
+) -> MatrixStrategy<ScalarStrategy, R, C>
+where
+    ScalarStrategy: Strategy + Clone + 'static,
+    ScalarStrategy::Value: Scalar,
+    R: Dim,
+    C: Dim,
+    DefaultAllocator: Allocator<ScalarStrategy::Value, R, C>,
+{
+    matrix_(value_strategy, rows.into(), cols.into())
+}
+
+/// Same as `matrix`, but without the additional anonymous generic types
+fn matrix_<R, C, ScalarStrategy>(
+    value_strategy: ScalarStrategy,
+    rows: DimRange<R>,
+    cols: DimRange<C>,
+) -> MatrixStrategy<ScalarStrategy, R, C>
+where
+    ScalarStrategy: Strategy + Clone + 'static,
+    ScalarStrategy::Value: Scalar,
+    R: Dim,
+    C: Dim,
+    DefaultAllocator: Allocator<ScalarStrategy::Value, R, C>,
+{
+    let nrows = rows.lower_bound().value()..=rows.upper_bound().value();
+    let ncols = cols.lower_bound().value()..=cols.upper_bound().value();
+
+    // Even though we can use this function to generate fixed-size matrices,
+    // we currently generate all matrices with heap allocated Vec data.
+    // TODO: Avoid heap allocation for fixed-size matrices.
+    // Doing this *properly* would probably require us to implement a custom
+    // strategy and valuetree with custom shrinking logic, which is not trivial
+
+    // Perhaps more problematic, however, is the poor shrinking behavior the current setup leads to.
+    // Shrinking in proptest basically happens in "reverse" of the combinators, so
+    // by first generating the dimensions and then the elements, we get shrinking that first
+    // tries to completely shrink the individual elements before trying to reduce the dimension.
+    // This is clearly the opposite of what we want. I can't find any good way around this
+    // short of writing our own custom value tree, which we should probably do at some point.
+    // TODO: Custom implementation of value tree for better shrinking behavior.
+
+    let strategy = nrows
+        .prop_flat_map(move |nrows| (Just(nrows), ncols.clone()))
+        .prop_flat_map(move |(nrows, ncols)| {
+            (
+                Just(nrows),
+                Just(ncols),
+                vec(value_strategy.clone(), nrows * ncols),
+            )
+        })
+        .prop_map(|(nrows, ncols, values)| {
+            // Note: R/C::from_usize will panic if nrows/ncols does not fit in the dimension type.
+            // However, this should never fail, because we should only be generating
+            // this stuff in the first place
+            MatrixMN::from_iterator_generic(R::from_usize(nrows), C::from_usize(ncols), values)
+        })
+        .boxed();
+
+    MatrixStrategy { strategy }
+}
+
+/// Create a strategy to generate column vectors containing values drawn from the given strategy,
+/// with length in the provided range.
+///
+/// This is a convenience function for calling
+/// [matrix(value_strategy, length, U1)](fn.matrix.html) and should
+/// be used when you only want to generate column vectors, as it's simpler and makes the intent
+/// clear.
+pub fn vector<D, ScalarStrategy>(
+    value_strategy: ScalarStrategy,
+    length: impl Into<DimRange<D>>,
+) -> MatrixStrategy<ScalarStrategy, D, U1>
+where
+    ScalarStrategy: Strategy + Clone + 'static,
+    ScalarStrategy::Value: Scalar,
+    D: Dim,
+    DefaultAllocator: Allocator<ScalarStrategy::Value, D>,
+{
+    matrix_(value_strategy, length.into(), U1.into())
+}
+
+impl<NParameters, R, C> Default for MatrixParameters<NParameters, R, C>
+where
+    NParameters: Default,
+    R: DimName,
+    C: DimName,
+{
+    fn default() -> Self {
+        Self {
+            rows: DimRange::from(R::name()),
+            cols: DimRange::from(C::name()),
+            value_parameters: NParameters::default(),
+        }
+    }
+}
+
+impl<NParameters, R> Default for MatrixParameters<NParameters, R, Dynamic>
+where
+    NParameters: Default,
+    R: DimName,
+{
+    fn default() -> Self {
+        Self {
+            rows: DimRange::from(R::name()),
+            cols: dynamic_dim_range(),
+            value_parameters: NParameters::default(),
+        }
+    }
+}
+
+impl<NParameters, C> Default for MatrixParameters<NParameters, Dynamic, C>
+where
+    NParameters: Default,
+    C: DimName,
+{
+    fn default() -> Self {
+        Self {
+            rows: dynamic_dim_range(),
+            cols: DimRange::from(C::name()),
+            value_parameters: NParameters::default(),
+        }
+    }
+}
+
+impl<NParameters> Default for MatrixParameters<NParameters, Dynamic, Dynamic>
+where
+    NParameters: Default,
+{
+    fn default() -> Self {
+        Self {
+            rows: dynamic_dim_range(),
+            cols: dynamic_dim_range(),
+            value_parameters: NParameters::default(),
+        }
+    }
+}
+
+impl<N, R, C> Arbitrary for MatrixMN<N, R, C>
+where
+    N: Scalar + Arbitrary,
+    <N as Arbitrary>::Strategy: Clone,
+    R: Dim,
+    C: Dim,
+    MatrixParameters<N::Parameters, R, C>: Default,
+    DefaultAllocator: Allocator<N, R, C>,
+{
+    type Parameters = MatrixParameters<N::Parameters, R, C>;
+
+    fn arbitrary_with(args: Self::Parameters) -> Self::Strategy {
+        let value_strategy = N::arbitrary_with(args.value_parameters);
+        matrix(value_strategy, args.rows, args.cols)
+    }
+
+    type Strategy = MatrixStrategy<N::Strategy, R, C>;
+}
+
+/// A strategy for generating matrices.
+#[derive(Debug)]
+pub struct MatrixStrategy<NStrategy, R: Dim, C: Dim>
+where
+    NStrategy: Strategy,
+    NStrategy::Value: Scalar,
+    DefaultAllocator: Allocator<NStrategy::Value, R, C>,
+{
+    // For now we only internally hold a boxed strategy. The reason for introducing this
+    // separate wrapper struct is so that we can replace the strategy logic with custom logic
+    // later down the road without introducing significant breaking changes
+    strategy: BoxedStrategy<MatrixMN<NStrategy::Value, R, C>>,
+}
+
+impl<NStrategy, R, C> Strategy for MatrixStrategy<NStrategy, R, C>
+where
+    NStrategy: Strategy,
+    NStrategy::Value: Scalar,
+    R: Dim,
+    C: Dim,
+    DefaultAllocator: Allocator<NStrategy::Value, R, C>,
+{
+    type Tree = MatrixValueTree<NStrategy::Value, R, C>;
+    type Value = MatrixMN<NStrategy::Value, R, C>;
+
+    fn new_tree(&self, runner: &mut TestRunner) -> NewTree<Self> {
+        let underlying_tree = self.strategy.new_tree(runner)?;
+        Ok(MatrixValueTree {
+            value_tree: underlying_tree,
+        })
+    }
+}
+
+/// A value tree for matrices.
+pub struct MatrixValueTree<N, R, C>
+where
+    N: Scalar,
+    R: Dim,
+    C: Dim,
+    DefaultAllocator: Allocator<N, R, C>,
+{
+    // For now we only wrap a boxed value tree. The reason for wrapping is that this allows us
+    // to swap out the value tree logic down the road without significant breaking changes.
+    value_tree: Box<dyn ValueTree<Value = MatrixMN<N, R, C>>>,
+}
+
+impl<N, R, C> ValueTree for MatrixValueTree<N, R, C>
+where
+    N: Scalar,
+    R: Dim,
+    C: Dim,
+    DefaultAllocator: Allocator<N, R, C>,
+{
+    type Value = MatrixMN<N, R, C>;
+
+    fn current(&self) -> Self::Value {
+        self.value_tree.current()
+    }
+
+    fn simplify(&mut self) -> bool {
+        self.value_tree.simplify()
+    }
+
+    fn complicate(&mut self) -> bool {
+        self.value_tree.complicate()
+    }
+}

tests/lib.rs

@@ -19,5 +19,9 @@ extern crate quickcheck;
 mod core;
 mod geometry;
 mod linalg;
+
+#[cfg(feature = "proptest")]
+mod proptest;
+
 //#[cfg(feature = "sparse")]
 //mod sparse;

tests/proptest/mod.rs

@@ -0,0 +1,184 @@
+//! Tests for proptest-related functionality.
+use nalgebra::base::dimension::*;
+use nalgebra::proptest::{matrix, DimRange, MatrixStrategy};
+use nalgebra::{DMatrix, DVector, Dim, Matrix3, Matrix4, MatrixMN, Vector3};
+use proptest::prelude::*;
+use proptest::strategy::ValueTree;
+use proptest::test_runner::TestRunner;
+
+/// Generate a proptest that tests that all matrices generated with the
+/// provided rows and columns conform to the constraints defined by the
+/// input.
+macro_rules! generate_matrix_sanity_test {
+    ($test_name:ident, $rows:expr, $cols:expr) => {
+        proptest! {
+            #[test]
+            fn $test_name(a in matrix(-5 ..= 5i32, $rows, $cols)) {
+                // let a: MatrixMN<_, $rows, $cols> = a;
+                let rows_range = DimRange::from($rows);
+                let cols_range = DimRange::from($cols);
+                prop_assert!(a.nrows() >= rows_range.lower_bound().value()
+                          && a.nrows() <= rows_range.upper_bound().value());
+                prop_assert!(a.ncols() >= cols_range.lower_bound().value()
+                          && a.ncols() <= cols_range.upper_bound().value());
+                prop_assert!(a.iter().all(|x_ij| *x_ij >= -5 && *x_ij <= 5));
+            }
+        }
+    };
+}
+
+// Test all fixed-size matrices with row/col dimensions up to 3
+generate_matrix_sanity_test!(test_matrix_u0_u0, U0, U0);
+generate_matrix_sanity_test!(test_matrix_u1_u0, U1, U0);
+generate_matrix_sanity_test!(test_matrix_u0_u1, U0, U1);
+generate_matrix_sanity_test!(test_matrix_u1_u1, U1, U1);
+generate_matrix_sanity_test!(test_matrix_u2_u1, U2, U1);
+generate_matrix_sanity_test!(test_matrix_u1_u2, U1, U2);
+generate_matrix_sanity_test!(test_matrix_u2_u2, U2, U2);
+generate_matrix_sanity_test!(test_matrix_u3_u2, U3, U2);
+generate_matrix_sanity_test!(test_matrix_u2_u3, U2, U3);
+generate_matrix_sanity_test!(test_matrix_u3_u3, U3, U3);
+
+// Similarly test all heap-allocated but fixed dim ranges
+generate_matrix_sanity_test!(test_matrix_0_0, 0, 0);
+generate_matrix_sanity_test!(test_matrix_0_1, 0, 1);
+generate_matrix_sanity_test!(test_matrix_1_0, 1, 0);
+generate_matrix_sanity_test!(test_matrix_1_1, 1, 1);
+generate_matrix_sanity_test!(test_matrix_2_1, 2, 1);
+generate_matrix_sanity_test!(test_matrix_1_2, 1, 2);
+generate_matrix_sanity_test!(test_matrix_2_2, 2, 2);
+generate_matrix_sanity_test!(test_matrix_3_2, 3, 2);
+generate_matrix_sanity_test!(test_matrix_2_3, 2, 3);
+generate_matrix_sanity_test!(test_matrix_3_3, 3, 3);
+
+// Test arbitrary inputs
+generate_matrix_sanity_test!(test_matrix_input_1, U5, 1..=5);
+generate_matrix_sanity_test!(test_matrix_input_2, 3..=4, 1..=5);
+generate_matrix_sanity_test!(test_matrix_input_3, 1..=2, U3);
+generate_matrix_sanity_test!(test_matrix_input_4, 3, U4);
+
+#[test]
+fn test_matrix_output_types() {
+    // Test that the dimension types are correct for the given inputs
+    let _: MatrixStrategy<_, U3, U4> = matrix(-5..5, U3, U4);
+    let _: MatrixStrategy<_, U3, U3> = matrix(-5..5, U3, U3);
+    let _: MatrixStrategy<_, U3, Dynamic> = matrix(-5..5, U3, 1..=5);
+    let _: MatrixStrategy<_, Dynamic, U3> = matrix(-5..5, 1..=5, U3);
+    let _: MatrixStrategy<_, Dynamic, Dynamic> = matrix(-5..5, 1..=5, 1..=5);
+}
+
+// Below we have some tests to ensure that specific instances of MatrixMN are usable
+// in a typical proptest scenario where we (implicitly) use the `Arbitrary` trait
+proptest! {
+    #[test]
+    fn ensure_arbitrary_test_compiles_matrix3(_: Matrix3<i32>) {}
+
+    #[test]
+    fn ensure_arbitrary_test_compiles_matrixmn_u3_dynamic(_: MatrixMN<i32, U3, Dynamic>) {}
+
+    #[test]
+    fn ensure_arbitrary_test_compiles_matrixmn_dynamic_u3(_: MatrixMN<i32, Dynamic, U3>) {}
+
+    #[test]
+    fn ensure_arbitrary_test_compiles_dmatrix(_: DMatrix<i32>) {}
+
+    #[test]
+    fn ensure_arbitrary_test_compiles_vector3(_: Vector3<i32>) {}
+
+    #[test]
+    fn ensure_arbitrary_test_compiles_dvector(_: DVector<i32>) {}
+}
+
+#[test]
+fn matrix_samples_all_possible_outputs() {
+    // Test that the proptest generation covers all possible outputs for a small space of inputs
+    // given enough samples.
+
+    // We use a deterministic test runner to make the test "stable".
+    let mut runner = TestRunner::deterministic();
+
+    let strategy = matrix(0..=2usize, 0..=3, 0..=3);
+
+    // We use flags to record whether values and combinations of dimensions were encountered.
+    // For example, if we encounter value 1, we set the value flag of 1 to true,
+    // and if we encounted matrix dimensions 4x3, we set the flag of [4, 3] to true.
+    let mut value_encountered = Vector3::new(false, false, false);
+    let mut dimensions_encountered = Matrix4::repeat(false);
+
+    for _ in 0..1000 {
+        let tree = strategy
+            .new_tree(&mut runner)
+            .expect("Tree generation should not fail");
+        let matrix = tree.current();
+
+        dimensions_encountered[(matrix.nrows(), matrix.ncols())] = true;
+
+        for &value in matrix.iter() {
+            value_encountered[value] = true;
+        }
+    }
+
+    assert!(
+        value_encountered.iter().all(|v| *v),
+        "Did not sample all possible values."
+    );
+    assert!(
+        dimensions_encountered.iter().all(|v| *v),
+        "Did not sample all possible matrix dimensions."
+    );
+}
+
+#[test]
+fn matrix_shrinking_satisfies_constraints() {
+    // We use a deterministic test runner to make the test "stable".
+    let mut runner = TestRunner::deterministic();
+
+    let strategy = matrix(-1..=2, 1..=3, 2..=4);
+
+    let num_matrices = 25;
+
+    macro_rules! maybeprintln {
+        ($($arg:tt)*) => {
+            // Uncomment the below line to enable printing of matrix sequences. This is handy
+            // for manually inspecting the sequences of simplified matrices.
+            // println!($($arg)*)
+        };
+    }
+
+    maybeprintln!("========================== (begin generation process)");
+
+    for _ in 0..num_matrices {
+        let mut tree = strategy
+            .new_tree(&mut runner)
+            .expect("Tree generation should not fail.");
+
+        let mut current = Some(tree.current());
+
+        maybeprintln!("------------------");
+
+        while let Some(matrix) = current {
+            maybeprintln!("{}", matrix);
+
+            assert!(
+                matrix.iter().all(|&v| v >= -1 && v <= 2),
+                "All matrix elements must satisfy constraints"
+            );
+            assert!(
+                matrix.nrows() >= 1 && matrix.nrows() <= 3,
+                "Number of rows in matrix must satisfy constraints."
+            );
+            assert!(
+                matrix.ncols() >= 2 && matrix.ncols() <= 4,
+                "Number of columns in matrix must satisfy constraints."
+            );
+
+            current = if tree.simplify() {
+                Some(tree.current())
+            } else {
+                None
+            }
+        }
+    }
+
+    maybeprintln!("========================== (end of generation process)");
+}