Basis trait now uses internal iterators to avoid allocations.

This commit is contained in:
Sébastien Crozet 2013-07-01 16:33:22 +00:00
parent 6fd9696253
commit 68d601a642
4 changed files with 52 additions and 56 deletions

View File

@ -79,15 +79,13 @@ macro_rules! test_basis_impl(
($t: ty) => ( ($t: ty) => (
for 10000.times for 10000.times
{ {
let basis = Basis::canonical_basis::<$t>(); do Basis::canonical_basis::<$t> |e1|
{
do Basis::canonical_basis::<$t> |e2|
{ assert!(e1 == e2 || e1.dot(&e2).approx_eq(&Zero::zero())) }
// check vectors form an ortogonal basis assert!(e1.norm().approx_eq(&One::one()));
assert!( }
do basis.iter().zip(basis.iter()).all
|(e1, e2)| { e1 == e2 || e1.dot(e2).approx_eq(&Zero::zero()) }
);
// check vectors form an orthonormal basis
assert!(basis.iter().all(|e| e.norm().approx_eq(&One::one())));
} }
); );
) )
@ -98,17 +96,17 @@ macro_rules! test_subspace_basis_impl(
{ {
let v : $t = random(); let v : $t = random();
let v1 = v.normalized(); let v1 = v.normalized();
let subbasis = v1.orthogonal_subspace_basis();
// check vectors are orthogonal to v1 do v1.orthonormal_subspace_basis() |e1|
assert!(subbasis.iter().all(|e| v1.dot(e).approx_eq(&Zero::zero()))); {
// check vectors form an ortogonal basis // check vectors are orthogonal to v1
assert!( assert!(v1.dot(&e1).approx_eq(&Zero::zero()));
do subbasis.iter().zip(subbasis.iter()).all // check vectors form an orthonormal basis
|(e1, e2)| { e1 == e2 || e1.dot(e2).approx_eq(&Zero::zero()) } assert!(e1.norm().approx_eq(&One::one()));
); // check vectors form an ortogonal basis
// check vectors form an orthonormal basis do v1.orthonormal_subspace_basis() |e2|
assert!(subbasis.iter().all(|e| e.norm().approx_eq(&One::one()))); { assert!(e1 == e2 || e1.dot(&e2).approx_eq(&Zero::zero())) }
}
} }
); );
) )

View File

@ -1,8 +1,7 @@
pub trait Basis pub trait Basis
{ {
/// Computes the canonical basis of the space in which this object lives. /// Computes the canonical basis of the space in which this object lives.
// FIXME: need type-associated values // FIXME: implement the for loop protocol?
// FIXME: this will make allocations… this is bad fn canonical_basis(&fn(Self));
fn canonical_basis() -> ~[Self]; fn orthonormal_subspace_basis(&self, &fn(Self));
fn orthogonal_subspace_basis(&self) -> ~[Self];
} }

View File

