Add tests and basis generation.

This commit is contained in:
Sébastien Crozet 2013-05-18 17:04:03 +00:00
parent 39707b42dc
commit 890cdb73f2
14 changed files with 617 additions and 28 deletions

View File

@ -1,7 +1,11 @@
nalgebra_lib_path=lib
all:
rust build src/nalgebra.rc --out-dir lib # rustpkg install
rust build src/nalgebra.rc --out-dir $(nalgebra_lib_path)
test:
rust test src/nalgebra.rc
doc:
rust doc src/nalgebra.rc --output-dir doc
rust test src/nalgebra.rc
.PHONY:doc
.PHONY:doc, test

View File

@ -1,8 +1,10 @@
use core::num::{Zero, Algebraic};
use core::num::{Zero, One, Algebraic};
use core::rand::{Rand, Rng, RngUtil};
use std::cmp::FuzzyEq;
use traits::dot::Dot;
use traits::dim::Dim;
use traits::basis::Basis;
use traits::norm::Norm;
#[deriving(Eq)]
pub struct Vec1<T>
@ -33,12 +35,28 @@ impl<T:Copy + Mul<T, T> + Add<T, T> + Algebraic> Dot<T> for Vec1<T>
{
fn dot(&self, other : &Vec1<T>) -> T
{ self.x * other.x }
}
impl<T:Copy + Mul<T, T> + Add<T, T> + Quot<T, T> + Algebraic>
Norm<T> for Vec1<T>
{
fn sqnorm(&self) -> T
{ self.dot(self) }
fn norm(&self) -> T
{ self.sqnorm().sqrt() }
fn normalized(&self) -> Vec1<T>
{ Vec1(self.x / self.norm()) }
fn normalize(&mut self) -> T
{
let l = self.norm();
self.x /= l;
l
}
}
impl<T:Copy + Neg<T>> Neg<Vec1<T>> for Vec1<T>
@ -59,6 +77,15 @@ impl<T:Copy + Zero> Zero for Vec1<T>
{ self.x.is_zero() }
}
impl<T: Copy + One> Basis for Vec1<T>
{
fn canonical_basis() -> ~[Vec1<T>]
{ ~[ Vec1(One::one()) ] } // FIXME: this should be static
fn orthogonal_subspace_basis(&self) -> ~[Vec1<T>]
{ ~[] }
}
impl<T:FuzzyEq<T>> FuzzyEq<T> for Vec1<T>
{
fn fuzzy_eq(&self, other: &Vec1<T>) -> bool

View File

@ -1,9 +1,11 @@
use core::num::{Zero, Algebraic};
use core::num::{Zero, One, Algebraic};
use core::rand::{Rand, Rng, RngUtil};
use std::cmp::FuzzyEq;
use traits::dot::Dot;
use traits::dim::Dim;
use traits::cross::Cross;
use traits::basis::Basis;
use traits::norm::Norm;
use dim1::vec1::Vec1;
#[deriving(Eq)]
@ -38,12 +40,33 @@ impl<T:Copy + Mul<T, T> + Add<T, T> + Algebraic> Dot<T> for Vec2<T>
{
fn dot(&self, other : &Vec2<T>) -> T
{ self.x * other.x + self.y * other.y }
}
impl<T:Copy + Mul<T, T> + Add<T, T> + Quot<T, T> + Algebraic>
Norm<T> for Vec2<T>
{
fn sqnorm(&self) -> T
{ self.dot(self) }
fn norm(&self) -> T
{ self.sqnorm().sqrt() }
fn normalized(&self) -> Vec2<T>
{
let l = self.norm();
Vec2(self.x / l, self.y / l)
}
fn normalize(&mut self) -> T
{
let l = self.norm();
self.x /= l;
self.y /= l;
l
}
}
impl<T:Copy + Mul<T, T> + Sub<T, T>> Cross<Vec1<T>> for Vec2<T>
@ -70,6 +93,19 @@ impl<T:Copy + Zero> Zero for Vec2<T>
{ self.x.is_zero() && self.y.is_zero() }
}
impl<T: Copy + One + Zero + Neg<T>> Basis for Vec2<T>
{
fn canonical_basis() -> ~[Vec2<T>]
{
// FIXME: this should be static
~[ Vec2(One::one(), Zero::zero()),
Vec2(Zero::zero(), One::one()) ]
}
fn orthogonal_subspace_basis(&self) -> ~[Vec2<T>]
{ ~[ Vec2(-self.y, self.x) ] }
}
impl<T:FuzzyEq<T>> FuzzyEq<T> for Vec2<T>
{
fn fuzzy_eq(&self, other: &Vec2<T>) -> bool

View File

@ -128,7 +128,7 @@ Inv for Mat3<T>
- self.m12 * minor_m21_m33
+ self.m13 * minor_m21_m32;
assert!(det.is_zero());
assert!(!det.is_zero());
*self = Mat3(
(minor_m22_m33 / det),

View File

@ -1,9 +1,11 @@
use core::num::{Zero, Algebraic};
use core::num::{Zero, One, Algebraic, abs};
use core::rand::{Rand, Rng, RngUtil};
use std::cmp::FuzzyEq;
use traits::dim::Dim;
use traits::dot::Dot;
use traits::cross::Cross;
use traits::basis::Basis;
use traits::norm::Norm;
#[deriving(Eq)]
pub struct Vec3<T>
@ -44,12 +46,34 @@ impl<T:Copy + Mul<T, T> + Add<T, T> + Algebraic> Dot<T> for Vec3<T>
{
fn dot(&self, other : &Vec3<T>) -> T
{ self.x * other.x + self.y * other.y + self.z * other.z }
}
impl<T:Copy + Mul<T, T> + Add<T, T> + Quot<T, T> + Algebraic>
Norm<T> for Vec3<T>
{
fn sqnorm(&self) -> T
{ self.dot(self) }
fn norm(&self) -> T
{ self.sqnorm().sqrt() }
fn normalized(&self) -> Vec3<T>
{
let l = self.norm();
Vec3(self.x / l, self.y / l, self.z / l)
}
fn normalize(&mut self) -> T
{
let l = self.norm();
self.x /= l;
self.y /= l;
self.z /= l;
l
}
}
impl<T:Copy + Mul<T, T> + Sub<T, T>> Cross<Vec3<T>> for Vec3<T>
@ -76,6 +100,30 @@ impl<T:Copy + Zero> Zero for Vec3<T>
{ self.x.is_zero() && self.y.is_zero() && self.z.is_zero() }
}
impl<T: Copy + One + Zero + Neg<T> + Ord + Mul<T, T> + Sub<T, T> + Add<T, T> +
Quot<T, T> + Algebraic>
Basis for Vec3<T>
{
fn canonical_basis() -> ~[Vec3<T>]
{
// FIXME: this should be static
~[ Vec3(One::one(), Zero::zero(), Zero::zero()),
Vec3(Zero::zero(), One::one(), Zero::zero()),
Vec3(Zero::zero(), Zero::zero(), One::one()) ]
}
fn orthogonal_subspace_basis(&self) -> ~[Vec3<T>]
{
let a =
if (abs(self.x) > abs(self.y))
{ Vec3(self.z, Zero::zero(), -self.x).normalized() }
else
{ Vec3(Zero::zero(), -self.z, self.y).normalized() };
~[ a, a.cross(self) ]
}
}
impl<T:FuzzyEq<T>> FuzzyEq<T> for Vec3<T>
{
fn fuzzy_eq(&self, other: &Vec3<T>) -> bool

View File

@ -27,9 +27,12 @@ pub use traits::cross::Cross;
pub use traits::dim::Dim;
pub use traits::inv::Inv;
pub use traits::transpose::Transpose;
pub use traits::basis::Basis;
pub use traits::norm::Norm;
pub use traits::workarounds::rlmul::{RMul, LMul};
pub use traits::workarounds::trigonometric::Trigonometric;
pub use traits::workarounds::scalar_op::ScalarOp;
mod dim2
{
@ -68,10 +71,19 @@ mod traits
mod inv;
mod transpose;
mod dim;
mod basis;
mod norm;
mod workarounds
{
mod rlmul;
mod trigonometric;
mod scalar_op;
}
}
mod tests
{
mod mat;
mod vec;
}

View File

@ -1,9 +1,12 @@
use core::vec::{map_zip, from_elem, map, all, all2};
use core::num::{Zero, One, Algebraic};
use core::rand::{Rand, Rng, RngUtil};
use core::num::{Zero, Algebraic};
use core::vec::{map_zip, from_elem, map, all, all2};
use std::cmp::FuzzyEq;
use traits::dim::Dim;
use traits::dot::Dot;
use traits::norm::Norm;
use traits::basis::Basis;
use traits::workarounds::scalar_op::ScalarOp;
// D is a phantom parameter, used only as a dimensional token.
// Its allows use to encode the vector dimension at the type-level.
@ -24,6 +27,12 @@ impl<D: Dim, T> Dim for NVec<D, T>
{ Dim::dim::<D>() }
}
impl<D, T: Clone> Clone for NVec<D, T>
{
fn clone(&self) -> NVec<D, T>
{ NVec{ at: self.at.clone() } }
}
impl<D, T: Copy + Add<T,T>> Add<NVec<D, T>, NVec<D, T>> for NVec<D, T>
{
fn add(&self, other: &NVec<D, T>) -> NVec<D, T>
@ -54,12 +63,131 @@ Dot<T> for NVec<D, T>
res
}
}
impl<D: Dim, T: Copy + Mul<T, T> + Quot<T, T> + Add<T, T> + Sub<T, T>>
ScalarOp<T> for NVec<D, T>
{
fn scalar_mul(&self, s: &T) -> NVec<D, T>
{ NVec { at: map(self.at, |a| a * *s) } }
fn scalar_div(&self, s: &T) -> NVec<D, T>
{ NVec { at: map(self.at, |a| a / *s) } }
fn scalar_add(&self, s: &T) -> NVec<D, T>
{ NVec { at: map(self.at, |a| a + *s) } }
fn scalar_sub(&self, s: &T) -> NVec<D, T>
{ NVec { at: map(self.at, |a| a - *s) } }
fn scalar_mul_inplace(&mut self, s: &T)
{
for uint::range(0u, Dim::dim::<D>()) |i|
{ self.at[i] *= *s; }
}
fn scalar_div_inplace(&mut self, s: &T)
{
for uint::range(0u, Dim::dim::<D>()) |i|
{ self.at[i] /= *s; }
}
fn scalar_add_inplace(&mut self, s: &T)
{
for uint::range(0u, Dim::dim::<D>()) |i|
{ self.at[i] += *s; }
}
fn scalar_sub_inplace(&mut self, s: &T)
{
for uint::range(0u, Dim::dim::<D>()) |i|
{ self.at[i] -= *s; }
}
}
impl<D: Dim, T: Copy + Mul<T, T> + Add<T, T> + Quot<T, T> + Algebraic + Zero +
Clone>
Norm<T> for NVec<D, T>
{
fn sqnorm(&self) -> T
{ self.dot(self) }
fn norm(&self) -> T
{ self.sqnorm().sqrt() }
fn normalized(&self) -> NVec<D, T>
{
let mut res : NVec<D, T> = self.clone();
res.normalize();
res
}
fn normalize(&mut self) -> T
{
let l = self.norm();
for uint::range(0u, Dim::dim::<D>()) |i|
{ self.at[i] /= l; }
l
}
}
impl<D: Dim,
T: Copy + One + Zero + Neg<T> + Ord + Mul<T, T> + Sub<T, T> + Add<T, T> +
Quot<T, T> + Algebraic + Clone + FuzzyEq<T>>
Basis for NVec<D, T>
{
fn canonical_basis() -> ~[NVec<D, T>]
{
let dim = Dim::dim::<D>();
let mut res : ~[NVec<D, T>] = ~[];
for uint::range(0u, dim) |i|
{
let mut basis_element : NVec<D, T> = Zero::zero();
basis_element.at[i] = One::one();
res.push(basis_element);
}
res
}
fn orthogonal_subspace_basis(&self) -> ~[NVec<D, T>]
{
// compute the basis of the orthogonal subspace using Gram-Schmidt
// orthogonalization algorithm
let dim = Dim::dim::<D>();
let mut res : ~[NVec<D, T>] = ~[];
for uint::range(0u, dim) |i|
{
let mut basis_element : NVec<D, T> = Zero::zero();
basis_element.at[i] = One::one();
if (res.len() == dim - 1)
{ break; }
let mut elt = basis_element.clone();
elt -= self.scalar_mul(&basis_element.dot(self));
for res.each |v|
{ elt -= v.scalar_mul(&elt.dot(v)) };
if (!elt.sqnorm().fuzzy_eq(&Zero::zero()))
{ res.push(elt.normalized()); }
}
assert!(res.len() == dim - 1);
res
}
}
// FIXME: I dont really know how te generalize the cross product int

