1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146

// Matrix properties checks.
use approx::RelativeEq;
use num::{One, Zero};

use simba::scalar::{ClosedAdd, ClosedMul, ComplexField, RealField};

use crate::base::allocator::Allocator;
use crate::base::dimension::{Dim, DimMin};
use crate::base::storage::Storage;
use crate::base::{DefaultAllocator, Matrix, Scalar, SquareMatrix};

impl<N: Scalar, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
    /// The total number of elements of this matrix.
    ///
    /// # Examples:
    ///
    /// ```
    /// # use nalgebra::Matrix3x4;
    /// let mat = Matrix3x4::<f32>::zeros();
    /// assert_eq!(mat.len(), 12);
    /// ```
    #[inline]
    pub fn len(&self) -> usize {
        let (nrows, ncols) = self.shape();
        nrows * ncols
    }

    /// Returns true if the matrix contains no elements.
    ///
    /// # Examples:
    ///
    /// ```
    /// # use nalgebra::Matrix3x4;
    /// let mat = Matrix3x4::<f32>::zeros();
    /// assert!(!mat.is_empty());
    /// ```
    #[inline]
    pub fn is_empty(&self) -> bool {
        self.len() == 0
    }

    /// Indicates if this is a square matrix.
    #[inline]
    pub fn is_square(&self) -> bool {
        let (nrows, ncols) = self.shape();
        nrows == ncols
    }

    // TODO: RelativeEq prevents us from using those methods on integer matrices…
    /// Indicated if this is the identity matrix within a relative error of `eps`.
    ///
    /// If the matrix is diagonal, this checks that diagonal elements (i.e. at coordinates `(i, i)`
    /// for i from `0` to `min(R, C)`) are equal one; and that all other elements are zero.
    #[inline]
    pub fn is_identity(&self, eps: N::Epsilon) -> bool
    where
        N: Zero + One + RelativeEq,
        N::Epsilon: Copy,
    {
        let (nrows, ncols) = self.shape();
        let d;

        if nrows > ncols {
            d = ncols;

            for i in d..nrows {
                for j in 0..ncols {
                    if !relative_eq!(self[(i, j)], N::zero(), epsilon = eps) {
                        return false;
                    }
                }
            }
        } else {
            // nrows <= ncols
            d = nrows;

            for i in 0..nrows {
                for j in d..ncols {
                    if !relative_eq!(self[(i, j)], N::zero(), epsilon = eps) {
                        return false;
                    }
                }
            }
        }

        // Off-diagonal elements of the sub-square matrix.
        for i in 1..d {
            for j in 0..i {
                // TODO: use unsafe indexing.
                if !relative_eq!(self[(i, j)], N::zero(), epsilon = eps)
                    || !relative_eq!(self[(j, i)], N::zero(), epsilon = eps)
                {
                    return false;
                }
            }
        }

        // Diagonal elements of the sub-square matrix.
        for i in 0..d {
            if !relative_eq!(self[(i, i)], N::one(), epsilon = eps) {
                return false;
            }
        }

        true
    }
}

impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
    /// Checks that `Mᵀ × M = Id`.
    ///
    /// In this definition `Id` is approximately equal to the identity matrix with a relative error
    /// equal to `eps`.
    #[inline]
    pub fn is_orthogonal(&self, eps: N::Epsilon) -> bool
    where
        N: Zero + One + ClosedAdd + ClosedMul + RelativeEq,
        S: Storage<N, R, C>,
        N::Epsilon: Copy,
        DefaultAllocator: Allocator<N, R, C> + Allocator<N, C, C>,
    {
        (self.ad_mul(self)).is_identity(eps)
    }
}

impl<N: RealField, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S>
where
    DefaultAllocator: Allocator<N, D, D>,
{
    /// Checks that this matrix is orthogonal and has a determinant equal to 1.
    #[inline]
    pub fn is_special_orthogonal(&self, eps: N) -> bool
    where
        D: DimMin<D, Output = D>,
        DefaultAllocator: Allocator<(usize, usize), D>,
    {
        self.is_square() && self.is_orthogonal(eps) && self.determinant() > N::zero()
    }

    /// Returns `true` if this matrix is invertible.
    #[inline]
    pub fn is_invertible(&self) -> bool {
        // TODO: improve this?
        self.clone_owned().try_inverse().is_some()
    }
}