using System;

using System.Collections;

using System.Collections.Generic;

using System.IO;

using System.Linq;

using System.Runtime.InteropServices;

using System.Text;

using Numpy.Models;

using Microsoft.VisualStudio.TestTools.UnitTesting;

using Assert = NUnit.Framework.Assert;

namespace Numpy.UnitTest

[TestClass]

public class NumPy_linalgTest : BaseTestCase

{

[TestMethod]

public void dotTest()

}

[TestMethod]

public void multi_dotTest()

}

[TestMethod]

public void vdotTest()

}

[TestMethod]

public void innerTest()

}

[TestMethod]

public void outerTest()

}

[TestMethod]

public void matmulTest()

{

// For 2-D arrays it is the matrix product:

// >>> a = np.array([[1, 0],

// ... [0, 1]])

// >>> b = np.array([[4, 1],

// ... [2, 2]]

// >>> np.matmul(a, b)

// array([[4, 1],

// [2, 2]])

//

#if TODO

var given= a = np.array({{1, 0},;

var expected=

"... [0, 1]])";

Assert.AreEqual(expected, given.repr);

given= b = np.array({{4, 1}, ;

expected=

"... [2, 2]]";

Assert.AreEqual(expected, given.repr);

given= np.matmul(a, b);

expected=

"array([[4, 1],\n" +

" [2, 2]])";

Assert.AreEqual(expected, given.repr);

#endif

// For 2-D mixed with 1-D, the result is the usual.

// >>> a = np.array([[1, 0],

// ... [0, 1]]

// >>> b = np.array([1, 2])

// >>> np.matmul(a, b)

// array([1, 2])

// >>> np.matmul(b, a)

// array([1, 2])

//

#if TODO

given= a = np.array({{1, 0},;

expected=

"... [0, 1]]";

Assert.AreEqual(expected, given.repr);

given= b = np.array({1, 2});

given= np.matmul(a, b);

expected=

"array([1, 2])";

Assert.AreEqual(expected, given.repr);

given= np.matmul(b, a);

expected=

"array([1, 2])";

Assert.AreEqual(expected, given.repr);

#endif

// Broadcasting is conventional for stacks of arrays

// >>> a = np.arange(2 * 2 * 4).reshape((2, 2, 4))

// >>> b = np.arange(2 * 2 * 4).reshape((2, 4, 2))

// >>> np.matmul(a,b).shape

// (2, 2, 2)

// >>> np.matmul(a, b)[0, 1, 1]

// 98

// >>> sum(a[0, 1, :] * b[0 , :, 1])

// 98

//

#if TODO

given= a = np.arange(2 * 2 * 4).reshape((2, 2, 4));

given= b = np.arange(2 * 2 * 4).reshape((2, 4, 2));

given= np.matmul(a,b).shape;

expected=

"(2, 2, 2)";

Assert.AreEqual(expected, given.repr);

given= np.matmul(a, b){0, 1, 1};

expected=

"98";

Assert.AreEqual(expected, given.repr);

given= sum(a[0, 1, :] * b[0 , :, 1]);

expected=

"98";

Assert.AreEqual(expected, given.repr);

#endif

// Vector, vector returns the scalar inner product, but neither argument

// is complex-conjugated:

// >>> np.matmul([2j, 3j], [2j, 3j])

// (-13+0j)

//

#if TODO

given= np.matmul({2j, 3j}, {2j, 3j});

expected=

"(-13+0j)";

Assert.AreEqual(expected, given.repr);

#endif

// Scalar multiplication raises an error.

// >>> np.matmul([1,2], 3)

// Traceback (most recent call last):

// ...

// ValueError: matmul: Input operand 1 does not have enough dimensions ...

//

#if TODO

given= np.matmul({1,2}, 3);

expected=

"Traceback (most recent call last):\n" +

"...\n" +

"ValueError: matmul: Input operand 1 does not have enough dimensions ...";

Assert.AreEqual(expected, given.repr);

#endif

}

[TestMethod]

public void tensordotTest()

{

// A “traditional” example:

// >>> a = np.arange(60.).reshape(3,4,5)

// >>> b = np.arange(24.).reshape(4,3,2)

// >>> c = np.tensordot(a,b, axes=([1,0],[0,1]))

// >>> c.shape

// (5, 2)

// >>> c

// array([[ 4400., 4730.],

// [ 4532., 4874.],

// [ 4664., 5018.],

// [ 4796., 5162.],

// [ 4928., 5306.]])

// >>> # A slower but equivalent way of computing the same...

// >>> d = np.zeros((5,2))

// >>> for i in range(5):

// ... for j in range(2):

// ... for k in range(3):

// ... for n in range(4):

// ... d[i,j] += a[k,n,i] * b[n,k,j]

// >>> c == d

// array([[ True, True],

// [ True, True],

// [ True, True],

// [ True, True],

// [ True, True]])

//

#if TODO

var given= a = np.arange(60.).reshape(3,4,5);

given= b = np.arange(24.).reshape(4,3,2);

given= c = np.tensordot(a,b, axes=({1,0},{0,1}));

given= c.shape;

var expected=

"(5, 2)";

Assert.AreEqual(expected, given.repr);

given= c;

expected=

"array([[ 4400., 4730.],\n" +

" [ 4532., 4874.],\n" +

" [ 4664., 5018.],\n" +

" [ 4796., 5162.],\n" +

" [ 4928., 5306.]])";

Assert.AreEqual(expected, given.repr);

given= # A slower but equivalent way of computing the same...;

given= d = np.zeros((5,2));

given= for i in range(5):;

expected=

"... for j in range(2):\n" +

"... for k in range(3):\n" +

"... for n in range(4):\n" +

"... d[i,j] += a[k,n,i] * b[n,k,j]";

Assert.AreEqual(expected, given.repr);

given= c == d;

expected=

"array([[ True, True],\n" +

" [ True, True],\n" +

" [ True, True],\n" +

" [ True, True],\n" +

" [ True, True]])";

Assert.AreEqual(expected, given.repr);

#endif

// An extended example taking advantage of the overloading of + and *:

// >>> a = np.array(range(1, 9))

// >>> a.shape = (2, 2, 2)

// >>> A = np.array(('a', 'b', 'c', 'd'), dtype=object)

// >>> A.shape = (2, 2)

// >>> a; A

// array([[[1, 2],

// [3, 4]],

// [[5, 6],

// [7, 8]]])

// array([[a, b],

// [c, d]], dtype=object)

//

#if TODO

given= a = np.array(range(1, 9));

given= a.shape = (2, 2, 2);

given= A = np.array(('a', 'b', 'c', 'd'), dtype=object);

given= A.shape = (2, 2);

given= a; A;

expected=

"array([[[1, 2],\n" +

" [3, 4]],\n" +

" [[5, 6],\n" +

" [7, 8]]])\n" +

"array([[a, b],\n" +

" [c, d]], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A) # third argument default is 2 for double-contraction

// array([abbcccdddd, aaaaabbbbbbcccccccdddddddd], dtype=object)

//

#if TODO

given= np.tensordot(a, A) # third argument default is 2 for double-contraction;

expected=

"array([abbcccdddd, aaaaabbbbbbcccccccdddddddd], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, 1)

// array([[[acc, bdd],

// [aaacccc, bbbdddd]],

// [[aaaaacccccc, bbbbbdddddd],

// [aaaaaaacccccccc, bbbbbbbdddddddd]]], dtype=object)

//

#if TODO

given= np.tensordot(a, A, 1);

expected=

"array([[[acc, bdd],\n" +

" [aaacccc, bbbdddd]],\n" +

" [[aaaaacccccc, bbbbbdddddd],\n" +

" [aaaaaaacccccccc, bbbbbbbdddddddd]]], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, 0) # tensor product (result too long to incl.)

// array([[[[[a, b],

// [c, d]],

// ...

//

#if TODO

given= np.tensordot(a, A, 0) # tensor product (result too long to incl.);

expected=

"array([[[[[a, b],\n" +

" [c, d]],\n" +

" ...";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, (0, 1))

// array([[[abbbbb, cddddd],

// [aabbbbbb, ccdddddd]],

// [[aaabbbbbbb, cccddddddd],

// [aaaabbbbbbbb, ccccdddddddd]]], dtype=object)

//

#if TODO

given= np.tensordot(a, A, (0, 1));

expected=

"array([[[abbbbb, cddddd],\n" +

" [aabbbbbb, ccdddddd]],\n" +

" [[aaabbbbbbb, cccddddddd],\n" +

" [aaaabbbbbbbb, ccccdddddddd]]], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, (2, 1))

// array([[[abb, cdd],

// [aaabbbb, cccdddd]],

// [[aaaaabbbbbb, cccccdddddd],

// [aaaaaaabbbbbbbb, cccccccdddddddd]]], dtype=object)

//

#if TODO

given= np.tensordot(a, A, (2, 1));

expected=

"array([[[abb, cdd],\n" +

" [aaabbbb, cccdddd]],\n" +

" [[aaaaabbbbbb, cccccdddddd],\n" +

" [aaaaaaabbbbbbbb, cccccccdddddddd]]], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, ((0, 1), (0, 1)))

// array([abbbcccccddddddd, aabbbbccccccdddddddd], dtype=object)

//

#if TODO

given= np.tensordot(a, A, ((0, 1), (0, 1)));

expected=

"array([abbbcccccddddddd, aabbbbccccccdddddddd], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

// >>> np.tensordot(a, A, ((2, 1), (1, 0)))

// array([acccbbdddd, aaaaacccccccbbbbbbdddddddd], dtype=object)

//

#if TODO

given= np.tensordot(a, A, ((2, 1), (1, 0)));

expected=

"array([acccbbdddd, aaaaacccccccbbbbbbdddddddd], dtype=object)";

Assert.AreEqual(expected, given.repr);

#endif

}

[TestMethod]

public void einsumTest()

{

// >>> a = np.arange(25).reshape(5,5)

// >>> b = np.arange(5)

// >>> c = np.arange(6).reshape(2,3)

//

var a = np.arange(25).reshape(5, 5);

var b = np.arange(5);

var c = np.arange(6).reshape(2, 3);

// Trace of a matrix:

// >>> np.einsum('ii', a)

// 60

// >>> np.einsum(a, [0,0])

// 60

// >>> np.trace(a)

// 60

//

var given = np.einsum("ii", a);

var expected =

"60";

#if TODO

Assert.AreEqual(expected, given.repr);

given = np.einsum(new[] { a, (NDarray)np.array(0, 0) });

expected =

"60";

Assert.AreEqual(expected, given.repr);

#endif

given = np.trace(a);

expected =

"60";

Assert.AreEqual(expected, given.repr);

// Extract the diagonal (requires explicit form):

// >>> np.einsum('ii->i', a)

// array([ 0, 6, 12, 18, 24])

// >>> np.einsum(a, [0,0], [0])

// array([ 0, 6, 12, 18, 24])

// >>> np.diag(a)

// array([ 0, 6, 12, 18, 24])

//

given = np.einsum("ii->i", a);

expected =

"array([ 0, 6, 12, 18, 24])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(a, {0,0}, {0});

expected=

"array([ 0, 6, 12, 18, 24])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.diag(a);

expected =

"array([ 0, 6, 12, 18, 24])";

Assert.AreEqual(expected, given.repr);

// Sum over an axis (requires explicit form):

// >>> np.einsum('ij->i', a)

// array([ 10, 35, 60, 85, 110])

// >>> np.einsum(a, [0,1], [0])

// array([ 10, 35, 60, 85, 110])

// >>> np.sum(a, axis=1)

// array([ 10, 35, 60, 85, 110])

//

given = np.einsum("ij->i", a);

expected =

"array([ 10, 35, 60, 85, 110])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(a, {0,1}, {0});

expected=

"array([ 10, 35, 60, 85, 110])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.sum(a, axis: 1);

expected =

"array([ 10, 35, 60, 85, 110])";

Assert.AreEqual(expected, given.repr);

// For higher dimensional arrays summing a single axis can be done with ellipsis:

// >>> np.einsum('...j->...', a)

// array([ 10, 35, 60, 85, 110])

// >>> np.einsum(a, [Ellipsis,1], [Ellipsis])

// array([ 10, 35, 60, 85, 110])

//

given = np.einsum("...j->...", a);

expected =

"array([ 10, 35, 60, 85, 110])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(a, {Ellipsis,1}, {Ellipsis});

expected=

"array([ 10, 35, 60, 85, 110])";

Assert.AreEqual(expected, given.repr);

#endif

// Compute a matrix transpose, or reorder any number of axes:

// >>> np.einsum('ji', c)

// array([[0, 3],

// [1, 4],

// [2, 5]])

// >>> np.einsum('ij->ji', c)

// array([[0, 3],

// [1, 4],

// [2, 5]])

// >>> np.einsum(c, [1,0])

// array([[0, 3],

// [1, 4],

// [2, 5]])

// >>> np.transpose(c)

// array([[0, 3],

// [1, 4],

// [2, 5]])

//

given = np.einsum("ji", c);

expected =

"array([[0, 3],\n" +

" [1, 4],\n" +

" [2, 5]])";

Assert.AreEqual(expected, given.repr);

given = np.einsum("ij->ji", c);

expected =

"array([[0, 3],\n" +

" [1, 4],\n" +

" [2, 5]])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(c, {1,0});

expected=

"array([[0, 3],\n" +

" [1, 4],\n" +

" [2, 5]])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.transpose(c);

expected =

"array([[0, 3],\n" +

" [1, 4],\n" +

" [2, 5]])";

Assert.AreEqual(expected, given.repr);

// Vector inner products:

// >>> np.einsum('i,i', b, b)

// 30

// >>> np.einsum(b, [0], b, [0])

// 30

// >>> np.inner(b,b)

// 30

//

given = np.einsum("i,i", b, b);

expected =

"30";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(b, {0}, b, {0});

expected=

"30";

Assert.AreEqual(expected, given.repr);

#endif

given = np.inner(b, b);

expected =

"30";

Assert.AreEqual(expected, given.repr);

// Matrix vector multiplication:

// >>> np.einsum('ij,j', a, b)

// array([ 30, 80, 130, 180, 230])

// >>> np.einsum(a, [0,1], b, [1])

// array([ 30, 80, 130, 180, 230])

// >>> np.dot(a, b)

// array([ 30, 80, 130, 180, 230])

// >>> np.einsum('...j,j', a, b)

// array([ 30, 80, 130, 180, 230])

//

given = np.einsum("ij,j", a, b);

expected =

"array([ 30, 80, 130, 180, 230])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(a, {0,1}, b, {1});

expected=

"array([ 30, 80, 130, 180, 230])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.dot(a, b);

expected =

"array([ 30, 80, 130, 180, 230])";

Assert.AreEqual(expected, given.repr);

given = np.einsum("...j,j", a, b);

expected =

"array([ 30, 80, 130, 180, 230])";

Assert.AreEqual(expected, given.repr);

// Broadcasting and scalar multiplication:

// >>> np.einsum('..., ...', 3, c)

// array([[ 0, 3, 6],

// [ 9, 12, 15]])

// >>> np.einsum(',ij', 3, c)

// array([[ 0, 3, 6],

// [ 9, 12, 15]])

// >>> np.einsum(3, [Ellipsis], c, [Ellipsis])

// array([[ 0, 3, 6],

// [ 9, 12, 15]])

// >>> np.multiply(3, c)

// array([[ 0, 3, 6],

// [ 9, 12, 15]])

//

given = np.einsum("..., ...", np.asarray(3), c);

expected =

"array([[ 0, 3, 6],\n" +

" [ 9, 12, 15]])";

Assert.AreEqual(expected, given.repr);

given = np.einsum(",ij", np.asarray(3), c);

expected =

"array([[ 0, 3, 6],\n" +

" [ 9, 12, 15]])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(3, {Ellipsis}, c, {Ellipsis});

expected=

"array([[ 0, 3, 6],\n" +

" [ 9, 12, 15]])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.multiply(np.asarray(3), c);

expected =

"array([[ 0, 3, 6],\n" +

" [ 9, 12, 15]])";

Assert.AreEqual(expected, given.repr);

// Vector outer product:

// >>> np.einsum('i,j', np.arange(2)+1, b)

// array([[0, 1, 2, 3, 4],

// [0, 2, 4, 6, 8]])

// >>> np.einsum(np.arange(2)+1, [0], b, [1])

// array([[0, 1, 2, 3, 4],

// [0, 2, 4, 6, 8]])

// >>> np.outer(np.arange(2)+1, b)

// array([[0, 1, 2, 3, 4],

// [0, 2, 4, 6, 8]])

//

given = np.einsum("i,j", np.arange(2) + 1, b);

expected =

"array([[0, 1, 2, 3, 4],\n" +

" [0, 2, 4, 6, 8]])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(np.arange(2)+1, {0}, b, {1});

expected=

"array([[0, 1, 2, 3, 4],\n" +

" [0, 2, 4, 6, 8]])";

Assert.AreEqual(expected, given.repr);

#endif

given = np.outer(np.arange(2) + 1, b);

expected =

"array([[0, 1, 2, 3, 4],\n" +

" [0, 2, 4, 6, 8]])";

Assert.AreEqual(expected, given.repr);

// Tensor contraction:

// >>> a = np.arange(60.).reshape(3,4,5)

// >>> b = np.arange(24.).reshape(4,3,2)

// >>> np.einsum('ijk,jil->kl', a, b)

// array([[ 4400., 4730.],

// [ 4532., 4874.],

// [ 4664., 5018.],

// [ 4796., 5162.],

// [ 4928., 5306.]])

// >>> np.einsum(a, [0,1,2], b, [1,0,3], [2,3])

// array([[ 4400., 4730.],

// [ 4532., 4874.],

// [ 4664., 5018.],

// [ 4796., 5162.],

// [ 4928., 5306.]])

// >>> np.tensordot(a,b, axes=([1,0],[0,1]))

// array([[ 4400., 4730.],

// [ 4532., 4874.],

// [ 4664., 5018.],

// [ 4796., 5162.],

// [ 4928., 5306.]])

//

given = a = np.arange(60.0).reshape(3, 4, 5);

given = b = np.arange(24.0).reshape(4, 3, 2);

given = np.einsum("ijk,jil->kl", a, b);

expected =

"array([[4400., 4730.],\n" +

" [4532., 4874.],\n" +

" [4664., 5018.],\n" +

" [4796., 5162.],\n" +

" [4928., 5306.]])";

Assert.AreEqual(expected, given.repr);

#if TODO

given= np.einsum(a, {0,1,2}, b, {1,0,3}, {2,3});

expected=

"array([[ 4400., 4730.],\n" +

" [ 4532., 4874.],\n" +

" [ 4664., 5018.],\n" +

" [ 4796., 5162.],\n" +

" [ 4928., 5306.]])";

Assert.AreEqual(expected, given.repr);

given = np.tensordot(a,b, axes:new int[,]{{1,0},{0,1}});

expected=

"array([[ 4400., 4730.],\n" +