Proposal: NTuple into the STL and/or runtime?

[HyperSphereStudio]

Note This proposal is inspired by the Julia Programming Language

One usage of the NTuple is for N dimensional abstraction. I have noticed in the C# STL the usage of Vector2, Vector3 etc. While this does work, a lot of abstraction power can come from writing code that can take an arbitrary ntuple size (using generics for performance).

Benefits:

• Smaller Source Code Size (Less Code Bloat from having to write the same thing over and over again for different dimensions)
• Greater algorithmic control (You can account for every dimensional case when you write a function, rather than specializing to a particular dimension)
• No runtime penalty (memory, speed, garbage collection etc.) for abstraction (Due to use of generics)
• Users dont have to do things like use stackalloc [N] or fixed buffers for similar outcomes

Downsides:

• Potential Code Generation Bloating

I believe that the addition of the NTuples into the .NET ecosystem can increase the abstractional power of C# and allow for far more advanced numerical algorithms that would traditionally be very annoying to implement and not necessarily COMPLETE. Actual implementation into the STL could simply a library similar to below, but if the runtime implemented NTuples, the approach could be simplified (add in SIMD?).

I was inspired to write the NTuple library to create a parallel to [GeometryBasics]. (https://github.com/JuliaGeometry/GeometryBasics.jl/tree/master) library

Performance Note: This library produces a lot of garbage because of the closures at the moment... work will needed in order to make the speed (which it should be!) equivalent to manual versions.

public struct Vec<T, NT> where NT: NTuple<T, NT> where T: IFloatingPointIeee754<T> {
    public NT Data;
    
    public Vec(T v) => Data = NT.Create(_ => v);
    public Vec(NT data) => Data = data;
    public Vec(Vec<T, NT> v) => Data = v.Data;
    
    public static Vec<T, NT> Unit(int n) => new(NT.Create(i => i == n ? T.One : T.Zero));

    public static Vec<T, NT> UnitX { get; } = Unit(0);
    public static Vec<T, NT> UnitY { get; } = Unit(1);
    public static Vec<T, NT> UnitZ { get; } = Unit(2);
    public static Vec<T, NT> Zero { get; } = new(T.Zero);
    
    [IgnoreDataMember] public T LengthSquared => Data.Reduce(T.One, (a, b) => a+b*b);
    [IgnoreDataMember] public T Length => T.Sqrt(LengthSquared);
    [IgnoreDataMember] public T ReciprocalSqrt => T.ReciprocalSqrtEstimate(LengthSquared);
    [IgnoreDataMember] public T X => Data[0];
    [IgnoreDataMember] public T Y => Data[1];
    [IgnoreDataMember] public T Z => Data[2];
    [IgnoreDataMember] public int Dimension => Data.Count;

    public T this[int index] {
        get => Data[index];
        set => Data[index] = value;
    }


    public static Vec<T, NT> operator +(Vec<T, NT> a, Vec<T, NT> b) => new(NT.Create(i => a[i] + b[i]));
    public static Vec<T, NT> operator -(Vec<T, NT> a, Vec<T, NT> b) => new(NT.Create(i => a[i] - b[i]));
    public static Vec<T, NT> operator +(Vec<T, NT> a, T b) => new(NT.Create(i => a[i] + b));
    public static Vec<T, NT> operator +(T a, Vec<T, NT> b) => new(NT.Create(i => a + b[i]));
    public static Vec<T, NT> operator -(T a, Vec<T, NT> b) => new(NT.Create(i => a - b[i]));
    public static Vec<T, NT> operator -(Vec<T, NT> a, T b) => new(NT.Create(i => a[i] - b));
    public static Vec<T, NT> operator *(T a, Vec<T, NT> b) => new(NT.Create(i => a * b[i]));
    public static Vec<T, NT> operator *(Vec<T, NT> a, T b) => new(NT.Create(i => a[i] * b));
    public static Vec<T, NT> operator /(Vec<T, NT> a, T b) => new(NT.Create(i => a[i] / b));
}

//https://github.com/JuliaGeometry/GeometryBasics.jl/blob/master/src/primitives/rectangles.jl
public struct HyperRect<T, NT>(Vec<T, NT> origin, Vec<T, NT> widths) where T : IFloatingPointIeee754<T> where NT : NTuple<T, NT> {
    public Vec<T, NT> Origin = origin, Widths = widths;

    [IgnoreDataMember] public T Volume => Widths.Data.Prod();
    [IgnoreDataMember] public Vec<T, NT> Minima => Origin;
    [IgnoreDataMember] public Vec<T, NT> Maxima => Origin + Widths;
    [IgnoreDataMember] public Vec<T, NT> Center => Origin + Widths/(T.One + T.One);
    
    /*
     * Check if Rect `this` is contained in `v`. This does not use
    strict inequality, so Rects may share faces and this will still
    return true.
     */
    public bool Inside(HyperRect<T, NT> container) {
        var maxThis = Maxima;
        var maxV = container.Maxima;
        for (var i = 0; i < Origin.Dimension; i++) {
            if (!(maxThis[i] <= maxV[i] && Minima[i] >= container.Minima[i]))
                return false;
        }
        return true;
    }
}

The NTuple Library

To perform efficient field lookup using the element reference function I made a rule that the static allocated tuples should have a sequential layout with no packing, maybe this should be improved?

Here is the base interface for the NTuple I made. It mostly just implements the list interface into all NTuple variants

public interface NTuple<T, NT> : IList<T> where NT: NTuple<T, NT>{
    public static abstract NT Create(Func<int, T> itr);
    public static abstract ref T GetElementReference(ref NT v, int index);
    
    public T Reduce(T identity, Func<T, T, T> reduction) {
        var v = identity;
        for (int i = 0; i < Count; i++)
            v = reduction(v, this[i]);
        return v;
    }
    
    IEnumerator<T> IEnumerable<T>.GetEnumerator() {
        var i = 0;
        while (i < Count)
            yield return this[i];
    }

    IEnumerator IEnumerable.GetEnumerator() => GetEnumerator();
    bool ICollection<T>.IsReadOnly => false;
    int IList<T>.IndexOf(T item) {
        for (var i = 0; i < Count; i++) {
            if (EqualityComparer<T>.Default.Equals(this[i], item))
                return i;
        }
        return -1;
    }
    
    bool ICollection<T>.Contains(T item) {
        foreach (var w in this) {
            if (EqualityComparer<T>.Default.Equals(w, item))
                return true;
        }
        return false;
    }
    void ICollection<T>.CopyTo(T[] array, int arrayIndex) => CopyTo(array.AsSpan(arrayIndex));
    public void CopyTo(Span<T> values) {
        for (int i = 0; i < Count; i++)
            values[i] = this[i];
    }

    void ICollection<T>.Add(T item) => throw new NotSupportedException();
    void ICollection<T>.Clear() => throw new NotSupportedException();
    bool ICollection<T>.Remove(T item) => throw new NotSupportedException();
    void IList<T>.Insert(int index, T item) => throw new NotSupportedException();
    void IList<T>.RemoveAt(int index) => throw new NotSupportedException();
}

Here are the first three ntuple's (I made up to seven)

[StructLayout(LayoutKind.Explicit, Size = 0)]
public struct NTuple0<T> : NTuple<T, NTuple0<T>> {
    [IgnoreDataMember] public int Count => 0;
    public NTuple0(){}
    
    public static NTuple0<T> Create(Func<int, T> itr) => default;
    public static unsafe ref T GetElementReference(ref NTuple0<T> v, int index) => ref Unsafe.AsRef<T>((byte*)Unsafe.AsPointer(ref v));

    public T Reduce(T identity, Func<T, T, T> reduction) => identity;

    public T this[int index] {
        get => throw new ArgumentOutOfRangeException(nameof(index));
        set => throw new ArgumentOutOfRangeException(nameof(index));
    }

    public override string ToString() => "()";
}

[StructLayout(LayoutKind.Sequential, Pack = 0)]
public struct NTuple1<T> : NTuple<T, NTuple1<T>> {
    [IgnoreDataMember] public int Count => 1;
    public T Item1;

    public NTuple1(){}
    public NTuple1(T item1) => Item1 = item1;

    public static NTuple1<T> Create(Func<int, T> itr) => new(itr(0));
    public static unsafe ref T GetElementReference(ref NTuple1<T> v, int index) => ref Unsafe.AsRef<T>((byte*)Unsafe.AsPointer(ref v) + Unsafe.SizeOf<T>()*index);

    public T Reduce(T identity, Func<T, T, T> reduction) => reduction(identity, Item1);
    
    public T this[int index] {
        get => index switch {
            0 => Item1,
            _ => throw new ArgumentOutOfRangeException(nameof(index))
        };
        set {
            Item1 = index switch {
                0 => value,
                _ => throw new ArgumentOutOfRangeException(nameof(index))
            };
        }
    }
    
    public override string ToString() => $"({Item1})";
}

[StructLayout(LayoutKind.Sequential, Pack = 0)]
public struct NTuple2<T> : NTuple<T, NTuple2<T>> {
    [IgnoreDataMember] public int Count => 2;
    public T Item1, Item2;
    
    public static NTuple2<T> Create(Func<int, T> itr) => new(itr(0), itr(1));
    public static unsafe ref T GetElementReference(ref NTuple2<T> v, int index) => ref Unsafe.AsRef<T>((byte*)Unsafe.AsPointer(ref v) + Unsafe.SizeOf<T>()*index);
    public T Reduce(T identity, Func<T, T, T> reduction) => reduction(reduction(identity, Item1), Item2);

    public NTuple2(){}
    public NTuple2(T item1, T item2) {
        Item1 = item1;
        Item2 = item2;
    }
    
    public T this[int index] {
        get => index switch {
            0 => Item1,
            1 => Item2,
            _ => throw new ArgumentOutOfRangeException(nameof(index))
        };
        set {
            switch(index) {
                case 0:
                    Item1 = value;
                    break;
                case 1:
                    Item2 = value;
                    break;
                default:
                    throw new ArgumentOutOfRangeException(nameof(index));
            };
        }
    }
    public override string ToString() => $"({Item1},{Item2})";
}

[StructLayout(LayoutKind.Sequential, Pack = 0)]
public struct NTuple3<T> : NTuple<T, NTuple3<T>> {
    [IgnoreDataMember] public int Count => 3;
    public T Item1, Item2, Item3;

    public NTuple3(T item1, T item2, T item3) {
        Item1 = item1;
        Item2 = item2;
        Item3 = item3;
    }
    
    public static NTuple3<T> Create(Func<int, T> itr) => new(itr(0), itr(1), itr(2));
    public static unsafe ref T GetElementReference(ref NTuple3<T> v, int index) => ref Unsafe.AsRef<T>((byte*)Unsafe.AsPointer(ref v) + Unsafe.SizeOf<T>()*index);
    public T Reduce(T identity, Func<T, T, T> reduction) => reduction(reduction(reduction(identity, Item1), Item2), Item3);
    
    public T this[int index] {
        get => index switch {
            0 => Item1,
            1 => Item2,
            2 => Item3,
            _ => throw new ArgumentOutOfRangeException(nameof(index))
        };
        set {
            switch(index) {
                case 0:
                    Item1 = value;
                    break;
                case 1:
                    Item2 = value;
                    break;
                case 2:
                    Item3 = value;
                    break;
                default:
                    throw new ArgumentOutOfRangeException(nameof(index));
            };
        }
    }
    public override string ToString() => $"({Item1},{Item2},{Item3})";
}

This breaks the memory layout of the static tuples but may still be useful for super high dimensional data (the classical approach to NTuples):

public struct NTupleArray<T, N> : NTuple<T, NTupleArray<T, N>> where N: IValue<int> {
    [IgnoreDataMember] public int Count => Items.Length;
    public readonly T[] Items;

    public NTupleArray(params T[] items) => Items = items;

    public static NTupleArray<T, N> Create(Func<int, T> itr) {
        var array = new T[N.Value];
        for (var i = 0; i < N.Value; i++)
            array[i] = itr(i);
        return new(array);
    }
    public static ref T GetElementReference(ref NTupleArray<T, N> v, int index) => ref v.Items[index];

    public T this[int index] {
        get => Items[index];
        set => Items[index] = value;
    }
    
    public override string ToString() => $"({string.Join(",", Items)})";
}

public interface IValue<T>{
    public static abstract T Value { get; }
}

public interface IOneThousand : IValue<int> {
    public static virtual int Value => 1000;
}

P.S You may have noticed I have been trying to do practices that make C# a little more julian (generics wise) ;).

[huoyaoyuan]

Indeed related to #111973 and #89730.

[HyperSphereStudio]

I optimized the ntuple library to make it obvious to the JIT where to inline function calls... Here are the results from a (very) simple and non rigourous benchmark

BenchmarkDotNet v0.14.0, Windows 11 (10.0.26100.3194)
Intel Core Ultra 7 155U, 1 CPU, 14 logical and 12 physical cores
.NET SDK 8.0.309
  [Host]     : .NET 8.0.13 (8.0.1325.6609), X64 RyuJIT AVX2
  DefaultJob : .NET 8.0.13 (8.0.1325.6609), X64 RyuJIT AVX2

| Method                   | Mean       | Error     | StdDev    |
|------------------------- |-----------:|----------:|----------:|
| Dot_Product_NTUPLE3f     |   3.251 ns | 0.0948 ns | 0.1360 ns |
| Dot_Product_VECTOR3f     |   3.333 ns | 0.0971 ns | 0.1541 ns |
| Cross_Product_NTUPLE3f   |  14.450 ns | 0.3009 ns | 0.2955 ns |
| Cross_Product_VECTOR3f   |  22.280 ns | 0.4661 ns | 0.4132 ns |
| Add_Mult_Length_NTUPLE3f |  94.948 ns | 1.0157 ns | 0.9501 ns |
| Add_Mult_Length_Vector3f |  79.510 ns | 1.0932 ns | 1.0226 ns |
| Add_Mult_Length_NTUPLE4f | 112.097 ns | 0.9206 ns | 0.8611 ns |
| Add_Mult_Length_VECTOR4f |  77.905 ns | 0.8184 ns | 0.6834 ns |
| Add_Mult_Length_NTUPLE2f |  70.128 ns | 0.6374 ns | 0.5962 ns |
| Add_Mult_Length_VECTOR2f |  81.153 ns | 0.4275 ns | 0.3998 ns |

The results are pretty neck and neck with the current Vector2, Vector3 & Vector4 on release mode.

Without release mode, the ntuples are about 8x slower (to be expected)

Here are what some of the functions from the test look like:

   [Benchmark]
    public Vec<float, NTuple4<float>> Add_Mult_Length_NTUPLE4f() {
        var w = new Vec<float, NTuple4<float>>(new NTuple4<float>(3, 8, 7, 12));
        var k = new Vec<float, NTuple4<float>>(new NTuple4<float>(65, 2, 8, -5));
        var c = w + k;
        for (int i = 0; i < 10; i++) {
            c += w + k * i;
            w += i * c.Length;
            k -= w.Length;
        }
        return c;
    }

  [Benchmark]
    public Vector4 Add_Mult_Length_VECTOR4f() {
        var w = new Vector4(3, 8, 7, 12);
        var k = new Vector4(65, 2, 8, -5);
        var c = w + k;
        for (int i = 0; i < 10; i++) {
            c += w + k * i;
            w += new Vector4(i * c.Length());
            k -= new Vector4(w.Length());
        }
        return c;
    }

[Benchmark]
    public float Dot_Product_NTUPLE3f() {
        var w = new Vec<float, NTuple3<float>>(new NTuple3<float>(3, 8, 7));
        var k = new Vec<float, NTuple3<float>>(new NTuple3<float>(65, 2, 8));
        float s = 0;
        for (var i = 0; i < 10; i++) {
            s += Vec.Dot(w, k);
        }
        return s;
    }
    
    [Benchmark]
    public float Dot_Product_VECTOR3f() {
        var w = new Vector3(3, 8, 7);
        var k = new Vector3(65, 2, 8);
        float s = 0;
        for (var i = 0; i < 10; i++) {
            s += Vector3.Dot(w, k);
        }
        return s;
    }

Here is what the hyperdimensional vec class looks like now:

public static class Vec {
    public static T Dot<T, NT>(Vec<T, NT> a, Vec<T, NT> b) where T : IFloatingPointIeee754<T> where NT : NTuple<T, NT> {
        return NT.BinaryOperation<Multiply<T>>(a.Data, b.Data).Sum<T, NT>();
    }
    
    /*
     * Use Item#s rather than X, Y & Z for less JIT optimization
     */
    public static Vec<T, NTuple3<T>> Cross<T>(Vec<T, NTuple3<T>> a, Vec<T, NTuple3<T>> b)
        where T : IFloatingPointIeee754<T> =>
        new(new NTuple3<T>(
            a.Data.Item2 * b.Data.Item3 - a.Data.Item3 * b.Data.Item2,
            a.Data.Item3 * b.Data.Item1 - a.Data.Item1 * b.Data.Item3,
            a.Data.Item1 * b.Data.Item2 - a.Data.Item2 * b.Data.Item1));
}

public struct Vec<T, NT> where NT : NTuple<T, NT> where T : IFloatingPointIeee754<T> {
    public NT Data;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public Vec(T v)
    {
        Data = NT.Uninitialized;
        NTuple<T, NT>.AsSpan(ref Data).Fill(v);
    }

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public Vec(NT data) => Data = data;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public Vec(Vec<T, NT> v) => Data = v.Data;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> Unit(int n)
    {
        var w = NT.Zeroed;
        w[n] = T.One;
        return new(w);
    }

    public static Vec<T, NT> UnitX { get; } = NT.Zeroed.Count > 0 ? Unit(0) : default;
    public static Vec<T, NT> UnitY { get; } = NT.Zeroed.Count > 0 ? Unit(1) : default;
    public static Vec<T, NT> UnitZ { get; } = NT.Zeroed.Count > 0 ? Unit(2) : default;
    public static Vec<T, NT> Zero { get; } = new(T.Zero);

    public static Span<T> AsSpan(ref Vec<T, NT> v) => NTuple<T, NT>.AsSpan(ref v.Data);
    [IgnoreDataMember] public T LengthSquared => NT.Reduce<EuclidLengthSquared<T>>(T.AdditiveIdentity, Data);
    [IgnoreDataMember] public T Length => T.Sqrt(LengthSquared);
    [IgnoreDataMember] public T ReciprocalSqrt => T.ReciprocalSqrtEstimate(LengthSquared);

    [IgnoreDataMember]
    public T X {
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        get => Data[0];
        
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        set => Data[0] = value;
    }

    [IgnoreDataMember]
    public T Y
    {
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        get => Data[1];
        
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        set => Data[1] = value;
    }

    [IgnoreDataMember]
    public T Z {
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        get => Data[2];
        
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        set => Data[2] = value;
    }

    [IgnoreDataMember]
    public int Dimension
    {
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        get => Data.Count;
    }

    public T this[int index]
    {
        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        get => Data[index];

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        set => Data[index] = value;
    }

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator +(Vec<T, NT> a, Vec<T, NT> b) =>
        new(NTuple.BinaryElementwiseAdd<T, NT>(a.Data, b.Data));

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator -(Vec<T, NT> v) => v * -T.One;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator -(Vec<T, NT> a, Vec<T, NT> b) => a + -b;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator +(Vec<T, NT> a, T b) => new(NTuple.ScalarAdd(a.Data, b));

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator +(T a, Vec<T, NT> b) => new(NTuple.ScalarAdd(b.Data, a));

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator -(T a, Vec<T, NT> b) => -b + a;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator -(Vec<T, NT> a, T b) => a + -b;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator *(T a, Vec<T, NT> b) => new(NTuple.ScalarMultiply(b.Data, a));

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator *(Vec<T, NT> a, T b) => b * a;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    public static Vec<T, NT> operator /(Vec<T, NT> a, T b) => a * T.ReciprocalEstimate(b);

    public override string? ToString() => Data.ToString();
}