You’re an expert programmer with a detailed understanding of high quality coding practices, a technical focus and numerical computing.You always review the code I give to you before you answer to check if there would be code improvements needed to apply first before proceeding with your answer, and when you make code improvements consider the change weight and how it would effect working with large Tensors<float> or Tensor<Double> for large neural networks in real world application. You will think and plan step by step carefully for best possible answer with best performance and least memory footprint,  use ParallelizationThreshold when preferred or required for best practices. you will produce fully complete and working code. 
Programming language : c#
Application Type: class library
Describe the application in full (what does it do?): Neural network library

 public partial class Tensor<T> : IDisposable, ICloneable, IEquatable<Tensor<T>>, IEnumerable<T> where T : unmanaged, INumber<T>, IExponentialFunctions<T>, ITrigonometricFunctions<T>, ILogarithmicFunctions<T>, IFloatingPointIeee754<T>, IRootFunctions<T> { #region ElementWiseOperation Operation public Tensor<T> ElementWiseOperation(Tensor<T> other, Func<T, T, T> operation) { EnsureNotDisposed(); other.EnsureNotDisposed(); var broadcastInfo = GetOrAddBroadcastCache(shape, other.shape, other.strides); int[] resultShape = broadcastInfo.ResultShape.ToArray(); int resultSize = resultShape.Aggregate(1, (a, b) => a * b); Tensor<T> result = new Tensor<T>(resultShape); var stridesA = broadcastInfo.BroadcastStridesA.ToArray(); var stridesB = broadcastInfo.BroadcastStridesB.ToArray(); var dataA = data; var dataB = other.data; var dataResult = result.data; var offsetA = offset; var offsetB = other.offset; var offsetResult = result.offset; if (resultSize >= ParallelizationThreshold) { Parallel.For(0, resultSize, i => { int indexA = 0; int indexB = 0; int idx = i; for (int dim = resultShape.Length - 1; dim >= 0; dim--) { int currentIndex = idx % resultShape[dim]; idx /= resultShape[dim]; indexA += stridesA[dim] * currentIndex; indexB += stridesB[dim] * currentIndex; } dataResult[offsetResult + i] = operation(dataA[offsetA + indexA], dataB[offsetB + indexB]); }); } else { for (int i = 0; i < resultSize; i++) { int indexA = 0; int indexB = 0; int idx = i; for (int dim = resultShape.Length - 1; dim >= 0; dim--) { int currentIndex = idx % resultShape[dim]; idx /= resultShape[dim]; indexA += stridesA[dim] * currentIndex; indexB += stridesB[dim] * currentIndex; } dataResult[offsetResult + i] = operation(dataA[offsetA + indexA], dataB[offsetB + indexB]); } } return result; } public Tensor<T> Add(Tensor<T> other) { return ElementWiseOperation(other, (x, y) => x + y); } public Tensor<T> Subtract(Tensor<T> other) { return ElementWiseOperation(other, (x, y) => x - y); } public Tensor<T> Multiply(Tensor<T> other) { return ElementWiseOperation(other, (x, y) => x * y); } public Tensor<T> Divide(Tensor<T> other) { return ElementWiseOperation(other, (x, y) => x / y); } public static Tensor<T> operator +(Tensor<T> a, Tensor<T> b) { return a.Add(b); } public static Tensor<T> operator -(Tensor<T> a, Tensor<T> b) { return a.Subtract(b); } public static Tensor<T> operator *(Tensor<T> a, Tensor<T> b) { return a.Multiply(b); } public static Tensor<T> operator /(Tensor<T> a, Tensor<T> b) { return a.Divide(b); } #endregion #region ScalarOperation public Tensor<T> ScalarOperation(Func<T, T> operation) { EnsureNotDisposed(); T[] resultData = ArrayPool<T>.Shared.Rent(totalSize); var offsetResult = 0; if (totalSize >= ParallelizationThreshold) { Parallel.For(0, totalSize, i => { resultData[offsetResult + i] = operation(data[offset + i]); }); } else { for (int i = 0; i < totalSize; i++) { resultData[offsetResult + i] = operation(data[offset + i]); } } return new Tensor<T>(resultData, (int[])shape.Clone(), (int[])strides.Clone(), totalSize, isRowMajor, disposed: false, offset: offsetResult); } public Tensor<T> AddScalar(T scalar) { return ScalarOperation((x) => x + scalar); } public Tensor<T> SubtractScalar(T scalar) { return ScalarOperation((x) => x - scalar); } public Tensor<T> MultiplyScalar(T scalar) { return ScalarOperation((x) => x * scalar); } public Tensor<T> DivideScalar(T scalar) { return ScalarOperation((x) => x / scalar); } public static Tensor<T> operator +(Tensor<T> a, T scalar) { return a.AddScalar(scalar); } public static Tensor<T> operator +(T scalar, Tensor<T> a) { return a.AddScalar(scalar); } public static Tensor<T> operator -(Tensor<T> a, T scalar) { return a.SubtractScalar(scalar); } public static Tensor<T> operator -(T scalar, Tensor<T> a) { return a.Negate().AddScalar(scalar); } public static Tensor<T> operator *(Tensor<T> a, T scalar) { return a.MultiplyScalar(scalar); } public static Tensor<T> operator *(T scalar, Tensor<T> a) { return a.MultiplyScalar(scalar); } public static Tensor<T> operator /(Tensor<T> a, T scalar) { return a.DivideScalar(scalar); } public static Tensor<T> operator /(T scalar, Tensor<T> a) { return a.Reciprocal().MultiplyScalar(scalar); } #endregion #region In Place Operations public void ElementWiseInPlaceOperation(Tensor<T> other, Func<T, T, T> operation) { EnsureNotDisposed(); other.EnsureNotDisposed(); if (!shape.SequenceEqual(other.shape)) throw new ArgumentException($"Shapes must match for ElementWiseInPlaceOperation."); if (totalSize >= ParallelizationThreshold) { Parallel.For(0, totalSize, i => { data[i] = operation(data[i], other.data[i]); }); } else { for (int i = 0; i < totalSize; i++) { data[i] = operation(data[i], other.data[i]); } } } public void AddInPlace(Tensor<T> other) { ElementWiseInPlaceOperation(other, (x, y) => x + y); } public void SubtractInPlace(Tensor<T> other) { ElementWiseInPlaceOperation(other, (x, y) => x - y); } public void MultiplyInPlace(Tensor<T> other) { ElementWiseInPlaceOperation(other, (x, y) => x * y); } public void DivideInPlace(Tensor<T> other) { ElementWiseInPlaceOperation(other, (x, y) => x / y); } public void ScalarInPlaceOperation(Func<T, T> operation) { EnsureNotDisposed(); if (totalSize >= ParallelizationThreshold) { Parallel.For(0, totalSize, i => { data[i] = operation(data[i]); }); } else { for (int i = 0; i < totalSize; i++) { data[i] = operation(data[i]); } } } public void MultiplyScalarInPlace(T scalar) { ScalarInPlaceOperation(x => x * scalar); } public void AddScalarInPlace(T scalar) { ScalarInPlaceOperation(x => x + scalar); } public void SubtractScalarInPlace(T scalar) { ScalarInPlaceOperation(x => x - scalar); } public void DivideScalarInPlace(T scalar) { ScalarInPlaceOperation(x => x / scalar); } #endregion #region MatMul public Tensor<T> MatMulSimple(Tensor<T> other) { EnsureNotDisposed(); other.EnsureNotDisposed(); if (Rank != 2 || other.Rank != 2) throw new InvalidOperationException("MatMulSimple requires both tensors to be 2D."); int m = shape[0]; int n = shape[1]; int p = other.shape[1]; if (n != other.shape[0]) throw new InvalidOperationException("Inner dimensions must match for MatMulSimple."); var result = new Tensor<T>(m, p); var dataA = data; var dataB = other.data; var dataResult = result.data; for (int i = 0; i < m; i++) { for (int j = 0; j < p; j++) { T sum = T.Zero; for (int k = 0; k < n; k++) { sum += dataA[i * n + k] * dataB[k * p + j]; } dataResult[i * p + j] = sum; } } return result; } public Tensor<T> MatMul(Tensor<T> other) { EnsureNotDisposed(); other.EnsureNotDisposed(); var aShape = shape; var bShape = other.shape; if (Rank == 1 && other.Rank == 1) { if (shape[0] != other.shape[0]) throw new ArgumentException("Shapes are not aligned for dot product."); T sum = T.Zero; for (int i = 0; i < shape[0]; i++) { sum += this[i] * other[i]; } var resultTensor = new Tensor<T>(); resultTensor.data[0] = sum; return resultTensor; } if (Rank == 1) { var aReshaped = Reshape(1, shape[0]); var partialResult = aReshaped.MatMul(other); return partialResult.Reshape(partialResult.shape.Skip(1).ToArray()); } if (other.Rank == 1) { var bReshaped = other.Reshape(other.shape[0], 1); var partialResult = MatMul(bReshaped); return partialResult.Reshape(partialResult.shape.Take(Rank - 1).ToArray()); } int n = shape[Rank - 1]; int k = other.shape[other.Rank - 2]; if (n != k) throw new ArgumentException("Inner dimensions do not match for MatMul."); int[] aBatchShape = shape.Take(Rank - 2).ToArray(); int[] bBatchShape = other.shape.Take(other.Rank - 2).ToArray(); int[] resultBatchShape = BroadcastShapes(aBatchShape, bBatchShape); int m = shape[Rank - 2]; int p = other.shape[other.Rank - 1]; int[] resultShape = resultBatchShape.Concat(new int[] { m, p }).ToArray(); Tensor<T> result = new Tensor<T>(resultShape); var dataA = data; var dataB = other.data; var dataResult = result.data; var offsetA = offset; var offsetB = other.offset; var offsetResult = result.offset; int batchSize = resultBatchShape.Aggregate(1, (a, b) => a * b); var aStrides = strides; var bStrides = other.strides; var resultStrides = result.strides; if (batchSize >= ParallelizationThreshold) { Parallel.For(0, batchSize, batchIndex => { var batchIndices = UnflattenIndex(batchIndex, resultBatchShape); int aBatchOffset = GetBatchOffset(batchIndices, aBatchShape, aStrides); int bBatchOffset = GetBatchOffset(batchIndices, bBatchShape, bStrides); int resultBatchOffset = GetBatchOffset(batchIndices, resultBatchShape, resultStrides); for (int i = 0; i < m; i++) { for (int j = 0; j < p; j++) { T sum = T.Zero; for (int l = 0; l < n; l++) { int aIndex = offsetA + aBatchOffset + i * aStrides[Rank - 2] + l * aStrides[Rank - 1]; int bIndex = offsetB + bBatchOffset + l * bStrides[other.Rank - 2] + j * bStrides[other.Rank - 1]; sum += dataA[aIndex] * dataB[bIndex]; } int resultIndex = offsetResult + resultBatchOffset + i * resultStrides[result.Rank - 2] + j * resultStrides[result.Rank - 1]; dataResult[resultIndex] = sum; } } }); } else { for (int batchIndex = 0; batchIndex < batchSize; batchIndex++) { var batchIndices = UnflattenIndex(batchIndex, resultBatchShape); int aBatchOffset = GetBatchOffset(batchIndices, aBatchShape, aStrides); int bBatchOffset = GetBatchOffset(batchIndices, bBatchShape, bStrides); int resultBatchOffset = GetBatchOffset(batchIndices, resultBatchShape, resultStrides); for (int i = 0; i < m; i++) { for (int j = 0; j < p; j++) { T sum = T.Zero; for (int l = 0; l < n; l++) { int aIndex = offsetA + aBatchOffset + i * aStrides[Rank - 2] + l * aStrides[Rank - 1]; int bIndex = offsetB + bBatchOffset + l * bStrides[other.Rank - 2] + j * bStrides[other.Rank - 1]; sum += dataA[aIndex] * dataB[bIndex]; } int resultIndex = offsetResult + resultBatchOffset + i * resultStrides[result.Rank - 2] + j * resultStrides[result.Rank - 1]; dataResult[resultIndex] = sum; } } } } return result; } public Tensor<T> MatMulBatchOffset(Tensor<T> other) { return MatMul(other); } #endregion public Tensor<T> Negate() { return ScalarOperation((x) => -x); } public Tensor<T> Reciprocal() { return ScalarOperation((x) => T.One / x); } public Tensor<T> Log() { return ScalarOperation(x => T.Log(x)); } public Tensor<T> OneSubtract() { EnsureNotDisposed(); return ScalarOperation(x => T.One - x); } public void Fill(T value) { EnsureNotDisposed(); if (totalSize >= ParallelizationThreshold) { Parallel.For(0, totalSize, i => { data[i + offset] = value; }); } else { for (int i = 0; i < totalSize; i++) { data[i + offset] = value; } } } }

Task1: create mstest for the current given methods in Tensor<T> class.