M-Labs
/
nalgebra
3
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Add the types: DVec1, .., DVec6. 

They are stack-allocated, dynamically sized vectors with a maximum size.
This is useful for slicing small matrices, without allocation.
This commit is contained in:
Sébastien Crozet 2014-08-16 12:16:26 +02:00
parent 7d6884c3df
commit 40c9915870
13 changed files with 693 additions and 432 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Cargo.toml
View File

@ -3,6 +3,6 @@ name = "nalgebra"
version = "0.1.0"
authors = [ "Sébastien Crozet <developer@crozet.re>" ] # FIXME: add the contributors.
[[lib]]
[lib]
name = "nalgebra"
path = "src/lib.rs"

6
src/linalg/decompositions.rs
View File

@ -34,11 +34,12 @@ pub fn householder_matrix<N: Float,
///
/// # Arguments
/// * `m` - matrix to decompose
pub fn decomp_qr<N: Float,
pub fn qr<N: Float,
V: Indexable<uint, N> + Norm<N>,
M: Clone + Eye + ColSlice<V> + Transpose
+ Indexable<(uint, uint), N> + Mul<M, M>>
(m: &M) -> (M, M) {
(m: &M)
-> (M, M) {
let (rows, cols) = m.shape();
assert!(rows >= cols);
let mut q : M = Eye::new_identity(rows);
@ -67,4 +68,3 @@ pub fn decomp_qr<N: Float,
(q, r)
}

2
src/linalg/mod.rs
View File

@ -1,4 +1,4 @@
pub use self::decompositions::{decomp_qr, householder_matrix};
pub use self::decompositions::{qr, eigen_qr, householder_matrix};
mod decompositions;

6
src/na.rs
View File

@ -48,7 +48,8 @@ pub use traits::{
pub use structs::{
Identity,
DMat, DVec,
DMat,
DVec, DVec1, DVec2, DVec3, DVec4, DVec5, DVec6,
Iso2, Iso3, Iso4,
Mat1, Mat2, Mat3, Mat4,
Mat5, Mat6,
@ -57,7 +58,8 @@ pub use structs::{
};
pub use linalg::{
decomp_qr,
qr,
eigen_qr,
householder_matrix
};

4
src/structs/dmat.rs
View File

@ -245,7 +245,7 @@ impl<N: Clone> Indexable<(uint, uint), N> for DMat<N> {
#[inline]
unsafe fn unsafe_at(&self, rowcol: (uint, uint)) -> N {
let (row, col) = rowcol;
(*self.mij.as_slice().unsafe_ref(self.offset(row, col))).clone()
(*self.mij.as_slice().unsafe_get(self.offset(row, col))).clone()
}
#[inline]
@ -524,7 +524,7 @@ impl<N: Clone> ColSlice<DVec<N>> for DMat<N> {
// we can init from slice thanks to the matrix being column major
let start= self.offset(row_start, col_id);
let stop = self.offset(row_end, col_id);
let slice = DVec::from_vec(
let slice = DVec::from_slice(
row_end - row_start, self.mij.slice(start, stop));
slice
}

406
src/structs/dvec.rs
View File

@ -11,97 +11,13 @@ use std::iter::FromIterator;
use traits::geometry::{Dot, Norm};
use traits::structure::{Iterable, IterableMut, Indexable};
/// Vector with a dimension unknown at compile-time.
/// Heap allocated, dynamically sized vector.
#[deriving(Eq, PartialEq, Show, Clone)]
pub struct DVec<N> {
/// Components of the vector. Contains as much elements as the vector dimension.
pub at: Vec<N>
}
double_dispatch_binop_decl_trait!(DVec, DVecMulRhs)
double_dispatch_binop_decl_trait!(DVec, DVecDivRhs)
double_dispatch_binop_decl_trait!(DVec, DVecAddRhs)
double_dispatch_binop_decl_trait!(DVec, DVecSubRhs)
mul_redispatch_impl!(DVec, DVecMulRhs)
div_redispatch_impl!(DVec, DVecDivRhs)
add_redispatch_impl!(DVec, DVecAddRhs)
sub_redispatch_impl!(DVec, DVecSubRhs)
impl<N: Zero + Clone> DVec<N> {
/// Builds a vector filled with zeros.
///
/// # Arguments
/// * `dim` - The dimension of the vector.
#[inline]
pub fn new_zeros(dim: uint) -> DVec<N> {
DVec::from_elem(dim, Zero::zero())
}
/// Tests if all components of the vector are zeroes.
#[inline]
pub fn is_zero(&self) -> bool {
self.at.iter().all(|e| e.is_zero())
}
}
impl<N: Clone> Indexable<uint, N> for DVec<N> {
fn at(&self, i: uint) -> N {
assert!(i < self.at.len());
unsafe {
self.unsafe_at(i)
}
}
fn set(&mut self, i: uint, val: N) {
assert!(i < self.at.len());
unsafe {
self.unsafe_set(i, val);
}
}
fn swap(&mut self, i: uint, j: uint) {
assert!(i < self.at.len());
assert!(j < self.at.len());
self.at.as_mut_slice().swap(i, j);
}
fn shape(&self) -> uint {
self.at.len()
}
#[inline]
unsafe fn unsafe_at(&self, i: uint) -> N {
(*self.at.as_slice().unsafe_ref(i)).clone()
}
#[inline]
unsafe fn unsafe_set(&mut self, i: uint, val: N) {
*self.at.as_mut_slice().unsafe_mut_ref(i) = val
}
}
impl<N: One + Clone> DVec<N> {
/// Builds a vector filled with ones.
///
/// # Arguments
/// * `dim` - The dimension of the vector.
#[inline]
pub fn new_ones(dim: uint) -> DVec<N> {
DVec::from_elem(dim, One::one())
}
}
impl<N: Rand> DVec<N> {
/// Builds a vector filled with random values.
#[inline]
pub fn new_random(dim: uint) -> DVec<N> {
DVec::from_fn(dim, |_| rand::random())
}
}
impl<N> DVec<N> {
/// Creates an uninitialized vec.
#[inline]
@ -113,24 +29,6 @@ impl<N> DVec<N> {
at: vec
}
}
/// Gets a reference to of this vector data.
#[inline]
pub fn as_vec<'r>(&'r self) -> &'r [N] {
self.at.as_slice()
}
/// Gets a mutable reference to of this vector data.
#[inline]
pub fn as_mut_vec<'r>(&'r mut self) -> &'r mut [N] {
self.at.as_mut_slice()
}
/// Extracts this vector data.
#[inline]
pub fn to_vec(self) -> Vec<N> {
self.at
}
}
impl<N: Clone> DVec<N> {
@ -144,7 +42,7 @@ impl<N: Clone> DVec<N> {
///
/// The vector must have at least `dim` elements.
#[inline]
pub fn from_vec(dim: uint, vec: &[N]) -> DVec<N> {
pub fn from_slice(dim: uint, vec: &[N]) -> DVec<N> {
assert!(dim <= vec.len());
DVec {
@ -161,27 +59,6 @@ impl<N> DVec<N> {
}
}
impl<N> Collection for DVec<N> {
#[inline]
fn len(&self) -> uint {
self.at.len()
}
}
impl<N> Iterable<N> for DVec<N> {
#[inline]
fn iter<'l>(&'l self) -> Items<'l, N> {
self.at.iter()
}
}
impl<N> IterableMut<N> for DVec<N> {
#[inline]
fn mut_iter<'l>(&'l mut self) -> MutItems<'l, N> {
self.at.mut_iter()
}
}
impl<N> FromIterator<N> for DVec<N> {
#[inline]
fn from_iter<I: Iterator<N>>(mut param: I) -> DVec<N> {
@ -195,258 +72,71 @@ impl<N> FromIterator<N> for DVec<N> {
}
}
impl<N: Clone + Float + ApproxEq<N> + DVecMulRhs<N, DVec<N>>> DVec<N> {
/// Computes the canonical basis for the given dimension. A canonical basis is a set of
/// vectors, mutually orthogonal, with all its component equal to 0.0 except one which is equal
/// to 1.0.
pub fn canonical_basis_with_dim(dim: uint) -> Vec<DVec<N>> {
let mut res : Vec<DVec<N>> = Vec::new();
for i in range(0u, dim) {
let mut basis_element : DVec<N> = DVec::new_zeros(dim);
*basis_element.at.get_mut(i) = One::one();
res.push(basis_element);
}
res
}
/// Computes a basis of the space orthogonal to the vector. If the input vector is of dimension
/// `n`, this will return `n - 1` vectors.
pub fn orthogonal_subspace_basis(&self) -> Vec<DVec<N>> {
// compute the basis of the orthogonal subspace using Gram-Schmidt
// orthogonalization algorithm
let dim = self.at.len();
let mut res : Vec<DVec<N>> = Vec::new();
for i in range(0u, dim) {
let mut basis_element : DVec<N> = DVec::new_zeros(self.at.len());
*basis_element.at.get_mut(i) = One::one();
if res.len() == dim - 1 {
break;
}
let mut elt = basis_element.clone();
elt = elt - self * Dot::dot(&basis_element, self);
for v in res.iter() {
elt = elt - v * Dot::dot(&elt, v)
};
if !ApproxEq::approx_eq(&Norm::sqnorm(&elt), &Zero::zero()) {
res.push(Norm::normalize_cpy(&elt));
}
}
assert!(res.len() == dim - 1);
res
impl<N> Collection for DVec<N> {
#[inline]
fn len(&self) -> uint {
self.at.len()
}
}
impl<N: Add<N, N>> DVecAddRhs<N, DVec<N>> for DVec<N> {
#[inline]
fn binop(left: &DVec<N>, right: &DVec<N>) -> DVec<N> {
assert!(left.at.len() == right.at.len());
DVec {
at: left.at.iter().zip(right.at.iter()).map(|(a, b)| *a + *b).collect()
}
}
dvec_impl!(DVec, DVecMulRhs, DVecDivRhs, DVecAddRhs, DVecSubRhs)
/// Stack-allocated, dynamically sized vector with a maximum size of 1.
pub struct DVec1<N> {
at: [N, ..1],
dim: uint
}
impl<N: Sub<N, N>> DVecSubRhs<N, DVec<N>> for DVec<N> {
#[inline]
fn binop(left: &DVec<N>, right: &DVec<N>) -> DVec<N> {
assert!(left.at.len() == right.at.len());
DVec {
at: left.at.iter().zip(right.at.iter()).map(|(a, b)| *a - *b).collect()
}
}
small_dvec_impl!(DVec1, 1, DVec1MulRhs, DVec1DivRhs, DVec1AddRhs, DVec1SubRhs, 0)
small_dvec_from_impl!(DVec1, 1, Zero::zero())
/// Stack-allocated, dynamically sized vector with a maximum size of 2.
pub struct DVec2<N> {
at: [N, ..2],
dim: uint
}
impl<N: Neg<N>> Neg<DVec<N>> for DVec<N> {
#[inline]
fn neg(&self) -> DVec<N> {
DVec { at: self.at.iter().map(|a| -*a).collect() }
}
small_dvec_impl!(DVec2, 2, DVec2MulRhs, DVec2DivRhs, DVec2AddRhs, DVec2SubRhs, 0, 1)
small_dvec_from_impl!(DVec2, 2, Zero::zero(), Zero::zero())
/// Stack-allocated, dynamically sized vector with a maximum size of 3.
pub struct DVec3<N> {
at: [N, ..3],
dim: uint
}
impl<N: Num + Clone> Dot<N> for DVec<N> {
#[inline]
fn dot(a: &DVec<N>, b: &DVec<N>) -> N {
assert!(a.at.len() == b.at.len());
small_dvec_impl!(DVec3, 3, DVec3MulRhs, DVec3DivRhs, DVec3AddRhs, DVec3SubRhs, 0, 1, 2)
small_dvec_from_impl!(DVec3, 3, Zero::zero(), Zero::zero(), Zero::zero())
let mut res: N = Zero::zero();
for i in range(0u, a.at.len()) {
res = res + unsafe { a.unsafe_at(i) * b.unsafe_at(i) };
}
res
}
#[inline]
fn sub_dot(a: &DVec<N>, b: &DVec<N>, c: &DVec<N>) -> N {
let mut res: N = Zero::zero();
for i in range(0u, a.at.len()) {
res = res + unsafe { (a.unsafe_at(i) - b.unsafe_at(i)) * c.unsafe_at(i) };
}
res
}
/// Stack-allocated, dynamically sized vector with a maximum size of 4.
pub struct DVec4<N> {
at: [N, ..4],
dim: uint
}
impl<N: Float + Clone> Norm<N> for DVec<N> {
#[inline]
fn sqnorm(v: &DVec<N>) -> N {
Dot::dot(v, v)
}
small_dvec_impl!(DVec4, 4, DVec4MulRhs, DVec4DivRhs, DVec4AddRhs, DVec4SubRhs, 0, 1, 2, 3)
small_dvec_from_impl!(DVec4, 4, Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero())
#[inline]
fn norm(v: &DVec<N>) -> N {
Norm::sqnorm(v).sqrt()
}
#[inline]
fn normalize_cpy(v: &DVec<N>) -> DVec<N> {
let mut res : DVec<N> = v.clone();
let _ = res.normalize();
res
}
#[inline]
fn normalize(&mut self) -> N {
let l = Norm::norm(self);
for i in range(0u, self.at.len()) {
*self.at.get_mut(i) = self.at[i] / l;
}
l
}
/// Stack-allocated, dynamically sized vector with a maximum size of 5.
pub struct DVec5<N> {
at: [N, ..5],
dim: uint
}
impl<N: ApproxEq<N>> ApproxEq<N> for DVec<N> {
#[inline]
fn approx_epsilon(_: Option<DVec<N>>) -> N {
ApproxEq::approx_epsilon(None::<N>)
}
small_dvec_impl!(DVec5, 5, DVec5MulRhs, DVec5DivRhs, DVec5AddRhs, DVec5SubRhs, 0, 1, 2, 3, 4)
small_dvec_from_impl!(DVec5, 5, Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero())
#[inline]
fn approx_eq(a: &DVec<N>, b: &DVec<N>) -> bool {
let mut zip = a.at.iter().zip(b.at.iter());
zip.all(|(a, b)| ApproxEq::approx_eq(a, b))
}
#[inline]
fn approx_eq_eps(a: &DVec<N>, b: &DVec<N>, epsilon: &N) -> bool {
let mut zip = a.at.iter().zip(b.at.iter());
zip.all(|(a, b)| ApproxEq::approx_eq_eps(a, b, epsilon))
}
/// Stack-allocated, dynamically sized vector with a maximum size of 6.
pub struct DVec6<N> {
at: [N, ..6],
dim: uint
}
macro_rules! scalar_mul_impl (
($n: ident) => (
impl DVecMulRhs<$n, DVec<$n>> for $n {
#[inline]
fn binop(left: &DVec<$n>, right: &$n) -> DVec<$n> {
DVec { at: left.at.iter().map(|a| a * *right).collect() }
}
}
)
)
macro_rules! scalar_div_impl (
($n: ident) => (
impl DVecDivRhs<$n, DVec<$n>> for $n {
#[inline]
fn binop(left: &DVec<$n>, right: &$n) -> DVec<$n> {
DVec { at: left.at.iter().map(|a| a / *right).collect() }
}
}
)
)
macro_rules! scalar_add_impl (
($n: ident) => (
impl DVecAddRhs<$n, DVec<$n>> for $n {
#[inline]
fn binop(left: &DVec<$n>, right: &$n) -> DVec<$n> {
DVec { at: left.at.iter().map(|a| a + *right).collect() }
}
}
)
)
macro_rules! scalar_sub_impl (
($n: ident) => (
impl DVecSubRhs<$n, DVec<$n>> for $n {
#[inline]
fn binop(left: &DVec<$n>, right: &$n) -> DVec<$n> {
DVec { at: left.at.iter().map(|a| a - *right).collect() }
}
}
)
)
scalar_mul_impl!(f64)
scalar_mul_impl!(f32)
scalar_mul_impl!(u64)
scalar_mul_impl!(u32)
scalar_mul_impl!(u16)
scalar_mul_impl!(u8)
scalar_mul_impl!(i64)
scalar_mul_impl!(i32)
scalar_mul_impl!(i16)
scalar_mul_impl!(i8)
scalar_mul_impl!(uint)
scalar_mul_impl!(int)
scalar_div_impl!(f64)
scalar_div_impl!(f32)
scalar_div_impl!(u64)
scalar_div_impl!(u32)
scalar_div_impl!(u16)
scalar_div_impl!(u8)
scalar_div_impl!(i64)
scalar_div_impl!(i32)
scalar_div_impl!(i16)
scalar_div_impl!(i8)
scalar_div_impl!(uint)
scalar_div_impl!(int)
scalar_add_impl!(f64)
scalar_add_impl!(f32)
scalar_add_impl!(u64)
scalar_add_impl!(u32)
scalar_add_impl!(u16)
scalar_add_impl!(u8)
scalar_add_impl!(i64)
scalar_add_impl!(i32)
scalar_add_impl!(i16)
scalar_add_impl!(i8)
scalar_add_impl!(uint)
scalar_add_impl!(int)
scalar_sub_impl!(f64)
scalar_sub_impl!(f32)
scalar_sub_impl!(u64)
scalar_sub_impl!(u32)
scalar_sub_impl!(u16)
scalar_sub_impl!(u8)
scalar_sub_impl!(i64)
scalar_sub_impl!(i32)
scalar_sub_impl!(i16)
scalar_sub_impl!(i8)
scalar_sub_impl!(uint)
scalar_sub_impl!(int)
small_dvec_impl!(DVec6, 6, DVec6MulRhs, DVec6DivRhs, DVec6AddRhs, DVec6SubRhs, 0, 1, 2, 3, 4, 5)
small_dvec_from_impl!(DVec6, 6, Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero(), Zero::zero())

502
src/structs/dvec_macros.rs Normal file
View File

@ -0,0 +1,502 @@
#![macro_escape]
macro_rules! dvec_impl(
($dvec: ident, $mul: ident, $div: ident, $add: ident, $sub: ident) => (
double_dispatch_binop_decl_trait!($dvec, $mul)
double_dispatch_binop_decl_trait!($dvec, $div)
double_dispatch_binop_decl_trait!($dvec, $add)
double_dispatch_binop_decl_trait!($dvec, $sub)
mul_redispatch_impl!($dvec, $mul)
div_redispatch_impl!($dvec, $div)
add_redispatch_impl!($dvec, $add)
sub_redispatch_impl!($dvec, $sub)
impl<N: Zero + Clone> $dvec<N> {
/// Builds a vector filled with zeros.
///
/// # Arguments
/// * `dim` - The dimension of the vector.
#[inline]
pub fn new_zeros(dim: uint) -> $dvec<N> {
$dvec::from_elem(dim, Zero::zero())
}
/// Tests if all components of the vector are zeroes.
#[inline]
pub fn is_zero(&self) -> bool {
self.as_slice().iter().all(|e| e.is_zero())
}
}
impl<N> $dvec<N> {
/// Slices this vector.
#[inline]
pub fn as_slice<'a>(&'a self) -> &'a [N] {
self.at.slice_to(self.len())
}
/// Mutably slices this vector.
#[inline]
pub fn as_mut_slice<'a>(&'a mut self) -> &'a mut [N] {
let len = self.len();
self.at.mut_slice_to(len)
}
}
impl<N: Clone> Indexable<uint, N> for $dvec<N> {
#[inline]
fn at(&self, i: uint) -> N {
assert!(i < self.len());
unsafe {
self.unsafe_at(i)
}
}
#[inline]
fn set(&mut self, i: uint, val: N) {
assert!(i < self.len());
unsafe {
self.unsafe_set(i, val);
}
}
#[inline]
fn swap(&mut self, i: uint, j: uint) {
assert!(i < self.len());
assert!(j < self.len());
self.as_mut_slice().swap(i, j);
}
#[inline]
fn shape(&self) -> uint {
self.len()
}
#[inline]
unsafe fn unsafe_at(&self, i: uint) -> N {
(*self.at.as_slice().unsafe_get(i)).clone()
}
#[inline]
unsafe fn unsafe_set(&mut self, i: uint, val: N) {
*self.at.as_mut_slice().unsafe_mut_ref(i) = val
}
}
impl<N: One + Zero + Clone> $dvec<N> {
/// Builds a vector filled with ones.
///
/// # Arguments
/// * `dim` - The dimension of the vector.
#[inline]
pub fn new_ones(dim: uint) -> $dvec<N> {
$dvec::from_elem(dim, One::one())
}
}
impl<N: Rand + Zero> $dvec<N> {
/// Builds a vector filled with random values.
#[inline]
pub fn new_random(dim: uint) -> $dvec<N> {
$dvec::from_fn(dim, |_| rand::random())
}
}
impl<N> Iterable<N> for $dvec<N> {
#[inline]
fn iter<'l>(&'l self) -> Items<'l, N> {
self.as_slice().iter()
}
}
impl<N> IterableMut<N> for $dvec<N> {
#[inline]
fn mut_iter<'l>(&'l mut self) -> MutItems<'l, N> {
self.as_mut_slice().mut_iter()
}
}
impl<N: Clone + Float + ApproxEq<N> + $mul<N, $dvec<N>>> $dvec<N> {
/// Computes the canonical basis for the given dimension. A canonical basis is a set of
/// vectors, mutually orthogonal, with all its component equal to 0.0 except one which is equal
/// to 1.0.
pub fn canonical_basis_with_dim(dim: uint) -> Vec<$dvec<N>> {
let mut res : Vec<$dvec<N>> = Vec::new();
for i in range(0u, dim) {
let mut basis_element : $dvec<N> = $dvec::new_zeros(dim);
basis_element.set(i, One::one());
res.push(basis_element);
}
res
}
/// Computes a basis of the space orthogonal to the vector. If the input vector is of dimension
/// `n`, this will return `n - 1` vectors.
pub fn orthogonal_subspace_basis(&self) -> Vec<$dvec<N>> {
// compute the basis of the orthogonal subspace using Gram-Schmidt
// orthogonalization algorithm
let dim = self.len();
let mut res : Vec<$dvec<N>> = Vec::new();
for i in range(0u, dim) {
let mut basis_element : $dvec<N> = $dvec::new_zeros(self.len());
basis_element.set(i, One::one());
if res.len() == dim - 1 {
break;
}
let mut elt = basis_element.clone();
elt = elt - self * Dot::dot(&basis_element, self);
for v in res.iter() {
elt = elt - v * Dot::dot(&elt, v)
};
if !ApproxEq::approx_eq(&Norm::sqnorm(&elt), &Zero::zero()) {
res.push(Norm::normalize_cpy(&elt));
}
}
assert!(res.len() == dim - 1);
res
}
}
impl<N: Add<N, N> + Zero> $add<N, $dvec<N>> for $dvec<N> {
#[inline]
fn binop(left: &$dvec<N>, right: &$dvec<N>) -> $dvec<N> {
assert!(left.len() == right.len());
FromIterator::from_iter(left.as_slice().iter().zip(right.as_slice().iter()).map(|(a, b)| *a + *b))
}
}
impl<N: Sub<N, N> + Zero> $sub<N, $dvec<N>> for $dvec<N> {
#[inline]
fn binop(left: &$dvec<N>, right: &$dvec<N>) -> $dvec<N> {
assert!(left.len() == right.len());
FromIterator::from_iter(left.as_slice().iter().zip(right.as_slice().iter()).map(|(a, b)| *a - *b))
}
}
impl<N: Neg<N> + Zero> Neg<$dvec<N>> for $dvec<N> {
#[inline]
fn neg(&self) -> $dvec<N> {
FromIterator::from_iter(self.as_slice().iter().map(|a| -*a))
}
}
impl<N: Num + Clone> Dot<N> for $dvec<N> {
#[inline]
fn dot(a: &$dvec<N>, b: &$dvec<N>) -> N {
assert!(a.len() == b.len());
let mut res: N = Zero::zero();
for i in range(0u, a.len()) {
res = res + unsafe { a.unsafe_at(i) * b.unsafe_at(i) };
}
res
}
#[inline]
fn sub_dot(a: &$dvec<N>, b: &$dvec<N>, c: &$dvec<N>) -> N {
let mut res: N = Zero::zero();
for i in range(0u, a.len()) {
res = res + unsafe { (a.unsafe_at(i) - b.unsafe_at(i)) * c.unsafe_at(i) };
}
res
}
}
impl<N: Float + Clone> Norm<N> for $dvec<N> {
#[inline]
fn sqnorm(v: &$dvec<N>) -> N {
Dot::dot(v, v)
}
#[inline]
fn norm(v: &$dvec<N>) -> N {
Norm::sqnorm(v).sqrt()
}
#[inline]
fn normalize_cpy(v: &$dvec<N>) -> $dvec<N> {
let mut res : $dvec<N> = v.clone();
let _ = res.normalize();
res
}
#[inline]
fn normalize(&mut self) -> N {
let l = Norm::norm(self);
for n in self.as_mut_slice().mut_iter() {
*n = *n / l;
}
l
}
}
impl<N: ApproxEq<N>> ApproxEq<N> for $dvec<N> {
#[inline]
fn approx_epsilon(_: Option<$dvec<N>>) -> N {
ApproxEq::approx_epsilon(None::<N>)
}
#[inline]
fn approx_eq(a: &$dvec<N>, b: &$dvec<N>) -> bool {
let mut zip = a.as_slice().iter().zip(b.as_slice().iter());
zip.all(|(a, b)| ApproxEq::approx_eq(a, b))
}
#[inline]
fn approx_eq_eps(a: &$dvec<N>, b: &$dvec<N>, epsilon: &N) -> bool {
let mut zip = a.as_slice().iter().zip(b.as_slice().iter());
zip.all(|(a, b)| ApproxEq::approx_eq_eps(a, b, epsilon))
}
}
dvec_scalar_mul_impl!($dvec, f64, $mul)
dvec_scalar_mul_impl!($dvec, f32, $mul)
dvec_scalar_mul_impl!($dvec, u64, $mul)
dvec_scalar_mul_impl!($dvec, u32, $mul)
dvec_scalar_mul_impl!($dvec, u16, $mul)
dvec_scalar_mul_impl!($dvec, u8, $mul)
dvec_scalar_mul_impl!($dvec, i64, $mul)
dvec_scalar_mul_impl!($dvec, i32, $mul)
dvec_scalar_mul_impl!($dvec, i16, $mul)
dvec_scalar_mul_impl!($dvec, i8, $mul)
dvec_scalar_mul_impl!($dvec, uint, $mul)
dvec_scalar_mul_impl!($dvec, int, $mul)
dvec_scalar_div_impl!($dvec, f64, $div)
dvec_scalar_div_impl!($dvec, f32, $div)
dvec_scalar_div_impl!($dvec, u64, $div)
dvec_scalar_div_impl!($dvec, u32, $div)
dvec_scalar_div_impl!($dvec, u16, $div)
dvec_scalar_div_impl!($dvec, u8, $div)
dvec_scalar_div_impl!($dvec, i64, $div)
dvec_scalar_div_impl!($dvec, i32, $div)
dvec_scalar_div_impl!($dvec, i16, $div)
dvec_scalar_div_impl!($dvec, i8, $div)
dvec_scalar_div_impl!($dvec, uint, $div)
dvec_scalar_div_impl!($dvec, int, $div)
dvec_scalar_add_impl!($dvec, f64, $add)
dvec_scalar_add_impl!($dvec, f32, $add)
dvec_scalar_add_impl!($dvec, u64, $add)
dvec_scalar_add_impl!($dvec, u32, $add)
dvec_scalar_add_impl!($dvec, u16, $add)
dvec_scalar_add_impl!($dvec, u8, $add)
dvec_scalar_add_impl!($dvec, i64, $add)
dvec_scalar_add_impl!($dvec, i32, $add)
dvec_scalar_add_impl!($dvec, i16, $add)
dvec_scalar_add_impl!($dvec, i8, $add)
dvec_scalar_add_impl!($dvec, uint, $add)
dvec_scalar_add_impl!($dvec, int, $add)
dvec_scalar_sub_impl!($dvec, f64, $sub)
dvec_scalar_sub_impl!($dvec, f32, $sub)
dvec_scalar_sub_impl!($dvec, u64, $sub)
dvec_scalar_sub_impl!($dvec, u32, $sub)
dvec_scalar_sub_impl!($dvec, u16, $sub)
dvec_scalar_sub_impl!($dvec, u8, $sub)
dvec_scalar_sub_impl!($dvec, i64, $sub)
dvec_scalar_sub_impl!($dvec, i32, $sub)
dvec_scalar_sub_impl!($dvec, i16, $sub)
dvec_scalar_sub_impl!($dvec, i8, $sub)
dvec_scalar_sub_impl!($dvec, uint, $sub)
dvec_scalar_sub_impl!($dvec, int, $sub)
)
)
macro_rules! dvec_scalar_mul_impl (
($dvec: ident, $n: ident, $mul: ident) => (
impl $mul<$n, $dvec<$n>> for $n {
#[inline]
fn binop(left: &$dvec<$n>, right: &$n) -> $dvec<$n> {
FromIterator::from_iter(left.as_slice().iter().map(|a| a * *right))
}
}
)
)
macro_rules! dvec_scalar_div_impl (
($dvec: ident, $n: ident, $div: ident) => (
impl $div<$n, $dvec<$n>> for $n {
#[inline]
fn binop(left: &$dvec<$n>, right: &$n) -> $dvec<$n> {
FromIterator::from_iter(left.as_slice().iter().map(|a| a / *right))
}
}
)
)
macro_rules! dvec_scalar_add_impl (
($dvec: ident, $n: ident, $add: ident) => (
impl $add<$n, $dvec<$n>> for $n {
#[inline]
fn binop(left: &$dvec<$n>, right: &$n) -> $dvec<$n> {
FromIterator::from_iter(left.as_slice().iter().map(|a| a + *right))
}
}
)
)
macro_rules! dvec_scalar_sub_impl (
($dvec: ident, $n: ident, $sub: ident) => (
impl $sub<$n, $dvec<$n>> for $n {
#[inline]
fn binop(left: &$dvec<$n>, right: &$n) -> $dvec<$n> {
FromIterator::from_iter(left.as_slice().iter().map(|a| a - *right))
}
}
)
)
macro_rules! small_dvec_impl (
($dvec: ident, $dim: expr, $mul: ident, $div: ident, $add: ident, $sub: ident $(,$idx: expr)*) => (
impl<N> Collection for $dvec<N> {
#[inline]
fn len(&self) -> uint {
self.dim
}
}
impl<N: PartialEq> PartialEq for $dvec<N> {
#[inline]
fn eq(&self, other: &$dvec<N>) -> bool {
if self.len() != other.len() {
return false; // FIXME: fail instead?
}
for (a, b) in self.as_slice().iter().zip(other.as_slice().iter()) {
if *a != *b {
return false;
}
}
true
}
}
impl<N: Clone> Clone for $dvec<N> {
fn clone(&self) -> $dvec<N> {
let at: [N, ..$dim] = [ $( self.at[$idx].clone(), )* ];
$dvec {
at: at,
dim: self.dim
}
}
}
dvec_impl!($dvec, $mul, $div, $add, $sub)
)
)
macro_rules! small_dvec_from_impl (
($dvec: ident, $dim: expr $(,$zeros: expr)*) => (
impl<N: Clone + Zero> $dvec<N> {
/// Builds a vector filled with a constant.
#[inline]
pub fn from_elem(dim: uint, elem: N) -> $dvec<N> {
assert!(dim <= $dim);
let mut at: [N, ..$dim] = [ $( $zeros, )* ];
for n in at.mut_slice_to(dim).mut_iter() {
*n = elem.clone();
}
$dvec {
at: at,
dim: dim
}
}
}
impl<N: Clone + Zero> $dvec<N> {
/// Builds a vector filled with the components provided by a vector.
///
/// The vector must have at least `dim` elements.
#[inline]
pub fn from_slice(dim: uint, vec: &[N]) -> $dvec<N> {
assert!(dim <= vec.len() && dim <= $dim);
// FIXME: not safe.
let mut at: [N, ..$dim] = [ $( $zeros, )* ];
for (curr, other) in vec.iter().zip(at.mut_iter()) {
*other = curr.clone();
}
$dvec {
at: at,
dim: dim
}
}
}
impl<N: Zero> $dvec<N> {
/// Builds a vector filled with the result of a function.
#[inline(always)]
pub fn from_fn(dim: uint, f: |uint| -> N) -> $dvec<N> {
assert!(dim <= $dim);
let mut at: [N, ..$dim] = [ $( $zeros, )* ];
for i in range(0, dim) {
at[i] = f(i);
}
$dvec {
at: at,
dim: dim
}
}
}
impl<N: Zero> FromIterator<N> for $dvec<N> {
#[inline]
fn from_iter<I: Iterator<N>>(mut param: I) -> $dvec<N> {
let mut at: [N, ..$dim] = [ $( $zeros, )* ];
let mut dim = 0;
for n in param {
if dim == $dim {
break;
}
at[dim] = n;
dim = dim + 1;
}
$dvec {
at: at,
dim: dim
}
}
}
)
)

30
src/structs/mat.rs
View File

@ -6,9 +6,9 @@ use std::mem;
use std::num::{One, Zero};
use traits::operations::ApproxEq;
use std::slice::{Items, MutItems};
use structs::vec::{Vec1, Vec2, Vec3, Vec4, Vec5, Vec6, Vec1MulRhs, Vec4MulRhs,
Vec5MulRhs, Vec6MulRhs};
use structs::dvec::DVec;
use structs::vec::{Vec1, Vec2, Vec3, Vec4, Vec5, Vec6,
Vec1MulRhs, Vec4MulRhs, Vec5MulRhs, Vec6MulRhs};
use structs::dvec::{DVec1, DVec2, DVec3, DVec4, DVec5, DVec6};
use traits::structure::{Cast, Row, Col, Iterable, IterableMut, Dim, Indexable,
Eye, ColSlice, RowSlice};
@ -118,8 +118,8 @@ transpose_impl!(Mat1, 1)
approx_eq_impl!(Mat1)
row_impl!(Mat1, Vec1, 1)
col_impl!(Mat1, Vec1, 1)
col_slice_impl!(Mat1, Vec1, 1)
row_slice_impl!(Mat1, Vec1, 1)
col_slice_impl!(Mat1, Vec1, DVec1, 1)
row_slice_impl!(Mat1, Vec1, DVec1, 1)
to_homogeneous_impl!(Mat1, Mat2, 1, 2)
from_homogeneous_impl!(Mat1, Mat2, 1, 2)
outer_impl!(Vec1, Mat1)
@ -219,8 +219,8 @@ transpose_impl!(Mat2, 2)
approx_eq_impl!(Mat2)
row_impl!(Mat2, Vec2, 2)
col_impl!(Mat2, Vec2, 2)
col_slice_impl!(Mat2, Vec2, 2)
row_slice_impl!(Mat2, Vec2, 2)
col_slice_impl!(Mat2, Vec2, DVec2, 2)
row_slice_impl!(Mat2, Vec2, DVec2, 2)
to_homogeneous_impl!(Mat2, Mat3, 2, 3)
from_homogeneous_impl!(Mat2, Mat3, 2, 3)
outer_impl!(Vec2, Mat2)
@ -334,8 +334,8 @@ transpose_impl!(Mat3, 3)
approx_eq_impl!(Mat3)
// (specialized) row_impl!(Mat3, Vec3, 3)
// (specialized) col_impl!(Mat3, Vec3, 3)
col_slice_impl!(Mat3, Vec3, 3)
row_slice_impl!(Mat3, Vec3, 3)
col_slice_impl!(Mat3, Vec3, DVec3, 3)
row_slice_impl!(Mat3, Vec3, DVec3, 3)
to_homogeneous_impl!(Mat3, Mat4, 3, 4)
from_homogeneous_impl!(Mat3, Mat4, 3, 4)
outer_impl!(Vec3, Mat3)
@ -501,8 +501,8 @@ transpose_impl!(Mat4, 4)
approx_eq_impl!(Mat4)
row_impl!(Mat4, Vec4, 4)
col_impl!(Mat4, Vec4, 4)
col_slice_impl!(Mat4, Vec4, 4)
row_slice_impl!(Mat4, Vec4, 4)
col_slice_impl!(Mat4, Vec4, DVec4, 4)
row_slice_impl!(Mat4, Vec4, DVec4, 4)
to_homogeneous_impl!(Mat4, Mat5, 4, 5)
from_homogeneous_impl!(Mat4, Mat5, 4, 5)
outer_impl!(Vec4, Mat4)
@ -684,8 +684,8 @@ transpose_impl!(Mat5, 5)
approx_eq_impl!(Mat5)
row_impl!(Mat5, Vec5, 5)
col_impl!(Mat5, Vec5, 5)
col_slice_impl!(Mat5, Vec5, 5)
row_slice_impl!(Mat5, Vec5, 5)
col_slice_impl!(Mat5, Vec5, DVec5, 5)
row_slice_impl!(Mat5, Vec5, DVec5, 5)
to_homogeneous_impl!(Mat5, Mat6, 5, 6)
from_homogeneous_impl!(Mat5, Mat6, 5, 6)
outer_impl!(Vec5, Mat5)
@ -919,6 +919,6 @@ transpose_impl!(Mat6, 6)
approx_eq_impl!(Mat6)
row_impl!(Mat6, Vec6, 6)
col_impl!(Mat6, Vec6, 6)
col_slice_impl!(Mat6, Vec6, 6)
row_slice_impl!(Mat6, Vec6, 6)
col_slice_impl!(Mat6, Vec6, DVec6, 6)
row_slice_impl!(Mat6, Vec6, DVec6, 6)
outer_impl!(Vec6, Mat6)

36
src/structs/mat_macros.rs
View File

@ -226,20 +226,12 @@ macro_rules! indexable_impl(
)
macro_rules! col_slice_impl(
($t: ident, $tv: ident, $dim: expr) => (
impl<N: Clone + Zero> ColSlice<DVec<N>> for $t<N> {
fn col_slice(& self, cid: uint, rstart: uint, rend: uint) -> DVec<N> {
assert!(cid < $dim);
assert!(rstart < rend);
assert!(rend <= $dim);
($t: ident, $tv: ident, $slice: ident, $dim: expr) => (
impl<N: Clone + Zero> ColSlice<$slice<N>> for $t<N> {
fn col_slice(&self, cid: uint, rstart: uint, rend: uint) -> $slice<N> {
let col = self.col(cid);
let res = DVec::from_vec(
rend - rstart,
unsafe {
mem::transmute::<&$tv<N>, & [N, ..$dim]>
(&col).slice(rstart, rend)
});
res
$slice::from_slice(rend - rstart, col.as_slice().slice(rstart, rend))
}
}
)
@ -275,20 +267,12 @@ macro_rules! row_impl(
)
macro_rules! row_slice_impl(
($t: ident, $tv: ident, $dim: expr) => (
impl<N: Clone + Zero> RowSlice<DVec<N>> for $t<N> {
fn row_slice(& self, rid: uint, cstart: uint, cend: uint) -> DVec<N> {
assert!(rid < $dim);
assert!(cstart < cend);
assert!(cend <= $dim);
($t: ident, $tv: ident, $slice: ident, $dim: expr) => (
impl<N: Clone + Zero> RowSlice<$slice<N>> for $t<N> {
fn row_slice(&self, rid: uint, cstart: uint, cend: uint) -> $slice<N> {
let row = self.row(rid);
let res = DVec::from_vec(
cend - cstart,
unsafe {
mem::transmute::<&$tv<N>, & [N, ..$dim]>
(&row).slice(cstart, cend)
});
res
$slice::from_slice(cend - cstart, row.as_slice().slice(cstart, cend))
}
}
)

3
src/structs/mod.rs
View File

@ -1,7 +1,7 @@
//! Data structures and implementations.
pub use self::dmat::DMat;
pub use self::dvec::DVec;
pub use self::dvec::{DVec, DVec1, DVec2, DVec3, DVec4, DVec5, DVec6};
pub use self::vec::{Vec0, Vec1, Vec2, Vec3, Vec4, Vec5, Vec6};
pub use self::mat::{Identity, Mat1, Mat2, Mat3, Mat4, Mat5, Mat6};
pub use self::rot::{Rot2, Rot3, Rot4};
@ -18,6 +18,7 @@ pub use self::mat::{Mat1MulRhs, Mat2MulRhs, Mat3MulRhs, Mat4MulRhs, Mat5MulRhs,
mod metal;
mod dmat;
mod dvec_macros;
mod dvec;
mod vec_macros;
mod vec;

6
src/structs/vec.rs
View File

@ -38,6 +38,7 @@ new_impl!(Vec1, x)
ord_impl!(Vec1, x)
vec_axis_impl!(Vec1, x)
vec_cast_impl!(Vec1, Vec1Cast, x)
as_slice_impl!(Vec1, 1)
indexable_impl!(Vec1, 1)
at_fast_impl!(Vec1, 1)
new_repeat_impl!(Vec1, val, x)
@ -135,6 +136,7 @@ new_impl!(Vec2, x, y)
ord_impl!(Vec2, x, y)
vec_axis_impl!(Vec2, x, y)
vec_cast_impl!(Vec2, Vec2Cast, x, y)
as_slice_impl!(Vec2, 2)
indexable_impl!(Vec2, 2)
at_fast_impl!(Vec2, 2)
new_repeat_impl!(Vec2, val, x, y)
@ -234,6 +236,7 @@ new_impl!(Vec3, x, y, z)
ord_impl!(Vec3, x, y, z)
vec_axis_impl!(Vec3, x, y, z)
vec_cast_impl!(Vec3, Vec3Cast, x, y, z)
as_slice_impl!(Vec3, 3)
indexable_impl!(Vec3, 3)
at_fast_impl!(Vec3, 3)
new_repeat_impl!(Vec3, val, x, y, z)
@ -339,6 +342,7 @@ new_impl!(Vec4, x, y, z, w)
ord_impl!(Vec4, x, y, z, w)
vec_axis_impl!(Vec4, x, y, z, w)
vec_cast_impl!(Vec4, Vec4Cast, x, y, z, w)
as_slice_impl!(Vec4, 4)
indexable_impl!(Vec4, 4)
at_fast_impl!(Vec4, 4)
new_repeat_impl!(Vec4, val, x, y, z, w)
@ -442,6 +446,7 @@ new_impl!(Vec5, x, y, z, w, a)
ord_impl!(Vec5, x, y, z, w, a)
vec_axis_impl!(Vec5, x, y, z, w, a)
vec_cast_impl!(Vec5, Vec5Cast, x, y, z, w, a)
as_slice_impl!(Vec5, 5)
indexable_impl!(Vec5, 5)
at_fast_impl!(Vec5, 5)
new_repeat_impl!(Vec5, val, x, y, z, w, a)
@ -547,6 +552,7 @@ new_impl!(Vec6, x, y, z, w, a, b)
ord_impl!(Vec6, x, y, z, w, a, b)
vec_axis_impl!(Vec6, x, y, z, w, a, b)
vec_cast_impl!(Vec6, Vec6Cast, x, y, z, w, a, b)
as_slice_impl!(Vec6, 6)
indexable_impl!(Vec6, 6)
at_fast_impl!(Vec6, 6)
new_repeat_impl!(Vec6, val, x, y, z, w, a, b)

20
src/structs/vec_macros.rs
View File

@ -15,6 +15,26 @@ macro_rules! new_impl(
)
)
macro_rules! as_slice_impl(
($t: ident, $dim: expr) => (
impl<N> $t<N> {
/// Slices this vector.
pub fn as_slice<'a>(&'a self) -> &'a [N] {
unsafe {
mem::transmute::<&$t<N>, &[N, ..$dim]>(self).as_slice()
}
}
/// Mutably slices this vector.
pub fn as_mut_slice<'a>(&'a mut self) -> &'a mut [N] {
unsafe {
mem::transmute::<&mut $t<N>, &mut [N, ..$dim]>(self).as_mut_slice()
}
}
}
)
)
macro_rules! at_fast_impl(
($t: ident, $dim: expr) => (
impl<N: Clone> $t<N> {

96
src/tests/mat.rs
View File

@ -2,7 +2,6 @@ use std::num::{Float, abs};
use std::rand::random;
use na::{Vec1, Vec3, Mat1, Mat2, Mat3, Mat4, Mat5, Mat6, Rot3, DMat, DVec, Indexable};
use na;
use na::decomp_qr;
use std::cmp::{min, max};
macro_rules! test_inv_mat_impl(
@ -25,12 +24,12 @@ macro_rules! test_transpose_mat_impl(
);
)
macro_rules! test_decomp_qr_impl(
macro_rules! test_qr_impl(
($t: ty) => (
for _ in range(0u, 10000) {
let randmat : $t = random();
let (q, r) = decomp_qr(&randmat);
let (q, r) = na::qr(&randmat);
let recomp = q * r;
assert!(na::approx_eq(&randmat, &recomp));
@ -38,6 +37,35 @@ macro_rules! test_decomp_qr_impl(
);
)
macro_rules! test_eigen_qr_impl(
($t: ty) => {
for _ in range(0u, 10000) {
let randmat : $t = random();
let (eigenvalues, eigenvectors) = na::eigen_qr(&randmat, &Float::epsilon(), 1000);
// FIXME: provide a method to initialize a matrix from its diagonal!
let diag: $t = na::zero();
for i in range(0, na::dim::<$t>()) {
diag.set((i, i), eigenvalues.at(i));
}
let recomp = na::transpose(&eigenvectors) * diag * eigenvectors;
println!("mat: {}", randmat);
println!("eigenvectors: {}", eigenvectors);
println!("eigenvalues: {}", eigenvalues);
println!("recomp: {}", recomp);
assert!(false);
fail!("what!");
assert!(na::approx_eq(&randmat, &recomp));
}
}
)
#[test]
fn test_transpose_mat1() {
test_transpose_mat_impl!(Mat1<f64>);
@ -139,7 +167,7 @@ fn test_mean_dmat() {
]
);
assert!(na::approx_eq(&na::mean(&mat), &DVec::from_vec(3, [4.0f64, 5.0, 6.0])));
assert!(na::approx_eq(&na::mean(&mat), &DVec::from_slice(3, [4.0f64, 5.0, 6.0])));
}
#[test]
@ -223,14 +251,14 @@ fn test_dmat_from_vec() {
}
#[test]
fn test_decomp_qr() {
fn test_qr() {
for _ in range(0u, 10) {
let dim1: uint = random();
let dim2: uint = random();
let rows = min(40, max(dim1, dim2));
let cols = min(40, min(dim1, dim2));
let randmat: DMat<f64> = DMat::new_random(rows, cols);
let (q, r) = decomp_qr(&randmat);
let (q, r) = na::qr(&randmat);
let recomp = q * r;
assert!(na::approx_eq(&randmat, &recomp));
@ -238,40 +266,68 @@ fn test_decomp_qr() {
}
#[test]
fn test_decomp_qr_mat1() {
test_decomp_qr_impl!(Mat1<f64>);
fn test_qr_mat1() {
test_qr_impl!(Mat1<f64>);
}
#[test]
fn test_decomp_qr_mat2() {
test_decomp_qr_impl!(Mat2<f64>);
fn test_qr_mat2() {
test_qr_impl!(Mat2<f64>);
}
#[test]
fn test_decomp_qr_mat3() {
test_decomp_qr_impl!(Mat3<f64>);
fn test_qr_mat3() {
test_qr_impl!(Mat3<f64>);
}
#[test]
fn test_decomp_qr_mat4() {
test_decomp_qr_impl!(Mat4<f64>);
fn test_qr_mat4() {
test_qr_impl!(Mat4<f64>);
}
#[test]
fn test_decomp_qr_mat5() {
test_decomp_qr_impl!(Mat5<f64>);
fn test_qr_mat5() {
test_qr_impl!(Mat5<f64>);
}
#[test]
fn test_decomp_qr_mat6() {
test_decomp_qr_impl!(Mat6<f64>);
fn test_qr_mat6() {
test_qr_impl!(Mat6<f64>);
}
// #[test]
// fn test_eigen_qr_mat1() {
// test_eigen_qr_impl!(Mat1<f64>);
// }
//
// #[test]
// fn test_eigen_qr_mat2() {
// test_eigen_qr_impl!(Mat2<f64>);
// }
//
// #[test]
// fn test_eigen_qr_mat3() {
// test_eigen_qr_impl!(Mat3<f64>);
// }
//
// #[test]
// fn test_eigen_qr_mat4() {
// test_eigen_qr_impl!(Mat4<f64>);
// }
//
// #[test]
// fn test_eigen_qr_mat5() {
// test_eigen_qr_impl!(Mat5<f64>);
// }
//
// #[test]
// fn test_eigen_qr_mat6() {
// test_eigen_qr_impl!(Mat6<f64>);
// }
#[test]
fn test_from_fn() {
let actual: DMat<uint> = DMat::from_fn( 3, 4,
|i,j| 10*i + j);
let actual: DMat<uint> = DMat::from_fn(3, 4, |i, j| 10 * i + j);
let expected: DMat<uint> = DMat::from_row_vec(3, 4,
[0_0, 0_1, 0_2, 0_3,
1_0, 1_1, 1_2, 1_3,