1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
//!
//! Traits for picking the types to store the internal data of the algebraic structs
//!
//! Ideally, this would be done using entirely `const` generics, but this is impossible for a
//! couple reasons:
//! 1. Unlike for simple vectors and matrices, the number of elements in blades, evens, etc
//!    is not simply the dimension or grade but instead a more involved function of them. Hence,
//!    since generic `const` expressions aren't currently supported, we
//!    cannot use them to construct the backing array.
//! 2. We need to *also* support non-array storage for [`Dynamic`] dimensions and grades
//!
//! Instead, we use a trait based system with associated types to select what the backing buffer
//! should be.
//!
//! The core of this system are the [`AllocBlade<N,G>`], [`AllocEven<N>`], [`AllocOdd<N>`], and
//! [`AllocMultivector<N>`] traits. These are all given blanket `impl`s on all types for all
//! supported of dimensions and grades `N` and `G`. Then, in order to retrieve the storage the
//! for a particular algebraic struct, we can use the corresponding `Alloc*` trait and its
//! associated type `Buffer`, adding the `Alloc*` bound to the current scalar being used.
//!
//! At the moment, for static array allocation, only dimensions and grades up to 16 are supported.
//! This is both due to the fact we don't have generic const expressions *and* out of practical
//! considerations regarding sizes in high numbers of dimensions.
//! (`f32` Multivectors in 16D already use ~**260KB**!!) However,
//! if a high number of dimensions are required, using a [`Dynamic`] dimension or grade will allow
//! that number of dimensions to be used, just using the heap instead
//!
//!

use std::mem::MaybeUninit;
use std::iter::Iterator;

use crate::base::*;

use std::ops::{ Add, BitOr };
use typenum::{
    IsLessOrEqual, Sum, LeEq, True, U1
};

#[cfg(doc)] use crate::algebra::*;
#[cfg(doc)] use crate::subspace::*;

/// Implements [`Alloc`] and makes the default choice for what storage types to use
pub struct DefaultAllocator;

/// The storage type to use for `M`
pub type Allocate<M> = <DefaultAllocator as Alloc<M>>::Buffer;

/// The backing buffer type for a [`Blade`] with scalar `T`, dimension `N`, and grade `G`
pub type AllocateBlade<T,N,G> = <T as AllocBlade<N,G>>::Buffer;

/// The backing buffer type for an [`Even`] with scalar `T` and dimension `N`
pub type AllocateEven<T,N> = <T as AllocEven<N>>::Buffer;

/// The backing buffer type for an [`Odd`] with scalar `T` and dimension `N`
pub type AllocateOdd<T,N> = <T as AllocOdd<N>>::Buffer;

/// The backing buffer type for a [`Multivector`] with scalar `T` and dimension `N`
pub type AllocateMultivector<T,N> = <T as AllocMultivector<N>>::Buffer;

///
/// Picks the buffer type to use for the generic type `M` and provides extra helper methods
///
/// The intent is that this is implemented onto a ZST for every supported algebraic struct
/// in order to represent a particular allocation system. At the moment, only one exists
/// ([`DefaultAllocator`]), but more may be added in the future.
///
/// Additionally, this allows us to make the [`Allocate`] type alias so that every algebraic
/// struct uses `Allocate<Self>` for it's internal buffer.
///
/// Finally, this trait includes some helper methods for actually creating and managing the backing
/// buffer.
///
pub unsafe trait Alloc<M> {

    ///The type to store in the backing bufferfor a simplification of some of the systems
    type Scalar: Sized;

    ///A type representing the dimension and/or grade of this structure
    type Shape: Copy;

    ///The type to store the data in
    type Buffer: Storage<Self::Scalar, Uninit=Self::Uninit>;

    ///The type to store uninitialized data in
    type Uninit: UninitStorage<Self::Scalar, Init=Self::Buffer>;

    ///Returns the dimension and potentially grade depending on what type `M` is
    fn shape(this: &M) -> Self::Shape;
    ///Makes an uninitialized buffer
    fn uninit(shape: Self::Shape) -> Self::Uninit;
    ///Creates an `M` from an uninitialized buffer assuming it is initialized
    unsafe fn assume_init(uninit: Self::Uninit) -> M;

}

///Picks the buffer to use for a [`Blade`] of dimension `N`, grade `G`, and using `Self` as the scalar
pub trait AllocBlade<N:Dim,G:Dim>: Sized {
    type Buffer: BladeStorage<Self,N,G>;
}