63
src/tests/mat.rs Normal file
View File

@ -0,0 +1,63 @@
#[test]
use core::num::{One};
#[test]
use core::rand::{random};
#[test]
use std::cmp::FuzzyEq;
// #[test]
// use ndim::nmat::NMat;
#[test]
use dim1::mat1::Mat1;
#[test]
use dim2::mat2::Mat2;
#[test]
use dim3::mat3::Mat3;
// #[test]
// use traits::dim::d7;
// FIXME: this one fails with an ICE: node_id_to_type: no type for node [...]
// #[test]
// fn test_inv_nmat()
// {
// let randmat : NMat<d7, f64> = random();
//
// assert!((randmat.inverse() * randmat).fuzzy_eq(&One::one()));
// }
#[test]
fn test_inv_mat1()
{
for uint::range(0u, 10000u) |_|
{
let randmat : Mat1<f64> = random();
assert!((randmat.inverse() * randmat).fuzzy_eq(&One::one()));
}
}
#[test]
fn test_inv_mat2()
{
for uint::range(0u, 10000u) |_|
{
let randmat : Mat2<f64> = random();
assert!((randmat.inverse() * randmat).fuzzy_eq(&One::one()));
}
}
#[test]
fn test_inv_mat3()
{
for uint::range(0u, 10000u) |_|
{
let randmat : Mat3<f64> = random();
assert!((randmat.inverse() * randmat).fuzzy_eq(&One::one()));
}
}
#[test]
fn test_rot2()
{
}

