import random from sympy.combinatorics import Permutation from sympy.combinatorics.permutations import _af_invert from sympy.testing.pytest import raises from sympy.core.function import diff from sympy.core.symbol import symbols from sympy.functions.elementary.complexes import (adjoint, conjugate, transpose) from sympy.functions.elementary.exponential import (exp, log) from sympy.functions.elementary.trigonometric import (cos, sin) from sympy.tensor.array import Array, ImmutableDenseNDimArray, ImmutableSparseNDimArray, MutableSparseNDimArray from sympy.tensor.array.arrayop import tensorproduct, tensorcontraction, derive_by_array, permutedims, Flatten, \ tensordiagonal def test_import_NDimArray(): from sympy.tensor.array import NDimArray del NDimArray def test_tensorproduct(): x,y,z,t = symbols('x y z t') from sympy.abc import a,b,c,d assert tensorproduct() == 1 assert tensorproduct([x]) == Array([x]) assert tensorproduct([x], [y]) == Array([[x*y]]) assert tensorproduct([x], [y], [z]) == Array([[[x*y*z]]]) assert tensorproduct([x], [y], [z], [t]) == Array([[[[x*y*z*t]]]]) assert tensorproduct(x) == x assert tensorproduct(x, y) == x*y assert tensorproduct(x, y, z) == x*y*z assert tensorproduct(x, y, z, t) == x*y*z*t for ArrayType in [ImmutableDenseNDimArray, ImmutableSparseNDimArray]: A = ArrayType([x, y]) B = ArrayType([1, 2, 3]) C = ArrayType([a, b, c, d]) assert tensorproduct(A, B, C) == ArrayType([[[a*x, b*x, c*x, d*x], [2*a*x, 2*b*x, 2*c*x, 2*d*x], [3*a*x, 3*b*x, 3*c*x, 3*d*x]], [[a*y, b*y, c*y, d*y], [2*a*y, 2*b*y, 2*c*y, 2*d*y], [3*a*y, 3*b*y, 3*c*y, 3*d*y]]]) assert tensorproduct([x, y], [1, 2, 3]) == tensorproduct(A, B) assert tensorproduct(A, 2) == ArrayType([2*x, 2*y]) assert tensorproduct(A, [2]) == ArrayType([[2*x], [2*y]]) assert tensorproduct([2], A) == ArrayType([[2*x, 2*y]]) assert tensorproduct(a, A) == ArrayType([a*x, a*y]) assert tensorproduct(a, A, B) == ArrayType([[a*x, 2*a*x, 3*a*x], [a*y, 2*a*y, 3*a*y]]) assert tensorproduct(A, B, a) == ArrayType([[a*x, 2*a*x, 3*a*x], [a*y, 2*a*y, 3*a*y]]) assert tensorproduct(B, a, A) == ArrayType([[a*x, a*y], [2*a*x, 2*a*y], [3*a*x, 3*a*y]]) # tests for large scale sparse array for SparseArrayType in [ImmutableSparseNDimArray, MutableSparseNDimArray]: a = SparseArrayType({1:2, 3:4},(1000, 2000)) b = SparseArrayType({1:2, 3:4},(1000, 2000)) assert tensorproduct(a, b) == ImmutableSparseNDimArray({2000001: 4, 2000003: 8, 6000001: 8, 6000003: 16}, (1000, 2000, 1000, 2000)) def test_tensorcontraction(): from sympy.abc import a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,q,r,s,t,u,v,w,x B = Array(range(18), (2, 3, 3)) assert tensorcontraction(B, (1, 2)) == Array([12, 39]) C1 = Array([a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x], (2, 3, 2, 2)) assert tensorcontraction(C1, (0, 2)) == Array([[a + o, b + p], [e + s, f + t], [i + w, j + x]]) assert tensorcontraction(C1, (0, 2, 3)) == Array([a + p, e + t, i + x]) assert tensorcontraction(C1, (2, 3)) == Array([[a + d, e + h, i + l], [m + p, q + t, u + x]]) def test_derivative_by_array(): from sympy.abc import i, j, t, x, y, z bexpr = x*y**2*exp(z)*log(t) sexpr = sin(bexpr) cexpr = cos(bexpr) a = Array([sexpr]) assert derive_by_array(sexpr, t) == x*y**2*exp(z)*cos(x*y**2*exp(z)*log(t))/t assert derive_by_array(sexpr, [x, y, z]) == Array([bexpr/x*cexpr, 2*y*bexpr/y**2*cexpr, bexpr*cexpr]) assert derive_by_array(a, [x, y, z]) == Array([[bexpr/x*cexpr], [2*y*bexpr/y**2*cexpr], [bexpr*cexpr]]) assert derive_by_array(sexpr, [[x, y], [z, t]]) == Array([[bexpr/x*cexpr, 2*y*bexpr/y**2*cexpr], [bexpr*cexpr, bexpr/log(t)/t*cexpr]]) assert derive_by_array(a, [[x, y], [z, t]]) == Array([[[bexpr/x*cexpr], [2*y*bexpr/y**2*cexpr]], [[bexpr*cexpr], [bexpr/log(t)/t*cexpr]]]) assert derive_by_array([[x, y], [z, t]], [x, y]) == Array([[[1, 0], [0, 0]], [[0, 1], [0, 0]]]) assert derive_by_array([[x, y], [z, t]], [[x, y], [z, t]]) == Array([[[[1, 0], [0, 0]], [[0, 1], [0, 0]]], [[[0, 0], [1, 0]], [[0, 0], [0, 1]]]]) assert diff(sexpr, t) == x*y**2*exp(z)*cos(x*y**2*exp(z)*log(t))/t assert diff(sexpr, Array([x, y, z])) == Array([bexpr/x*cexpr, 2*y*bexpr/y**2*cexpr, bexpr*cexpr]) assert diff(a, Array([x, y, z])) == Array([[bexpr/x*cexpr], [2*y*bexpr/y**2*cexpr], [bexpr*cexpr]]) assert diff(sexpr, Array([[x, y], [z, t]])) == Array([[bexpr/x*cexpr, 2*y*bexpr/y**2*cexpr], [bexpr*cexpr, bexpr/log(t)/t*cexpr]]) assert diff(a, Array([[x, y], [z, t]])) == Array([[[bexpr/x*cexpr], [2*y*bexpr/y**2*cexpr]], [[bexpr*cexpr], [bexpr/log(t)/t*cexpr]]]) assert diff(Array([[x, y], [z, t]]), Array([x, y])) == Array([[[1, 0], [0, 0]], [[0, 1], [0, 0]]]) assert diff(Array([[x, y], [z, t]]), Array([[x, y], [z, t]])) == Array([[[[1, 0], [0, 0]], [[0, 1], [0, 0]]], [[[0, 0], [1, 0]], [[0, 0], [0, 1]]]]) # test for large scale sparse array for SparseArrayType in [ImmutableSparseNDimArray, MutableSparseNDimArray]: b = MutableSparseNDimArray({0:i, 1:j}, (10000, 20000)) assert derive_by_array(b, i) == ImmutableSparseNDimArray({0: 1}, (10000, 20000)) assert derive_by_array(b, (i, j)) == ImmutableSparseNDimArray({0: 1, 200000001: 1}, (2, 10000, 20000)) #https://github.com/sympy/sympy/issues/20655 U = Array([x, y, z]) E = 2 assert derive_by_array(E, U) == ImmutableDenseNDimArray([0, 0, 0]) def test_issue_emerged_while_discussing_10972(): ua = Array([-1,0]) Fa = Array([[0, 1], [-1, 0]]) po = tensorproduct(Fa, ua, Fa, ua) assert tensorcontraction(po, (1, 2), (4, 5)) == Array([[0, 0], [0, 1]]) sa = symbols('a0:144') po = Array(sa, [2, 2, 3, 3, 2, 2]) assert tensorcontraction(po, (0, 1), (2, 3), (4, 5)) == sa[0] + sa[108] + sa[111] + sa[124] + sa[127] + sa[140] + sa[143] + sa[16] + sa[19] + sa[3] + sa[32] + sa[35] assert tensorcontraction(po, (0, 1, 4, 5), (2, 3)) == sa[0] + sa[111] + sa[127] + sa[143] + sa[16] + sa[32] assert tensorcontraction(po, (0, 1), (4, 5)) == Array([[sa[0] + sa[108] + sa[111] + sa[3], sa[112] + sa[115] + sa[4] + sa[7], sa[11] + sa[116] + sa[119] + sa[8]], [sa[12] + sa[120] + sa[123] + sa[15], sa[124] + sa[127] + sa[16] + sa[19], sa[128] + sa[131] + sa[20] + sa[23]], [sa[132] + sa[135] + sa[24] + sa[27], sa[136] + sa[139] + sa[28] + sa[31], sa[140] + sa[143] + sa[32] + sa[35]]]) assert tensorcontraction(po, (0, 1), (2, 3)) == Array([[sa[0] + sa[108] + sa[124] + sa[140] + sa[16] + sa[32], sa[1] + sa[109] + sa[125] + sa[141] + sa[17] + sa[33]], [sa[110] + sa[126] + sa[142] + sa[18] + sa[2] + sa[34], sa[111] + sa[127] + sa[143] + sa[19] + sa[3] + sa[35]]]) def test_array_permutedims(): sa = symbols('a0:144') for ArrayType in [ImmutableDenseNDimArray, ImmutableSparseNDimArray]: m1 = ArrayType(sa[:6], (2, 3)) assert permutedims(m1, (1, 0)) == transpose(m1) assert m1.tomatrix().T == permutedims(m1, (1, 0)).tomatrix() assert m1.tomatrix().T == transpose(m1).tomatrix() assert m1.tomatrix().C == conjugate(m1).tomatrix() assert m1.tomatrix().H == adjoint(m1).tomatrix() assert m1.tomatrix().T == m1.transpose().tomatrix() assert m1.tomatrix().C == m1.conjugate().tomatrix() assert m1.tomatrix().H == m1.adjoint().tomatrix() raises(ValueError, lambda: permutedims(m1, (0,))) raises(ValueError, lambda: permutedims(m1, (0, 0))) raises(ValueError, lambda: permutedims(m1, (1, 2, 0))) # Some tests with random arrays: dims = 6 shape = [random.randint(1,5) for i in range(dims)] elems = [random.random() for i in range(tensorproduct(*shape))] ra = ArrayType(elems, shape) perm = list(range(dims)) # Randomize the permutation: random.shuffle(perm) # Test inverse permutation: assert permutedims(permutedims(ra, perm), _af_invert(perm)) == ra # Test that permuted shape corresponds to action by `Permutation`: assert permutedims(ra, perm).shape == tuple(Permutation(perm)(shape)) z = ArrayType.zeros(4,5,6,7) assert permutedims(z, (2, 3, 1, 0)).shape == (6, 7, 5, 4) assert permutedims(z, [2, 3, 1, 0]).shape == (6, 7, 5, 4) assert permutedims(z, Permutation([2, 3, 1, 0])).shape == (6, 7, 5, 4) po = ArrayType(sa, [2, 2, 3, 3, 2, 2]) raises(ValueError, lambda: permutedims(po, (1, 1))) raises(ValueError, lambda: po.transpose()) raises(ValueError, lambda: po.adjoint()) assert permutedims(po, reversed(range(po.rank()))) == ArrayType( [[[[[[sa[0], sa[72]], [sa[36], sa[108]]], [[sa[12], sa[84]], [sa[48], sa[120]]], [[sa[24], sa[96]], [sa[60], sa[132]]]], [[[sa[4], sa[76]], [sa[40], sa[112]]], [[sa[16], sa[88]], [sa[52], sa[124]]], [[sa[28], sa[100]], [sa[64], sa[136]]]], [[[sa[8], sa[80]], [sa[44], sa[116]]], [[sa[20], sa[92]], [sa[56], sa[128]]], [[sa[32], sa[104]], [sa[68], sa[140]]]]], [[[[sa[2], sa[74]], [sa[38], sa[110]]], [[sa[14], sa[86]], [sa[50], sa[122]]], [[sa[26], sa[98]], [sa[62], sa[134]]]], [[[sa[6], sa[78]], [sa[42], sa[114]]], [[sa[18], sa[90]], [sa[54], sa[126]]], [[sa[30], sa[102]], [sa[66], sa[138]]]], [[[sa[10], sa[82]], [sa[46], sa[118]]], [[sa[22], sa[94]], [sa[58], sa[130]]], [[sa[34], sa[106]], [sa[70], sa[142]]]]]], [[[[[sa[1], sa[73]], [sa[37], sa[109]]], [[sa[13], sa[85]], [sa[49], sa[121]]], [[sa[25], sa[97]], [sa[61], sa[133]]]], [[[sa[5], sa[77]], [sa[41], sa[113]]], [[sa[17], sa[89]], [sa[53], sa[125]]], [[sa[29], sa[101]], [sa[65], sa[137]]]], [[[sa[9], sa[81]], [sa[45], sa[117]]], [[sa[21], sa[93]], [sa[57], sa[129]]], [[sa[33], sa[105]], [sa[69], sa[141]]]]], [[[[sa[3], sa[75]], [sa[39], sa[111]]], [[sa[15], sa[87]], [sa[51], sa[123]]], [[sa[27], sa[99]], [sa[63], sa[135]]]], [[[sa[7], sa[79]], [sa[43], sa[115]]], [[sa[19], sa[91]], [sa[55], sa[127]]], [[sa[31], sa[103]], [sa[67], sa[139]]]], [[[sa[11], sa[83]], [sa[47], sa[119]]], [[sa[23], sa[95]], [sa[59], sa[131]]], [[sa[35], sa[107]], [sa[71], sa[143]]]]]]]) assert permutedims(po, (1, 0, 2, 3, 4, 5)) == ArrayType( [[[[[[sa[0], sa[1]], [sa[2], sa[3]]], [[sa[4], sa[5]], [sa[6], sa[7]]], [[sa[8], sa[9]], [sa[10], sa[11]]]], [[[sa[12], sa[13]], [sa[14], sa[15]]], [[sa[16], sa[17]], [sa[18], sa[19]]], [[sa[20], sa[21]], [sa[22], sa[23]]]], [[[sa[24], sa[25]], [sa[26], sa[27]]], [[sa[28], sa[29]], [sa[30], sa[31]]], [[sa[32], sa[33]], [sa[34], sa[35]]]]], [[[[sa[72], sa[73]], [sa[74], sa[75]]], [[sa[76], sa[77]], [sa[78], sa[79]]], [[sa[80], sa[81]], [sa[82], sa[83]]]], [[[sa[84], sa[85]], [sa[86], sa[87]]], [[sa[88], sa[89]], [sa[90], sa[91]]], [[sa[92], sa[93]], [sa[94], sa[95]]]], [[[sa[96], sa[97]], [sa[98], sa[99]]], [[sa[100], sa[101]], [sa[102], sa[103]]], [[sa[104], sa[105]], [sa[106], sa[107]]]]]], [[[[[sa[36], sa[37]], [sa[38], sa[39]]], [[sa[40], sa[41]], [sa[42], sa[43]]], [[sa[44], sa[45]], [sa[46], sa[47]]]], [[[sa[48], sa[49]], [sa[50], sa[51]]], [[sa[52], sa[53]], [sa[54], sa[55]]], [[sa[56], sa[57]], [sa[58], sa[59]]]], [[[sa[60], sa[61]], [sa[62], sa[63]]], [[sa[64], sa[65]], [sa[66], sa[67]]], [[sa[68], sa[69]], [sa[70], sa[71]]]]], [ [[[sa[108], sa[109]], [sa[110], sa[111]]], [[sa[112], sa[113]], [sa[114], sa[115]]], [[sa[116], sa[117]], [sa[118], sa[119]]]], [[[sa[120], sa[121]], [sa[122], sa[123]]], [[sa[124], sa[125]], [sa[126], sa[127]]], [[sa[128], sa[129]], [sa[130], sa[131]]]], [[[sa[132], sa[133]], [sa[134], sa[135]]], [[sa[136], sa[137]], [sa[138], sa[139]]], [[sa[140], sa[141]], [sa[142], sa[143]]]]]]]) assert permutedims(po, (0, 2, 1, 4, 3, 5)) == ArrayType( [[[[[[sa[0], sa[1]], [sa[4], sa[5]], [sa[8], sa[9]]], [[sa[2], sa[3]], [sa[6], sa[7]], [sa[10], sa[11]]]], [[[sa[36], sa[37]], [sa[40], sa[41]], [sa[44], sa[45]]], [[sa[38], sa[39]], [sa[42], sa[43]], [sa[46], sa[47]]]]], [[[[sa[12], sa[13]], [sa[16], sa[17]], [sa[20], sa[21]]], [[sa[14], sa[15]], [sa[18], sa[19]], [sa[22], sa[23]]]], [[[sa[48], sa[49]], [sa[52], sa[53]], [sa[56], sa[57]]], [[sa[50], sa[51]], [sa[54], sa[55]], [sa[58], sa[59]]]]], [[[[sa[24], sa[25]], [sa[28], sa[29]], [sa[32], sa[33]]], [[sa[26], sa[27]], [sa[30], sa[31]], [sa[34], sa[35]]]], [[[sa[60], sa[61]], [sa[64], sa[65]], [sa[68], sa[69]]], [[sa[62], sa[63]], [sa[66], sa[67]], [sa[70], sa[71]]]]]], [[[[[sa[72], sa[73]], [sa[76], sa[77]], [sa[80], sa[81]]], [[sa[74], sa[75]], [sa[78], sa[79]], [sa[82], sa[83]]]], [[[sa[108], sa[109]], [sa[112], sa[113]], [sa[116], sa[117]]], [[sa[110], sa[111]], [sa[114], sa[115]], [sa[118], sa[119]]]]], [[[[sa[84], sa[85]], [sa[88], sa[89]], [sa[92], sa[93]]], [[sa[86], sa[87]], [sa[90], sa[91]], [sa[94], sa[95]]]], [[[sa[120], sa[121]], [sa[124], sa[125]], [sa[128], sa[129]]], [[sa[122], sa[123]], [sa[126], sa[127]], [sa[130], sa[131]]]]], [[[[sa[96], sa[97]], [sa[100], sa[101]], [sa[104], sa[105]]], [[sa[98], sa[99]], [sa[102], sa[103]], [sa[106], sa[107]]]], [[[sa[132], sa[133]], [sa[136], sa[137]], [sa[140], sa[141]]], [[sa[134], sa[135]], [sa[138], sa[139]], [sa[142], sa[143]]]]]]]) po2 = po.reshape(4, 9, 2, 2) assert po2 == ArrayType([[[[sa[0], sa[1]], [sa[2], sa[3]]], [[sa[4], sa[5]], [sa[6], sa[7]]], [[sa[8], sa[9]], [sa[10], sa[11]]], [[sa[12], sa[13]], [sa[14], sa[15]]], [[sa[16], sa[17]], [sa[18], sa[19]]], [[sa[20], sa[21]], [sa[22], sa[23]]], [[sa[24], sa[25]], [sa[26], sa[27]]], [[sa[28], sa[29]], [sa[30], sa[31]]], [[sa[32], sa[33]], [sa[34], sa[35]]]], [[[sa[36], sa[37]], [sa[38], sa[39]]], [[sa[40], sa[41]], [sa[42], sa[43]]], [[sa[44], sa[45]], [sa[46], sa[47]]], [[sa[48], sa[49]], [sa[50], sa[51]]], [[sa[52], sa[53]], [sa[54], sa[55]]], [[sa[56], sa[57]], [sa[58], sa[59]]], [[sa[60], sa[61]], [sa[62], sa[63]]], [[sa[64], sa[65]], [sa[66], sa[67]]], [[sa[68], sa[69]], [sa[70], sa[71]]]], [[[sa[72], sa[73]], [sa[74], sa[75]]], [[sa[76], sa[77]], [sa[78], sa[79]]], [[sa[80], sa[81]], [sa[82], sa[83]]], [[sa[84], sa[85]], [sa[86], sa[87]]], [[sa[88], sa[89]], [sa[90], sa[91]]], [[sa[92], sa[93]], [sa[94], sa[95]]], [[sa[96], sa[97]], [sa[98], sa[99]]], [[sa[100], sa[101]], [sa[102], sa[103]]], [[sa[104], sa[105]], [sa[106], sa[107]]]], [[[sa[108], sa[109]], [sa[110], sa[111]]], [[sa[112], sa[113]], [sa[114], sa[115]]], [[sa[116], sa[117]], [sa[118], sa[119]]], [[sa[120], sa[121]], [sa[122], sa[123]]], [[sa[124], sa[125]], [sa[126], sa[127]]], [[sa[128], sa[129]], [sa[130], sa[131]]], [[sa[132], sa[133]], [sa[134], sa[135]]], [[sa[136], sa[137]], [sa[138], sa[139]]], [[sa[140], sa[141]], [sa[142], sa[143]]]]]) assert permutedims(po2, (3, 2, 0, 1)) == ArrayType([[[[sa[0], sa[4], sa[8], sa[12], sa[16], sa[20], sa[24], sa[28], sa[32]], [sa[36], sa[40], sa[44], sa[48], sa[52], sa[56], sa[60], sa[64], sa[68]], [sa[72], sa[76], sa[80], sa[84], sa[88], sa[92], sa[96], sa[100], sa[104]], [sa[108], sa[112], sa[116], sa[120], sa[124], sa[128], sa[132], sa[136], sa[140]]], [[sa[2], sa[6], sa[10], sa[14], sa[18], sa[22], sa[26], sa[30], sa[34]], [sa[38], sa[42], sa[46], sa[50], sa[54], sa[58], sa[62], sa[66], sa[70]], [sa[74], sa[78], sa[82], sa[86], sa[90], sa[94], sa[98], sa[102], sa[106]], [sa[110], sa[114], sa[118], sa[122], sa[126], sa[130], sa[134], sa[138], sa[142]]]], [[[sa[1], sa[5], sa[9], sa[13], sa[17], sa[21], sa[25], sa[29], sa[33]], [sa[37], sa[41], sa[45], sa[49], sa[53], sa[57], sa[61], sa[65], sa[69]], [sa[73], sa[77], sa[81], sa[85], sa[89], sa[93], sa[97], sa[101], sa[105]], [sa[109], sa[113], sa[117], sa[121], sa[125], sa[129], sa[133], sa[137], sa[141]]], [[sa[3], sa[7], sa[11], sa[15], sa[19], sa[23], sa[27], sa[31], sa[35]], [sa[39], sa[43], sa[47], sa[51], sa[55], sa[59], sa[63], sa[67], sa[71]], [sa[75], sa[79], sa[83], sa[87], sa[91], sa[95], sa[99], sa[103], sa[107]], [sa[111], sa[115], sa[119], sa[123], sa[127], sa[131], sa[135], sa[139], sa[143]]]]]) # test for large scale sparse array for SparseArrayType in [ImmutableSparseNDimArray, MutableSparseNDimArray]: A = SparseArrayType({1:1, 10000:2}, (10000, 20000, 10000)) assert permutedims(A, (0, 1, 2)) == A assert permutedims(A, (1, 0, 2)) == SparseArrayType({1: 1, 100000000: 2}, (20000, 10000, 10000)) B = SparseArrayType({1:1, 20000:2}, (10000, 20000)) assert B.transpose() == SparseArrayType({10000: 1, 1: 2}, (20000, 10000)) def test_flatten(): from sympy.matrices.dense import Matrix for ArrayType in [ImmutableDenseNDimArray, ImmutableSparseNDimArray, Matrix]: A = ArrayType(range(24)).reshape(4, 6) assert [i for i in Flatten(A)] == [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23] for i, v in enumerate(Flatten(A)): assert i == v def test_tensordiagonal(): from sympy.matrices.dense import eye expr = Array(range(9)).reshape(3, 3) raises(ValueError, lambda: tensordiagonal(expr, [0], [1])) raises(ValueError, lambda: tensordiagonal(expr, [0, 0])) assert tensordiagonal(eye(3), [0, 1]) == Array([1, 1, 1]) assert tensordiagonal(expr, [0, 1]) == Array([0, 4, 8]) x, y, z = symbols("x y z") expr2 = tensorproduct([x, y, z], expr) assert tensordiagonal(expr2, [1, 2]) == Array([[0, 4*x, 8*x], [0, 4*y, 8*y], [0, 4*z, 8*z]]) assert tensordiagonal(expr2, [0, 1]) == Array([[0, 3*y, 6*z], [x, 4*y, 7*z], [2*x, 5*y, 8*z]]) assert tensordiagonal(expr2, [0, 1, 2]) == Array([0, 4*y, 8*z]) # assert tensordiagonal(expr2, [0]) == permutedims(expr2, [1, 2, 0]) # assert tensordiagonal(expr2, [1]) == permutedims(expr2, [0, 2, 1]) # assert tensordiagonal(expr2, [2]) == expr2 # assert tensordiagonal(expr2, [1], [2]) == expr2 # assert tensordiagonal(expr2, [0], [1]) == permutedims(expr2, [2, 0, 1]) a, b, c, X, Y, Z = symbols("a b c X Y Z") expr3 = tensorproduct([x, y, z], [1, 2, 3], [a, b, c], [X, Y, Z]) assert tensordiagonal(expr3, [0, 1, 2, 3]) == Array([x*a*X, 2*y*b*Y, 3*z*c*Z]) assert tensordiagonal(expr3, [0, 1], [2, 3]) == tensorproduct([x, 2*y, 3*z], [a*X, b*Y, c*Z]) # assert tensordiagonal(expr3, [0], [1, 2], [3]) == tensorproduct([x, y, z], [a, 2*b, 3*c], [X, Y, Z]) assert tensordiagonal(tensordiagonal(expr3, [2, 3]), [0, 1]) == tensorproduct([a*X, b*Y, c*Z], [x, 2*y, 3*z]) raises(ValueError, lambda: tensordiagonal([[1, 2, 3], [4, 5, 6]], [0, 1])) raises(ValueError, lambda: tensordiagonal(expr3.reshape(3, 3, 9), [1, 2]))