///Picks the buffer to use for an [`Even`] of dimension `N` and using `Self` as the scalar
pub trait AllocEven<N:Dim>: Sized {
    type Buffer: EvenStorage<Self,N>;
}

///Picks the buffer to use for an [`Odd`] of dimension `N` and using `Self` as the scalar
pub trait AllocOdd<N:Dim>: Sized {
    type Buffer: OddStorage<Self,N>;
}

///A marker trait for when a dimension `N` is supported for both [`Even`]s and [`Odd`]s
pub trait AllocVersor<N:Dim>: AllocEven<N> + AllocOdd<N> {}
impl<T:AllocEven<N>+AllocOdd<N>, N:Dim> AllocVersor<N> for T {}

///Picks the buffer to use for an [`Multivector`] of dimension `N`, grade `G`, and using `Self` as the scalar
pub trait AllocMultivector<N:Dim>: Sized {
    type Buffer: MultivectorStorage<Self,N>;
}

///
/// Marks if a [Blade] of a given dimension and grade is guaranteed to always be [simple](SimpleBlade)
///
/// This trait is implemented for scalars, vectors, psuedovectors, and psuedoscalars by bounding the
/// grade to be 0, 1, N, or N-1.
///
/// It is used in a number of locations to allow for certain operations that are only possible with
/// simple blades, ie mutating `SimpleBlade`s, [projecting](Blade::project) onto a Blade, etc.
///
/// At the moment, the implementation uses `typenum` expressions to determine type eligibility,
/// which may be lead to some provability issues when using generic types, especially
/// for psuedoscalars and psuedovectors. However, this will be simplified dramatically once the
/// `#[feature(marker_trait_attr)]` is stabilized,
///
pub trait AllocSimpleBlade<N:Dim,G:Dim>: AllocBlade<N,G> {}

impl<T:AllocBlade<Dynamic,Const<0>>> AllocSimpleBlade<Dynamic,Const<0>> for T {}
impl<T:AllocBlade<Dynamic,Const<1>>> AllocSimpleBlade<Dynamic,Const<1>> for T {}
impl<T:AllocBlade<N,G>,N:Dim,G:Dim> AllocSimpleBlade<N,G> for T where
    //establish that `N` and `G` are constant numbers
    N: DimName+ToTypenum,
    G: DimName+ToTypenum,

    //establish that we can compare `G` with `1` and `N` with `G+1`
    G::Typenum: Add<U1> + IsLessOrEqual<U1>,
    N::Typenum: IsLessOrEqual<Sum<G::Typenum,U1>>,

    //`or` together the two comparisons
    LeEq<G::Typenum, U1>: BitOr<LeEq<N::Typenum, Sum<G::Typenum,U1>>, Output=True>
{}

impl<T, const N: usize> AllocBlade<Const<N>, Dynamic> for T {
    type Buffer = DynBladeStorage<T, Const<N>, Dynamic>;
}

impl<T, const G: usize> AllocBlade<Dynamic, Const<G>> for T {
    type Buffer = DynBladeStorage<T, Dynamic, Const<G>>;
}

impl<T> AllocBlade<Dynamic, Dynamic> for T {
    type Buffer = DynBladeStorage<T, Dynamic, Dynamic>;
}

impl<T> AllocEven<Dynamic> for T {
    type Buffer = DynEvenStorage<T, Dynamic>;
}

impl<T> AllocOdd<Dynamic> for T {
    type Buffer = DynOddStorage<T, Dynamic>;
}

impl<T> AllocMultivector<Dynamic> for T {
    type Buffer = DynMultivectorStorage<T, Dynamic>;
}

#[inline(always)]
fn uninit_array<T, const L: usize>() -> [MaybeUninit<T>; L] {
    //TODO: use `MaybeUninit::uninit_array()` when stabilized
    unsafe {
        //the outer MaybeUninit wraps the [MaybeUninit<T>; L] array
        MaybeUninit::uninit().assume_init()
    }
}

#[inline(always)]
fn array_from_iter<T, I: IntoIterator<Item=T>, const L: usize>(iter:I, kind:&str) -> [T;L] {
    //TODO: use `MaybeUninit::uninit_array()` when stabilized
    let mut uninit: [MaybeUninit<T>;L] = uninit_array();

    //fill the array and count how many spaces we actually fill
    let mut count = 0;
    for (i, x) in (0..L).zip(iter) {
        uninit[i] = MaybeUninit::new(x);
        count = i+1;
    }

    //panic if not enough elements
    if count!=L {
        panic!("Not enough elements to fill {}", kind);
    }

    unsafe { <[MaybeUninit<T>;L] as UninitStorage<T>>::assume_init(uninit) }
}

