commit
012923d275
|
@ -48,9 +48,7 @@ pub fn qr<N, V, M>(m: &M) -> (M, M)
|
||||||
let mut q : M = Eye::new_identity(rows);
|
let mut q : M = Eye::new_identity(rows);
|
||||||
let mut r = *m;
|
let mut r = *m;
|
||||||
|
|
||||||
let iterations = min(rows - 1, cols);
|
for ite in 0..min(rows - 1, cols) {
|
||||||
|
|
||||||
for ite in (0us .. iterations) {
|
|
||||||
let mut v = r.col_slice(ite, ite, rows);
|
let mut v = r.col_slice(ite, ite, rows);
|
||||||
let alpha =
|
let alpha =
|
||||||
if unsafe { v.unsafe_at(ite) } >= ::zero() {
|
if unsafe { v.unsafe_at(ite) } >= ::zero() {
|
||||||
|
@ -84,19 +82,19 @@ pub fn eigen_qr<N, V, VS, M>(m: &M, eps: &N, niter: usize) -> (M, V)
|
||||||
let mut eigenvalues = *m;
|
let mut eigenvalues = *m;
|
||||||
// let mut shifter: M = Eye::new_identity(rows);
|
// let mut shifter: M = Eye::new_identity(rows);
|
||||||
|
|
||||||
let mut iter = 0us;
|
let mut iter = 0;
|
||||||
for _ in (0us .. niter) {
|
for _ in 0..niter {
|
||||||
let mut stop = true;
|
let mut stop = true;
|
||||||
|
|
||||||
for j in (0us .. ::dim::<M>()) {
|
for j in 0..::dim::<M>() {
|
||||||
for i in (0us .. j) {
|
for i in 0..j {
|
||||||
if unsafe { eigenvalues.unsafe_at((i, j)) }.abs() >= *eps {
|
if unsafe { eigenvalues.unsafe_at((i, j)) }.abs() >= *eps {
|
||||||
stop = false;
|
stop = false;
|
||||||
break;
|
break;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
for i in (j + 1 .. ::dim::<M>()) {
|
for i in j + 1..::dim::<M>() {
|
||||||
if unsafe { eigenvalues.unsafe_at((i, j)) }.abs() >= *eps {
|
if unsafe { eigenvalues.unsafe_at((i, j)) }.abs() >= *eps {
|
||||||
stop = false;
|
stop = false;
|
||||||
break;
|
break;
|
||||||
|
|
|
@ -130,7 +130,7 @@ impl<N> DMat<N> {
|
||||||
DMat {
|
DMat {
|
||||||
nrows: nrows,
|
nrows: nrows,
|
||||||
ncols: ncols,
|
ncols: ncols,
|
||||||
mij: (0us .. nrows * ncols).map(|i| { let m = i / nrows; f(i - m * nrows, m) }).collect()
|
mij: (0..nrows * ncols).map(|i| { let m = i / nrows; f(i - m * nrows, m) }).collect()
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
@ -157,7 +157,7 @@ impl<N> DMat<N> {
|
||||||
/// The returned vector contains the matrix data in column-major order.
|
/// The returned vector contains the matrix data in column-major order.
|
||||||
#[inline]
|
#[inline]
|
||||||
pub fn as_vec(&self) -> &[N] {
|
pub fn as_vec(&self) -> &[N] {
|
||||||
&self.mij[]
|
&self.mij
|
||||||
}
|
}
|
||||||
|
|
||||||
/// Gets a mutable reference to this matrix data.
|
/// Gets a mutable reference to this matrix data.
|
||||||
|
@ -180,7 +180,7 @@ impl<N: One + Zero + Clone + Copy> Eye for DMat<N> {
|
||||||
fn new_identity(dim: usize) -> DMat<N> {
|
fn new_identity(dim: usize) -> DMat<N> {
|
||||||
let mut res = DMat::new_zeros(dim, dim);
|
let mut res = DMat::new_zeros(dim, dim);
|
||||||
|
|
||||||
for i in (0us .. dim) {
|
for i in 0..dim {
|
||||||
let _1: N = ::one();
|
let _1: N = ::one();
|
||||||
res[(i, i)] = _1;
|
res[(i, i)] = _1;
|
||||||
}
|
}
|
||||||
|
@ -237,7 +237,7 @@ impl<N: Copy> Indexable<(usize, usize), N> for DMat<N> {
|
||||||
unsafe fn unsafe_at(&self, rowcol: (usize, usize)) -> N {
|
unsafe fn unsafe_at(&self, rowcol: (usize, usize)) -> N {
|
||||||
let (row, col) = rowcol;
|
let (row, col) = rowcol;
|
||||||
|
|
||||||
*self.mij[].get_unchecked(self.offset(row, col))
|
*self.mij.get_unchecked(self.offset(row, col))
|
||||||
}
|
}
|
||||||
|
|
||||||
#[inline]
|
#[inline]
|
||||||
|
@ -269,7 +269,7 @@ impl<N> Index<(usize, usize)> for DMat<N> {
|
||||||
assert!(j < self.ncols);
|
assert!(j < self.ncols);
|
||||||
|
|
||||||
unsafe {
|
unsafe {
|
||||||
self.mij[].get_unchecked(self.offset(i, j))
|
self.mij.get_unchecked(self.offset(i, j))
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -295,12 +295,12 @@ impl<N: Copy + Mul<N, Output = N> + Add<N, Output = N> + Zero> Mul<DMat<N>> for
|
||||||
|
|
||||||
let mut res = unsafe { DMat::new_uninitialized(self.nrows, right.ncols) };
|
let mut res = unsafe { DMat::new_uninitialized(self.nrows, right.ncols) };
|
||||||
|
|
||||||
for i in (0us .. self.nrows) {
|
for i in 0..self.nrows {
|
||||||
for j in (0us .. right.ncols) {
|
for j in 0..right.ncols {
|
||||||
let mut acc: N = ::zero();
|
let mut acc: N = ::zero();
|
||||||
|
|
||||||
unsafe {
|
unsafe {
|
||||||
for k in (0us .. self.ncols) {
|
for k in 0..self.ncols {
|
||||||
acc = acc
|
acc = acc
|
||||||
+ self.unsafe_at((i, k)) * right.unsafe_at((k, j));
|
+ self.unsafe_at((i, k)) * right.unsafe_at((k, j));
|
||||||
}
|
}
|
||||||
|
@ -322,10 +322,10 @@ impl<N: Copy + Add<N, Output = N> + Mul<N, Output = N> + Zero> Mul<DVec<N>> for
|
||||||
|
|
||||||
let mut res : DVec<N> = unsafe { DVec::new_uninitialized(self.nrows) };
|
let mut res : DVec<N> = unsafe { DVec::new_uninitialized(self.nrows) };
|
||||||
|
|
||||||
for i in (0us .. self.nrows) {
|
for i in 0..self.nrows {
|
||||||
let mut acc: N = ::zero();
|
let mut acc: N = ::zero();
|
||||||
|
|
||||||
for j in (0us .. self.ncols) {
|
for j in 0..self.ncols {
|
||||||
unsafe {
|
unsafe {
|
||||||
acc = acc + self.unsafe_at((i, j)) * right.unsafe_at(j);
|
acc = acc + self.unsafe_at((i, j)) * right.unsafe_at(j);
|
||||||
}
|
}
|
||||||
|
@ -347,10 +347,10 @@ impl<N: Copy + Add<N, Output = N> + Mul<N, Output = N> + Zero> Mul<DMat<N>> for
|
||||||
|
|
||||||
let mut res : DVec<N> = unsafe { DVec::new_uninitialized(right.ncols) };
|
let mut res : DVec<N> = unsafe { DVec::new_uninitialized(right.ncols) };
|
||||||
|
|
||||||
for i in (0us .. right.ncols) {
|
for i in 0..right.ncols {
|
||||||
let mut acc: N = ::zero();
|
let mut acc: N = ::zero();
|
||||||
|
|
||||||
for j in (0us .. right.nrows) {
|
for j in 0..right.nrows {
|
||||||
unsafe {
|
unsafe {
|
||||||
acc = acc + self.unsafe_at(j) * right.unsafe_at((j, i));
|
acc = acc + self.unsafe_at(j) * right.unsafe_at((j, i));
|
||||||
}
|
}
|
||||||
|
@ -382,7 +382,7 @@ impl<N: BaseNum + Clone> Inv for DMat<N> {
|
||||||
let mut res: DMat<N> = Eye::new_identity(dim);
|
let mut res: DMat<N> = Eye::new_identity(dim);
|
||||||
|
|
||||||
// inversion using Gauss-Jordan elimination
|
// inversion using Gauss-Jordan elimination
|
||||||
for k in (0us .. dim) {
|
for k in 0..dim {
|
||||||
// search a non-zero value on the k-th column
|
// search a non-zero value on the k-th column
|
||||||
// FIXME: would it be worth it to spend some more time searching for the
|
// FIXME: would it be worth it to spend some more time searching for the
|
||||||
// max instead?
|
// max instead?
|
||||||
|
@ -403,7 +403,7 @@ impl<N: BaseNum + Clone> Inv for DMat<N> {
|
||||||
|
|
||||||
// swap pivot line
|
// swap pivot line
|
||||||
if n0 != k {
|
if n0 != k {
|
||||||
for j in (0us .. dim) {
|
for j in 0..dim {
|
||||||
let off_n0_j = self.offset(n0, j);
|
let off_n0_j = self.offset(n0, j);
|
||||||
let off_k_j = self.offset(k, j);
|
let off_k_j = self.offset(k, j);
|
||||||
|
|
||||||
|
@ -415,26 +415,26 @@ impl<N: BaseNum + Clone> Inv for DMat<N> {
|
||||||
unsafe {
|
unsafe {
|
||||||
let pivot = self.unsafe_at((k, k));
|
let pivot = self.unsafe_at((k, k));
|
||||||
|
|
||||||
for j in (k .. dim) {
|
for j in k..dim {
|
||||||
let selfval = self.unsafe_at((k, j)) / pivot;
|
let selfval = self.unsafe_at((k, j)) / pivot;
|
||||||
self.unsafe_set((k, j), selfval);
|
self.unsafe_set((k, j), selfval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for j in (0us .. dim) {
|
for j in 0..dim {
|
||||||
let resval = res.unsafe_at((k, j)) / pivot;
|
let resval = res.unsafe_at((k, j)) / pivot;
|
||||||
res.unsafe_set((k, j), resval);
|
res.unsafe_set((k, j), resval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for l in (0us .. dim) {
|
for l in 0..dim {
|
||||||
if l != k {
|
if l != k {
|
||||||
let normalizer = self.unsafe_at((l, k));
|
let normalizer = self.unsafe_at((l, k));
|
||||||
|
|
||||||
for j in (k .. dim) {
|
for j in k..dim {
|
||||||
let selfval = self.unsafe_at((l, j)) - self.unsafe_at((k, j)) * normalizer;
|
let selfval = self.unsafe_at((l, j)) - self.unsafe_at((k, j)) * normalizer;
|
||||||
self.unsafe_set((l, j), selfval);
|
self.unsafe_set((l, j), selfval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for j in (0us .. dim) {
|
for j in 0..dim {
|
||||||
let resval = res.unsafe_at((l, j)) - res.unsafe_at((k, j)) * normalizer;
|
let resval = res.unsafe_at((l, j)) - res.unsafe_at((k, j)) * normalizer;
|
||||||
res.unsafe_set((l, j), resval);
|
res.unsafe_set((l, j), resval);
|
||||||
}
|
}
|
||||||
|
@ -462,8 +462,8 @@ impl<N: Clone + Copy> Transpose for DMat<N> {
|
||||||
else {
|
else {
|
||||||
let mut res = unsafe { DMat::new_uninitialized(self.ncols, self.nrows) };
|
let mut res = unsafe { DMat::new_uninitialized(self.ncols, self.nrows) };
|
||||||
|
|
||||||
for i in (0us .. self.nrows) {
|
for i in 0..self.nrows {
|
||||||
for j in (0us .. self.ncols) {
|
for j in 0..self.ncols {
|
||||||
unsafe {
|
unsafe {
|
||||||
res.unsafe_set((j, i), self.unsafe_at((i, j)))
|
res.unsafe_set((j, i), self.unsafe_at((i, j)))
|
||||||
}
|
}
|
||||||
|
@ -477,8 +477,8 @@ impl<N: Clone + Copy> Transpose for DMat<N> {
|
||||||
#[inline]
|
#[inline]
|
||||||
fn transpose_mut(&mut self) {
|
fn transpose_mut(&mut self) {
|
||||||
if self.nrows == self.ncols {
|
if self.nrows == self.ncols {
|
||||||
for i in (1us .. self.nrows) {
|
for i in 1..self.nrows {
|
||||||
for j in (0us .. self.ncols - 1) {
|
for j in 0..self.ncols - 1 {
|
||||||
let off_i_j = self.offset(i, j);
|
let off_i_j = self.offset(i, j);
|
||||||
let off_j_i = self.offset(j, i);
|
let off_j_i = self.offset(j, i);
|
||||||
|
|
||||||
|
@ -500,8 +500,8 @@ impl<N: BaseNum + Cast<f64> + Clone> Mean<DVec<N>> for DMat<N> {
|
||||||
let mut res: DVec<N> = DVec::new_zeros(self.ncols);
|
let mut res: DVec<N> = DVec::new_zeros(self.ncols);
|
||||||
let normalizer: N = Cast::from(1.0f64 / Cast::from(self.nrows));
|
let normalizer: N = Cast::from(1.0f64 / Cast::from(self.nrows));
|
||||||
|
|
||||||
for i in (0us .. self.nrows) {
|
for i in 0..self.nrows {
|
||||||
for j in (0us .. self.ncols) {
|
for j in 0..self.ncols {
|
||||||
unsafe {
|
unsafe {
|
||||||
let acc = res.unsafe_at(j) + self.unsafe_at((i, j)) * normalizer;
|
let acc = res.unsafe_at(j) + self.unsafe_at((i, j)) * normalizer;
|
||||||
res.unsafe_set(j, acc);
|
res.unsafe_set(j, acc);
|
||||||
|
@ -522,8 +522,8 @@ impl<N: BaseNum + Cast<f64> + Clone> Cov<DMat<N>> for DMat<N> {
|
||||||
let mean = self.mean();
|
let mean = self.mean();
|
||||||
|
|
||||||
// FIXME: use the rows iterator when available
|
// FIXME: use the rows iterator when available
|
||||||
for i in (0us .. self.nrows) {
|
for i in 0..self.nrows {
|
||||||
for j in (0us .. self.ncols) {
|
for j in 0..self.ncols {
|
||||||
unsafe {
|
unsafe {
|
||||||
centered.unsafe_set((i, j), self.unsafe_at((i, j)) - mean.unsafe_at(j));
|
centered.unsafe_set((i, j), self.unsafe_at((i, j)) - mean.unsafe_at(j));
|
||||||
}
|
}
|
||||||
|
@ -559,8 +559,8 @@ impl<N: Copy> RowSlice<DVec<N>> for DMat<N> {
|
||||||
let mut slice : DVec<N> = unsafe {
|
let mut slice : DVec<N> = unsafe {
|
||||||
DVec::new_uninitialized(self.nrows)
|
DVec::new_uninitialized(self.nrows)
|
||||||
};
|
};
|
||||||
let mut slice_idx = 0us;
|
let mut slice_idx = 0;
|
||||||
for col_id in (col_start .. col_end) {
|
for col_id in col_start..col_end {
|
||||||
unsafe {
|
unsafe {
|
||||||
slice.unsafe_set(slice_idx, self.unsafe_at((row_id, col_id)));
|
slice.unsafe_set(slice_idx, self.unsafe_at((row_id, col_id)));
|
||||||
}
|
}
|
||||||
|
@ -586,7 +586,7 @@ impl<N: Copy + Clone + Zero> Diag<DVec<N>> for DMat<N> {
|
||||||
|
|
||||||
assert!(diag.len() == smallest_dim);
|
assert!(diag.len() == smallest_dim);
|
||||||
|
|
||||||
for i in (0us .. smallest_dim) {
|
for i in 0..smallest_dim {
|
||||||
unsafe { self.unsafe_set((i, i), diag.unsafe_at(i)) }
|
unsafe { self.unsafe_set((i, i), diag.unsafe_at(i)) }
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -597,7 +597,7 @@ impl<N: Copy + Clone + Zero> Diag<DVec<N>> for DMat<N> {
|
||||||
|
|
||||||
let mut diag: DVec<N> = DVec::new_zeros(smallest_dim);
|
let mut diag: DVec<N> = DVec::new_zeros(smallest_dim);
|
||||||
|
|
||||||
for i in (0us .. smallest_dim) {
|
for i in 0..smallest_dim {
|
||||||
unsafe { diag.unsafe_set(i, self.unsafe_at((i, i))) }
|
unsafe { diag.unsafe_set(i, self.unsafe_at((i, i))) }
|
||||||
}
|
}
|
||||||
|
|
||||||
|
@ -631,8 +631,8 @@ impl<N: ApproxEq<N>> ApproxEq<N> for DMat<N> {
|
||||||
|
|
||||||
impl<N: Debug + Copy + Display> Debug for DMat<N> {
|
impl<N: Debug + Copy + Display> Debug for DMat<N> {
|
||||||
fn fmt(&self, form:&mut Formatter) -> Result {
|
fn fmt(&self, form:&mut Formatter) -> Result {
|
||||||
for i in (0us .. self.nrows()) {
|
for i in 0..self.nrows() {
|
||||||
for j in (0us .. self.ncols()) {
|
for j in 0..self.ncols() {
|
||||||
let _ = write!(form, "{} ", self[(i, j)]);
|
let _ = write!(form, "{} ", self[(i, j)]);
|
||||||
}
|
}
|
||||||
let _ = write!(form, "\n");
|
let _ = write!(form, "\n");
|
||||||
|
|
|
@ -57,7 +57,7 @@ impl<N> DVec<N> {
|
||||||
/// Builds a vector filled with the result of a function.
|
/// Builds a vector filled with the result of a function.
|
||||||
#[inline(always)]
|
#[inline(always)]
|
||||||
pub fn from_fn<F: FnMut(usize) -> N>(dim: usize, mut f: F) -> DVec<N> {
|
pub fn from_fn<F: FnMut(usize) -> N>(dim: usize, mut f: F) -> DVec<N> {
|
||||||
DVec { at: (0us .. dim).map(|i| f(i)).collect() }
|
DVec { at: (0..dim).map(|i| f(i)).collect() }
|
||||||
}
|
}
|
||||||
|
|
||||||
#[inline]
|
#[inline]
|
||||||
|
|
|
@ -128,7 +128,7 @@ macro_rules! dvec_impl(
|
||||||
fn axpy(&mut self, a: &N, x: &$dvec<N>) {
|
fn axpy(&mut self, a: &N, x: &$dvec<N>) {
|
||||||
assert!(self.len() == x.len());
|
assert!(self.len() == x.len());
|
||||||
|
|
||||||
for i in (0us .. x.len()) {
|
for i in 0..x.len() {
|
||||||
unsafe {
|
unsafe {
|
||||||
let self_i = self.unsafe_at(i);
|
let self_i = self.unsafe_at(i);
|
||||||
self.unsafe_set(i, self_i + *a * x.unsafe_at(i))
|
self.unsafe_set(i, self_i + *a * x.unsafe_at(i))
|
||||||
|
@ -144,7 +144,7 @@ macro_rules! dvec_impl(
|
||||||
pub fn canonical_basis_with_dim(dim: usize) -> Vec<$dvec<N>> {
|
pub fn canonical_basis_with_dim(dim: usize) -> Vec<$dvec<N>> {
|
||||||
let mut res : Vec<$dvec<N>> = Vec::new();
|
let mut res : Vec<$dvec<N>> = Vec::new();
|
||||||
|
|
||||||
for i in (0us .. dim) {
|
for i in 0..dim {
|
||||||
let mut basis_element : $dvec<N> = $dvec::new_zeros(dim);
|
let mut basis_element : $dvec<N> = $dvec::new_zeros(dim);
|
||||||
|
|
||||||
basis_element.set(i, ::one());
|
basis_element.set(i, ::one());
|
||||||
|
@ -163,7 +163,7 @@ macro_rules! dvec_impl(
|
||||||
let dim = self.len();
|
let dim = self.len();
|
||||||
let mut res : Vec<$dvec<N>> = Vec::new();
|
let mut res : Vec<$dvec<N>> = Vec::new();
|
||||||
|
|
||||||
for i in (0us .. dim) {
|
for i in 0..dim {
|
||||||
let mut basis_element : $dvec<N> = $dvec::new_zeros(self.len());
|
let mut basis_element : $dvec<N> = $dvec::new_zeros(self.len());
|
||||||
|
|
||||||
basis_element.set(i, ::one());
|
basis_element.set(i, ::one());
|
||||||
|
@ -274,7 +274,7 @@ macro_rules! dvec_impl(
|
||||||
fn dot(&self, other: &$dvec<N>) -> N {
|
fn dot(&self, other: &$dvec<N>) -> N {
|
||||||
assert!(self.len() == other.len());
|
assert!(self.len() == other.len());
|
||||||
let mut res: N = ::zero();
|
let mut res: N = ::zero();
|
||||||
for i in (0us .. self.len()) {
|
for i in 0..self.len() {
|
||||||
res = res + unsafe { self.unsafe_at(i) * other.unsafe_at(i) };
|
res = res + unsafe { self.unsafe_at(i) * other.unsafe_at(i) };
|
||||||
}
|
}
|
||||||
res
|
res
|
||||||
|
@ -484,7 +484,7 @@ macro_rules! small_dvec_from_impl (
|
||||||
|
|
||||||
let mut at: [N; $dim] = [ $( $zeros, )* ];
|
let mut at: [N; $dim] = [ $( $zeros, )* ];
|
||||||
|
|
||||||
for i in (0us .. dim) {
|
for i in 0..dim {
|
||||||
at[i] = f(i);
|
at[i] = f(i);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
|
@ -417,7 +417,7 @@ macro_rules! diag_impl(
|
||||||
|
|
||||||
#[inline]
|
#[inline]
|
||||||
fn set_diag(&mut self, diag: &$tv<N>) {
|
fn set_diag(&mut self, diag: &$tv<N>) {
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
unsafe { self.unsafe_set((i, i), diag.unsafe_at(i)) }
|
unsafe { self.unsafe_set((i, i), diag.unsafe_at(i)) }
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -426,7 +426,7 @@ macro_rules! diag_impl(
|
||||||
fn diag(&self) -> $tv<N> {
|
fn diag(&self) -> $tv<N> {
|
||||||
let mut diag: $tv<N> = ::zero();
|
let mut diag: $tv<N> = ::zero();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
unsafe { diag.unsafe_set(i, self.unsafe_at((i, i))) }
|
unsafe { diag.unsafe_set(i, self.unsafe_at((i, i))) }
|
||||||
}
|
}
|
||||||
|
|
||||||
|
@ -445,12 +445,12 @@ macro_rules! mat_mul_mat_impl(
|
||||||
// careful! we need to comute other * self here (self is the rhs).
|
// careful! we need to comute other * self here (self is the rhs).
|
||||||
let mut res: $t<N> = ::zero();
|
let mut res: $t<N> = ::zero();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
let mut acc: N = ::zero();
|
let mut acc: N = ::zero();
|
||||||
|
|
||||||
unsafe {
|
unsafe {
|
||||||
for k in (0us .. $dim) {
|
for k in 0..$dim {
|
||||||
acc = acc + self.at_fast((i, k)) * right.at_fast((k, j));
|
acc = acc + self.at_fast((i, k)) * right.at_fast((k, j));
|
||||||
}
|
}
|
||||||
|
|
||||||
|
@ -474,8 +474,8 @@ macro_rules! vec_mul_mat_impl(
|
||||||
fn mul(self, right: $t<N>) -> $v<N> {
|
fn mul(self, right: $t<N>) -> $v<N> {
|
||||||
let mut res : $v<N> = $zero();
|
let mut res : $v<N> = $zero();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
unsafe {
|
unsafe {
|
||||||
let val = res.at_fast(i) + self.at_fast(j) * right.at_fast((j, i));
|
let val = res.at_fast(i) + self.at_fast(j) * right.at_fast((j, i));
|
||||||
res.set_fast(i, val)
|
res.set_fast(i, val)
|
||||||
|
@ -498,8 +498,8 @@ macro_rules! mat_mul_vec_impl(
|
||||||
fn mul(self, right: $v<N>) -> $v<N> {
|
fn mul(self, right: $v<N>) -> $v<N> {
|
||||||
let mut res : $v<N> = $zero();
|
let mut res : $v<N> = $zero();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
unsafe {
|
unsafe {
|
||||||
let val = res.at_fast(i) + self.at_fast((i, j)) * right.at_fast(j);
|
let val = res.at_fast(i) + self.at_fast((i, j)) * right.at_fast(j);
|
||||||
res.set_fast(i, val)
|
res.set_fast(i, val)
|
||||||
|
@ -544,7 +544,7 @@ macro_rules! inv_impl(
|
||||||
let mut res: $t<N> = ::one();
|
let mut res: $t<N> = ::one();
|
||||||
|
|
||||||
// inversion using Gauss-Jordan elimination
|
// inversion using Gauss-Jordan elimination
|
||||||
for k in (0us .. $dim) {
|
for k in 0..$dim {
|
||||||
// search a non-zero value on the k-th column
|
// search a non-zero value on the k-th column
|
||||||
// FIXME: would it be worth it to spend some more time searching for the
|
// FIXME: would it be worth it to spend some more time searching for the
|
||||||
// max instead?
|
// max instead?
|
||||||
|
@ -565,7 +565,7 @@ macro_rules! inv_impl(
|
||||||
|
|
||||||
// swap pivot line
|
// swap pivot line
|
||||||
if n0 != k {
|
if n0 != k {
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
self.swap((n0, j), (k, j));
|
self.swap((n0, j), (k, j));
|
||||||
res.swap((n0, j), (k, j));
|
res.swap((n0, j), (k, j));
|
||||||
}
|
}
|
||||||
|
@ -573,26 +573,26 @@ macro_rules! inv_impl(
|
||||||
|
|
||||||
let pivot = self.at((k, k));
|
let pivot = self.at((k, k));
|
||||||
|
|
||||||
for j in (k .. $dim) {
|
for j in k..$dim {
|
||||||
let selfval = self.at((k, j)) / pivot;
|
let selfval = self.at((k, j)) / pivot;
|
||||||
self.set((k, j), selfval);
|
self.set((k, j), selfval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
let resval = res.at((k, j)) / pivot;
|
let resval = res.at((k, j)) / pivot;
|
||||||
res.set((k, j), resval);
|
res.set((k, j), resval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for l in (0us .. $dim) {
|
for l in 0..$dim {
|
||||||
if l != k {
|
if l != k {
|
||||||
let normalizer = self.at((l, k));
|
let normalizer = self.at((l, k));
|
||||||
|
|
||||||
for j in (k .. $dim) {
|
for j in k..$dim {
|
||||||
let selfval = self.at((l, j)) - self.at((k, j)) * normalizer;
|
let selfval = self.at((l, j)) - self.at((k, j)) * normalizer;
|
||||||
self.set((l, j), selfval);
|
self.set((l, j), selfval);
|
||||||
}
|
}
|
||||||
|
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
let resval = res.at((l, j)) - res.at((k, j)) * normalizer;
|
let resval = res.at((l, j)) - res.at((k, j)) * normalizer;
|
||||||
res.set((l, j), resval);
|
res.set((l, j), resval);
|
||||||
}
|
}
|
||||||
|
@ -621,8 +621,8 @@ macro_rules! transpose_impl(
|
||||||
|
|
||||||
#[inline]
|
#[inline]
|
||||||
fn transpose_mut(&mut self) {
|
fn transpose_mut(&mut self) {
|
||||||
for i in (1us .. $dim) {
|
for i in 1..$dim {
|
||||||
for j in (0us .. i) {
|
for j in 0..i {
|
||||||
self.swap((i, j), (j, i))
|
self.swap((i, j), (j, i))
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -666,8 +666,8 @@ macro_rules! to_homogeneous_impl(
|
||||||
fn to_homogeneous(&self) -> $t2<N> {
|
fn to_homogeneous(&self) -> $t2<N> {
|
||||||
let mut res: $t2<N> = ::one();
|
let mut res: $t2<N> = ::one();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
for j in (0us .. $dim) {
|
for j in 0..$dim {
|
||||||
res.set((i, j), self.at((i, j)))
|
res.set((i, j), self.at((i, j)))
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -685,8 +685,8 @@ macro_rules! from_homogeneous_impl(
|
||||||
fn from(m: &$t2<N>) -> $t<N> {
|
fn from(m: &$t2<N>) -> $t<N> {
|
||||||
let mut res: $t<N> = ::one();
|
let mut res: $t<N> = ::one();
|
||||||
|
|
||||||
for i in (0us .. $dim2) {
|
for i in 0..$dim2 {
|
||||||
for j in (0us .. $dim2) {
|
for j in 0..$dim2 {
|
||||||
res.set((i, j), m.at((i, j)))
|
res.set((i, j), m.at((i, j)))
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -706,8 +706,8 @@ macro_rules! outer_impl(
|
||||||
#[inline]
|
#[inline]
|
||||||
fn outer(&self, other: &$t<N>) -> $m<N> {
|
fn outer(&self, other: &$t<N>) -> $m<N> {
|
||||||
let mut res: $m<N> = ::zero();
|
let mut res: $m<N> = ::zero();
|
||||||
for i in (0us .. Dim::dim(None::<$t<N>>)) {
|
for i in 0..Dim::dim(None::<$t<N>>) {
|
||||||
for j in (0us .. Dim::dim(None::<$t<N>>)) {
|
for j in 0..Dim::dim(None::<$t<N>>) {
|
||||||
res.set((i, j), self.at(i) * other.at(j))
|
res.set((i, j), self.at(i) * other.at(j))
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
|
@ -308,7 +308,7 @@ macro_rules! basis_impl(
|
||||||
impl<N: Copy + BaseFloat + ApproxEq<N>> Basis for $t<N> {
|
impl<N: Copy + BaseFloat + ApproxEq<N>> Basis for $t<N> {
|
||||||
#[inline]
|
#[inline]
|
||||||
fn canonical_basis<F: FnMut($t<N>) -> bool>(mut f: F) {
|
fn canonical_basis<F: FnMut($t<N>) -> bool>(mut f: F) {
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
if !f(Basis::canonical_basis_element(i).unwrap()) { return }
|
if !f(Basis::canonical_basis_element(i).unwrap()) { return }
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
@ -319,7 +319,7 @@ macro_rules! basis_impl(
|
||||||
// orthogonalization algorithm
|
// orthogonalization algorithm
|
||||||
let mut basis: Vec<$t<N>> = Vec::new();
|
let mut basis: Vec<$t<N>> = Vec::new();
|
||||||
|
|
||||||
for i in (0us .. $dim) {
|
for i in 0..$dim {
|
||||||
let mut basis_element : $t<N> = ::zero();
|
let mut basis_element : $t<N> = ::zero();
|
||||||
|
|
||||||
unsafe {
|
unsafe {
|
||||||
|
|
Loading…
Reference in New Issue