ndarray_linalg/
generate.rs

1//! Generator functions for matrices
2
3use ndarray::*;
4use rand::prelude::*;
5
6use super::convert::*;
7use super::error::*;
8use super::qr::*;
9use super::types::*;
10
11/// Hermite conjugate matrix
12pub fn conjugate<A, Si, So>(a: &ArrayBase<Si, Ix2>) -> ArrayBase<So, Ix2>
13where
14    A: Scalar,
15    Si: Data<Elem = A>,
16    So: DataOwned<Elem = A> + DataMut,
17{
18    let mut a: ArrayBase<So, Ix2> = replicate(&a.t());
19    for val in a.iter_mut() {
20        *val = val.conj();
21    }
22    a
23}
24
25/// Generate random array with given shape
26///
27/// - This function uses [rand::thread_rng].
28///   See [random_using] for using another RNG
29pub fn random<A, S, Sh, D>(sh: Sh) -> ArrayBase<S, D>
30where
31    A: Scalar,
32    S: DataOwned<Elem = A>,
33    D: Dimension,
34    Sh: ShapeBuilder<Dim = D>,
35{
36    let mut rng = thread_rng();
37    random_using(sh, &mut rng)
38}
39
40/// Generate random array with given RNG
41///
42/// - See [random] for using default RNG
43pub fn random_using<A, S, Sh, D, R>(sh: Sh, rng: &mut R) -> ArrayBase<S, D>
44where
45    A: Scalar,
46    S: DataOwned<Elem = A>,
47    D: Dimension,
48    Sh: ShapeBuilder<Dim = D>,
49    R: Rng,
50{
51    ArrayBase::from_shape_fn(sh, |_| A::rand(rng))
52}
53
54/// Generate random unitary matrix using QR decomposition
55///
56/// - Be sure that this it **NOT** a uniform distribution.
57///   Use it only for test purpose.
58/// - This function uses [rand::thread_rng].
59///   See [random_unitary_using] for using another RNG.
60pub fn random_unitary<A>(n: usize) -> Array2<A>
61where
62    A: Scalar + Lapack,
63{
64    let mut rng = thread_rng();
65    random_unitary_using(n, &mut rng)
66}
67
68/// Generate random unitary matrix using QR decomposition with given RNG
69///
70/// - Be sure that this it **NOT** a uniform distribution.
71///   Use it only for test purpose.
72/// - See [random_unitary] for using default RNG.
73pub fn random_unitary_using<A, R>(n: usize, rng: &mut R) -> Array2<A>
74where
75    A: Scalar + Lapack,
76    R: Rng,
77{
78    let a: Array2<A> = random_using((n, n), rng);
79    let (q, _r) = a.qr_into().unwrap();
80    q
81}
82
83/// Generate random regular matrix
84///
85/// - Be sure that this it **NOT** a uniform distribution.
86///   Use it only for test purpose.
87/// - This function uses [rand::thread_rng].
88///   See [random_regular_using] for using another RNG.
89pub fn random_regular<A>(n: usize) -> Array2<A>
90where
91    A: Scalar + Lapack,
92{
93    let mut rng = rand::thread_rng();
94    random_regular_using(n, &mut rng)
95}
96
97/// Generate random regular matrix with given RNG
98///
99/// - Be sure that this it **NOT** a uniform distribution.
100///   Use it only for test purpose.
101/// - See [random_regular] for using default RNG.
102pub fn random_regular_using<A, R>(n: usize, rng: &mut R) -> Array2<A>
103where
104    A: Scalar + Lapack,
105    R: Rng,
106{
107    let a: Array2<A> = random_using((n, n), rng);
108    let (q, mut r) = a.qr_into().unwrap();
109    for i in 0..n {
110        r[(i, i)] = A::one() + A::from_real(r[(i, i)].abs());
111    }
112    q.dot(&r)
113}
114
115/// Random Hermite matrix
116///
117/// - This function uses [rand::thread_rng].
118///   See [random_hermite_using] for using another RNG.
119pub fn random_hermite<A, S>(n: usize) -> ArrayBase<S, Ix2>
120where
121    A: Scalar,
122    S: DataOwned<Elem = A> + DataMut,
123{
124    let mut rng = rand::thread_rng();
125    random_hermite_using(n, &mut rng)
126}
127
128/// Random Hermite matrix with given RNG
129///
130/// - See [random_hermite] for using default RNG.
131pub fn random_hermite_using<A, S, R>(n: usize, rng: &mut R) -> ArrayBase<S, Ix2>
132where
133    A: Scalar,
134    S: DataOwned<Elem = A> + DataMut,
135    R: Rng,
136{
137    let mut a: ArrayBase<S, Ix2> = random_using((n, n), rng);
138    for i in 0..n {
139        a[(i, i)] = a[(i, i)] + a[(i, i)].conj();
140        for j in (i + 1)..n {
141            a[(i, j)] = a[(j, i)].conj();
142        }
143    }
144    a
145}
146
147/// Random Hermite Positive-definite matrix
148///
149/// - Eigenvalue of matrix must be larger than 1 (thus non-singular)
150/// - This function uses [rand::thread_rng].
151///   See [random_hpd_using] for using another RNG.
152///
153pub fn random_hpd<A, S>(n: usize) -> ArrayBase<S, Ix2>
154where
155    A: Scalar,
156    S: DataOwned<Elem = A> + DataMut,
157{
158    let mut rng = rand::thread_rng();
159    random_hpd_using(n, &mut rng)
160}
161
162/// Random Hermite Positive-definite matrix with given RNG
163///
164/// - Eigenvalue of matrix must be larger than 1 (thus non-singular)
165/// - See [random_hpd] for using default RNG.
166///
167pub fn random_hpd_using<A, S, R>(n: usize, rng: &mut R) -> ArrayBase<S, Ix2>
168where
169    A: Scalar,
170    S: DataOwned<Elem = A> + DataMut,
171    R: Rng,
172{
173    let a: Array2<A> = random_using((n, n), rng);
174    let ah: Array2<A> = conjugate(&a);
175    ArrayBase::eye(n) + &ah.dot(&a)
176}
177
178/// construct matrix from diag
179pub fn from_diag<A>(d: &[A]) -> Array2<A>
180where
181    A: Scalar,
182{
183    let n = d.len();
184    let mut e = Array::zeros((n, n));
185    for i in 0..n {
186        e[(i, i)] = d[i];
187    }
188    e
189}
190
191/// stack vectors into matrix horizontally
192pub fn hstack<A, S>(xs: &[ArrayBase<S, Ix1>]) -> Result<Array<A, Ix2>>
193where
194    A: Scalar,
195    S: Data<Elem = A>,
196{
197    let views: Vec<_> = xs.iter().map(|x| x.view()).collect();
198    stack(Axis(1), &views).map_err(Into::into)
199}
200
201/// stack vectors into matrix vertically
202pub fn vstack<A, S>(xs: &[ArrayBase<S, Ix1>]) -> Result<Array<A, Ix2>>
203where
204    A: Scalar,
205    S: Data<Elem = A>,
206{
207    let views: Vec<_> = xs.iter().map(|x| x.view()).collect();
208    stack(Axis(0), &views).map_err(Into::into)
209}