macro_rules! impl_alloc{

    ($n:literal $($N:literal)*; $($G:literal)*; @$cmd:ident $($pairs:tt)*) => {
        impl_alloc!($($N)*; $($G)*; @$cmd $($pairs)* $(($n, $G))*);

        impl_alloc!($n @$cmd);

    };

    //reuse the macro loop to produce tests for the given dimension
    ($N:literal @tests) => {
        assert_eq!(
            std::mem::size_of::<AllocateEven<f32, Const<$N>>>(),
            //this has some weird behavior
            std::mem::size_of::<f32>() * even_elements($N)
        );

        assert_eq!(
            std::mem::size_of::<AllocateOdd<f32, Const<$N>>>(),
            //this has some weird behavior
            std::mem::size_of::<f32>() * odd_elements($N)
        );

        assert_eq!(
            std::mem::size_of::<AllocateMultivector<f32, Const<$N>>>(),
            std::mem::size_of::<f32>() * 2usize.pow($N)
        );
    };

    ($N:literal @impl) => {
        impl<T> AllocEven<Const<$N>> for T {
            type Buffer = [T; even_elements($N)];
        }

        unsafe impl<T> EvenStorage<T, Const<$N>> for [T; even_elements($N) ] {
            fn dim(&self) -> Const<$N> { Const::<$N> }
            fn uninit(_: Const<$N>,) -> Self::Uninit { uninit_array() }
            fn from_iterator<I:IntoIterator<Item=T>>(_: Const<$N>, iter: I) -> Self {
                array_from_iter(iter, "value")
            }
        }

        impl<T> AllocOdd<Const<$N>> for T {
            type Buffer = [T; odd_elements($N)];
        }

        unsafe impl<T> OddStorage<T, Const<$N>> for [T; odd_elements($N) ] {
            fn dim(&self) -> Const<$N> { Const::<$N> }
            fn uninit(_: Const<$N>,) -> Self::Uninit { uninit_array() }
            fn from_iterator<I:IntoIterator<Item=T>>(_: Const<$N>, iter: I) -> Self {
                array_from_iter(iter, "value")
            }
        }

        impl<T> AllocMultivector<Const<$N>> for T {
            type Buffer = [T; 2usize.pow($N)];
        }

        unsafe impl<T> MultivectorStorage<T, Const<$N>> for [T; 2usize.pow($N)] {
            fn dim(&self) -> Const<$N> { Const::<$N> }
            fn uninit(_: Const<$N>,) -> Self::Uninit { uninit_array() }
            fn from_iterator<I:IntoIterator<Item=T>>(_: Const<$N>, iter: I) -> Self {
                array_from_iter(iter, "multivector")
            }
        }

    };

    (; $($_G:literal)*; @impl $(($N:literal, $G:literal))*) => {
        $(
            impl<T> AllocBlade<Const<$N>, Const<$G>> for T {
                type Buffer = [T; binom($N, $G)];
            }

            unsafe impl<T> BladeStorage<T, Const<$N>, Const<$G>> for [T; binom($N, $G)] {
                fn dim(&self) -> Const<$N> { Const::<$N> }
                fn grade(&self) -> Const<$G> { Const::<$G> }

                fn uninit(_: Const<$N>, _: Const<$G>) -> Self::Uninit { uninit_array() }

                fn from_iterator<I:IntoIterator<Item=T>>(_: Const<$N>, _: Const<$G>, iter: I) -> Self {
                    array_from_iter(iter, "blade")
                }
            }

        )*
    };

    //reuse the macro loop to produce tests for the given dimension and grade
    (; $($_G:literal)*; @tests $(($N:literal, $G:literal))*) => {
        $(
            assert_eq!(
                std::mem::size_of::<AllocateBlade<f32, Const<$N>, Const<$G>>>(),
                std::mem::size_of::<f32>() * binom($N, $G)
            );
        )*
    };
}

impl_alloc!(
    0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16; 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16; @impl
);

#[test]
#[cfg(test)]
fn buffer_sizes() {
    impl_alloc!(
        0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16; 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16; @tests
    );
}