212
src/tests/vec.rs Normal file
View File

@ -0,0 +1,212 @@
#[test]
use core::num::{Zero, One};
#[test]
use core::rand::{random};
#[test]
use core::vec::{all, all2};
#[test]
use std::cmp::FuzzyEq;
#[test]
use dim3::vec3::Vec3;
#[test]
use dim2::vec2::Vec2;
#[test]
use dim1::vec1::Vec1;
#[test]
use ndim::nvec::NVec;
#[test]
use traits::dim::d7;
#[test]
use traits::basis::Basis;
#[test]
fn test_cross_vec3()
{
for uint::range(0u, 10000u) |_|
{
let v1 : Vec3<f64> = random();
let v2 : Vec3<f64> = random();
let v3 : Vec3<f64> = v1.cross(&v2);
assert!(v3.dot(&v2).fuzzy_eq(&Zero::zero()));
assert!(v3.dot(&v1).fuzzy_eq(&Zero::zero()));
}
}
#[test]
fn test_dot_nvec()
{
for uint::range(0u, 10000u) |_|
{
let v1 : NVec<d7, f64> = random();
let v2 : NVec<d7, f64> = random();
assert!(v1.dot(&v2).fuzzy_eq(&v2.dot(&v1)));
}
}
#[test]
fn test_commut_dot_vec3()
{
for uint::range(0u, 10000u) |_|
{
let v1 : Vec3<f64> = random();
let v2 : Vec3<f64> = random();
assert!(v1.dot(&v2).fuzzy_eq(&v2.dot(&v1)));
}
}
#[test]
fn test_commut_dot_vec2()
{
for uint::range(0u, 10000u) |_|
{
let v1 : Vec2<f64> = random();
let v2 : Vec2<f64> = random();
assert!(v1.dot(&v2).fuzzy_eq(&v2.dot(&v1)));
}
}
#[test]
fn test_commut_dot_vec1()
{
for uint::range(0u, 10000u) |_|
{
let v1 : Vec1<f64> = random();
let v2 : Vec1<f64> = random();
assert!(v1.dot(&v2).fuzzy_eq(&v2.dot(&v1)));
}
}
#[test]
fn test_basis_vec1()
{
for uint::range(0u, 10000u) |_|
{
let basis = Basis::canonical_basis::<Vec1<f64>>();
// check vectors form an ortogonal basis
assert!(all2(basis, basis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(basis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_basis_vec2()
{
for uint::range(0u, 10000u) |_|
{
let basis = Basis::canonical_basis::<Vec2<f64>>();
// check vectors form an ortogonal basis
assert!(all2(basis, basis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(basis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_basis_vec3()
{
for uint::range(0u, 10000u) |_|
{
let basis = Basis::canonical_basis::<Vec3<f64>>();
// check vectors form an ortogonal basis
assert!(all2(basis, basis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(basis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_basis_nvec()
{
for uint::range(0u, 10000u) |_|
{
let basis = Basis::canonical_basis::<NVec<d7, f64>>();
// check vectors form an ortogonal basis
assert!(all2(basis, basis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(basis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_subspace_basis_vec1()
{
for uint::range(0u, 10000u) |_|
{
let v : Vec1<f64> = random();
let v1 = v.normalized();
let subbasis = v1.orthogonal_subspace_basis();
// check vectors are orthogonal to v1
assert!(all(subbasis, |e| v1.dot(e).fuzzy_eq(&Zero::zero())));
// check vectors form an ortogonal basis
assert!(all2(subbasis, subbasis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(subbasis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_subspace_basis_vec2()
{
for uint::range(0u, 10000u) |_|
{
let v : Vec2<f64> = random();
let v1 = v.normalized();
let subbasis = v1.orthogonal_subspace_basis();
// check vectors are orthogonal to v1
assert!(all(subbasis, |e| v1.dot(e).fuzzy_eq(&Zero::zero())));
// check vectors form an ortogonal basis
assert!(all2(subbasis, subbasis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(subbasis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
#[test]
fn test_subspace_basis_vec3()
{
for uint::range(0u, 10000u) |_|
{
let v : Vec3<f64> = random();
let v1 = v.normalized();
let subbasis = v1.orthogonal_subspace_basis();
// check vectors are orthogonal to v1
assert!(all(subbasis, |e| v1.dot(e).fuzzy_eq(&Zero::zero())));
// check vectors form an ortogonal basis
assert!(all2(subbasis, subbasis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// check vectors form an orthonormal basis
assert!(all(subbasis, |e| e.norm().fuzzy_eq(&One::one())));
}
}
// ICE
//
// #[test]
// fn test_subspace_basis_vecn()
// {
// for uint::range(0u, 10000u) |_|
// {
// let v : NVec<d7, f64> = random();
// let v1 = v.normalized();
// let subbasis = v1.orthogonal_subspace_basis();
//
// // check vectors are orthogonal to v1
// assert!(all(subbasis, |e| v1.dot(e).fuzzy_eq(&Zero::zero())));
// // check vectors form an ortogonal basis
// assert!(all2(subbasis, subbasis, |e1, e2| e1 == e2 || e1.dot(e2).fuzzy_eq(&Zero::zero())));
// // check vectors form an orthonormal basis
// assert!(all(subbasis, |e| e.norm().fuzzy_eq(&One::one())));
// }
// }

5
src/traits/basis.rs Normal file
View File

@ -0,0 +1,5 @@
pub trait Basis
{
fn canonical_basis() -> ~[Self]; // FIXME: is it the right pointer?
fn orthogonal_subspace_basis(&self) -> ~[Self];
}

View File

@ -1,5 +1,4 @@
pub trait Dim
{
pub trait Dim {
/// The dimension of the object.
fn dim() -> uint;
}
@ -8,20 +7,44 @@ pub trait Dim
// object at the type-level.
/// Dimensional token for 0-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d0;
/// Dimensional token for 1-dimension. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d1;
/// Dimensional token for 2-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d2;
/// Dimensional token for 3-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d3;
/// Dimensional token for 4-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d4;
/// Dimensional token for 5-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d5;
/// Dimensional token for 6-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d6;
/// Dimensional token for 7-dimensions. Dimensional tokens are the preferred
/// way to specify at the type level the dimension of n-dimensional objects.
#[deriving(Eq)]
pub struct d7;
impl Dim for d0
{ fn dim() -> uint { 0 } }
@ -36,3 +59,12 @@ impl Dim for d3
impl Dim for d4
{ fn dim() -> uint { 4 } }
impl Dim for d5
{ fn dim() -> uint { 5 } }
impl Dim for d6
{ fn dim() -> uint { 6 } }
impl Dim for d7
{ fn dim() -> uint { 7 } }

View File

@ -2,13 +2,4 @@ pub trait Dot<T>
{
/// Computes the dot (inner) product of two objects.
fn dot(&self, &Self) -> T;
/// Computes the norm a an object.
fn norm(&self) -> T;
/**
* Computes the squared norm of an object.
*
* Computes the squared norm of an object. Computation of the squared norm
* is usually faster than the norm itself.
*/
fn sqnorm(&self) -> T; // { self.dot(self); }
}

19
src/traits/norm.rs Normal file
View File

@ -0,0 +1,19 @@
pub trait Norm<T>
{
/// Computes the norm a an object.
fn norm(&self) -> T;
/**
* Computes the squared norm of an object.
*
* Computes the squared norm of an object. Computation of the squared norm
* is usually faster than the norm itself.
*/
fn sqnorm(&self) -> T;
/// Returns the normalized version of the argument.
fn normalized(&self) -> Self;
/// Inplace version of `normalized`.
fn normalize(&mut self) -> T;
}

View File

@ -0,0 +1,12 @@
pub trait ScalarOp<T>
{
fn scalar_mul(&self, &T) -> Self;
fn scalar_div(&self, &T) -> Self;
fn scalar_add(&self, &T) -> Self;
fn scalar_sub(&self, &T) -> Self;
fn scalar_mul_inplace(&mut self, &T);
fn scalar_div_inplace(&mut self, &T);
fn scalar_add_inplace(&mut self, &T);
fn scalar_sub_inplace(&mut self, &T);
}