@ -88,27 +88,23 @@ macro_rules! basis_impl(
($t: ident, $dim: expr) => ( ($t: ident, $dim: expr) => (
impl<N: Copy + DivisionRing + Algebraic + ApproxEq<N>> Basis for $t<N> impl<N: Copy + DivisionRing + Algebraic + ApproxEq<N>> Basis for $t<N>
{ {
pub fn canonical_basis() -> ~[$t<N>] pub fn canonical_basis(f: &fn($t<N>))
{ {
let mut res : ~[$t<N>] = ~[];
for iterate(0u, $dim) |i| for iterate(0u, $dim) |i|
{ {
let mut basis_element : $t<N> = Zero::zero(); let mut basis_element : $t<N> = Zero::zero();
basis_element.at[i] = One::one(); basis_element.at[i] = One::one();
res.push(basis_element); f(basis_element);
} }
res
} }
pub fn orthogonal_subspace_basis(&self) -> ~[$t<N>] pub fn orthonormal_subspace_basis(&self, f: &fn($t<N>))
{ {
// compute the basis of the orthogonal subspace using Gram-Schmidt // compute the basis of the orthogonal subspace using Gram-Schmidt
// orthogonalization algorithm // orthogonalization algorithm
let mut res : ~[$t<N>] = ~[]; let mut basis: ~[$t<N>] = ~[];
for iterate(0u, $dim) |i| for iterate(0u, $dim) |i|
{ {
@ -116,21 +112,25 @@ macro_rules! basis_impl(
basis_element.at[i] = One::one(); basis_element.at[i] = One::one();
if res.len() == $dim - 1 if basis.len() == $dim - 1
{ break; } { break; }
let mut elt = copy basis_element; let mut elt = copy basis_element;
elt = elt - self.scalar_mul(&basis_element.dot(self)); elt = elt - self.scalar_mul(&basis_element.dot(self));
for res.iter().advance |v| for basis.iter().advance |v|
{ elt = elt - v.scalar_mul(&elt.dot(v)) }; { elt = elt - v.scalar_mul(&elt.dot(v)) };
if !elt.sqnorm().approx_eq(&Zero::zero()) if !elt.sqnorm().approx_eq(&Zero::zero())
{ res.push(elt.normalized()); } {
} let new_element = elt.normalized();
res f(copy new_element);
basis.push(new_element);
}
}
} }
} }
) )

View File

@ -27,44 +27,42 @@ impl<N: Mul<N, N> + Sub<N, N>> Cross<Vec3<N>> for Vec3<N>
impl<N: One> Basis for Vec1<N> impl<N: One> Basis for Vec1<N>
{ {
#[inline] #[inline(always)]
fn canonical_basis() -> ~[Vec1<N>] fn canonical_basis(f: &fn(Vec1<N>))
{ ~[ Vec1::new([One::one()]) ] } // FIXME: this should be static { f(Vec1::new([One::one()])) }
#[inline] #[inline(always)]
fn orthogonal_subspace_basis(&self) -> ~[Vec1<N>] fn orthonormal_subspace_basis(&self, _: &fn(Vec1<N>))
{ ~[] } { }
} }
impl<N: Copy + One + Zero + Neg<N>> Basis for Vec2<N> impl<N: Copy + One + Zero + Neg<N>> Basis for Vec2<N>
{ {
#[inline] #[inline]
fn canonical_basis() -> ~[Vec2<N>] fn canonical_basis(f: &fn(Vec2<N>))
{ {
// FIXME: this should be static f(Vec2::new([One::one(), Zero::zero()]));
~[ Vec2::new([One::one(), Zero::zero()]), f(Vec2::new([Zero::zero(), One::one()]));
Vec2::new([Zero::zero(), One::one()]) ]
} }
#[inline] #[inline]
fn orthogonal_subspace_basis(&self) -> ~[Vec2<N>] fn orthonormal_subspace_basis(&self, f: &fn(Vec2<N>))
{ ~[ Vec2::new([-self.at[1], copy self.at[0]]) ] } { f(Vec2::new([-self.at[1], copy self.at[0]])) }
} }
impl<N: Copy + DivisionRing + Ord + Algebraic> impl<N: Copy + DivisionRing + Ord + Algebraic>
Basis for Vec3<N> Basis for Vec3<N>
{ {
#[inline] #[inline(always)]
fn canonical_basis() -> ~[Vec3<N>] fn canonical_basis(f: &fn(Vec3<N>))
{ {
// FIXME: this should be static f(Vec3::new([One::one(), Zero::zero(), Zero::zero()]));
~[ Vec3::new([One::one(), Zero::zero(), Zero::zero()]), f(Vec3::new([Zero::zero(), One::one(), Zero::zero()]));
Vec3::new([Zero::zero(), One::one(), Zero::zero()]), f(Vec3::new([Zero::zero(), Zero::zero(), One::one()]));
Vec3::new([Zero::zero(), Zero::zero(), One::one()]) ]
} }
#[inline] #[inline(always)]
fn orthogonal_subspace_basis(&self) -> ~[Vec3<N>] fn orthonormal_subspace_basis(&self, f: &fn(Vec3<N>))
{ {
let a = let a =
if abs(copy self.at[0]) > abs(copy self.at[1]) if abs(copy self.at[0]) > abs(copy self.at[1])
@ -72,6 +70,7 @@ Basis for Vec3<N>
else else
{ Vec3::new([Zero::zero(), -self.at[2], copy self.at[1]]).normalized() }; { Vec3::new([Zero::zero(), -self.at[2], copy self.at[1]]).normalized() };
~[ a.cross(self), a ] f(a.cross(self));
f(a);
} }
} }