N维数组¶

SymPy 的 N 维数组模块。

提供了四个类来处理 N 维数组，由密集/稀疏（即是否在内存中存储所有元素或仅存储非零元素）和可变/不可变（不可变类是 SymPy 对象，但创建后不能更改）的组合给出。

示例

以下示例展示了 Array 的使用。这是 ImmutableDenseNDimArray 的缩写，即一个不可变且密集的 N 维数组，其他类与此类似。对于可变类，也可以在对象构造后更改元素值。

数组构造可以检测嵌套列表和元组的形状：

>>> from sympy import Array
>>> a1 = Array([[1, 2], [3, 4], [5, 6]])
>>> a1
[[1, 2], [3, 4], [5, 6]]
>>> a1.shape
(3, 2)
>>> a1.rank()
2
>>> from sympy.abc import x, y, z
>>> a2 = Array([[[x, y], [z, x*z]], [[1, x*y], [1/x, x/y]]])
>>> a2
[[[x, y], [z, x*z]], [[1, x*y], [1/x, x/y]]]
>>> a2.shape
(2, 2, 2)
>>> a2.rank()
3

否则可以传递一个一维数组，后面跟着一个形状元组：

>>> m1 = Array(range(12), (3, 4))
>>> m1
[[0, 1, 2, 3], [4, 5, 6, 7], [8, 9, 10, 11]]
>>> m2 = Array(range(12), (3, 2, 2))
>>> m2
[[[0, 1], [2, 3]], [[4, 5], [6, 7]], [[8, 9], [10, 11]]]
>>> m2[1,1,1]
7
>>> m2.reshape(4, 3)
[[0, 1, 2], [3, 4, 5], [6, 7, 8], [9, 10, 11]]

切片支持：

>>> m2[:, 1, 1]
[3, 7, 11]

逐元素导数：

>>> from sympy.abc import x, y, z
>>> m3 = Array([x**3, x*y, z])
>>> m3.diff(x)
[3*x**2, y, 0]
>>> m3.diff(z)
[0, 0, 1]

与其他 SymPy 表达式的乘法是按元素应用的：

>>> (1+x)*m3
[x**3*(x + 1), x*y*(x + 1), z*(x + 1)]

要对 N 维数组的每个元素应用函数，请使用 applyfunc：

>>> m3.applyfunc(lambda x: x/2)
[x**3/2, x*y/2, z/2]

N 维数组可以通过 tolist() 方法转换为嵌套列表：

>>> m2.tolist()
[[[0, 1], [2, 3]], [[4, 5], [6, 7]], [[8, 9], [10, 11]]]
>>> isinstance(m2.tolist(), list)
True

如果秩为2，可以使用 tomatrix() 将其转换为矩阵：

>>> m1.tomatrix()
Matrix([
[0, 1,  2,  3],
[4, 5,  6,  7],
[8, 9, 10, 11]])

产品和收缩¶

数组 \(A_{i_1,\ldots,i_n}\) 和 \(B_{j_1,\ldots,j_m}\) 的张量积生成组合数组 \(P = A \otimes B\)，定义为

\(P_{i_1,\ldots,i_n,j_1,\ldots,j_m} := A_{i_1,\ldots,i_n}\cdot B_{j_1,\ldots,j_m}.\)

它可以通过 tensorproduct(...) 获得：

>>> from sympy import Array, tensorproduct
>>> from sympy.abc import x,y,z,t
>>> A = Array([x, y, z, t])
>>> B = Array([1, 2, 3, 4])
>>> tensorproduct(A, B)
[[x, 2*x, 3*x, 4*x], [y, 2*y, 3*y, 4*y], [z, 2*z, 3*z, 4*z], [t, 2*t, 3*t, 4*t]]

如果你不想立即计算张量积，你可以使用 ArrayTensorProduct，它会创建一个未计算的张量积表达式：

>>> from sympy.tensor.array.expressions import ArrayTensorProduct
>>> ArrayTensorProduct(A, B)
ArrayTensorProduct([x, y, z, t], [1, 2, 3, 4])

在 ArrayTensorProduct 上调用 .as_explicit() 等同于直接调用 tensorproduct(...)：

>>> ArrayTensorProduct(A, B).as_explicit()
[[x, 2*x, 3*x, 4*x], [y, 2*y, 3*y, 4*y], [z, 2*z, 3*z, 4*z], [t, 2*t, 3*t, 4*t]]

秩为1的数组与矩阵的张量积生成一个秩为3的数组：

>>> from sympy import eye
>>> p1 = tensorproduct(A, eye(4))
>>> p1
[[[x, 0, 0, 0], [0, x, 0, 0], [0, 0, x, 0], [0, 0, 0, x]], [[y, 0, 0, 0], [0, y, 0, 0], [0, 0, y, 0], [0, 0, 0, y]], [[z, 0, 0, 0], [0, z, 0, 0], [0, 0, z, 0], [0, 0, 0, z]], [[t, 0, 0, 0], [0, t, 0, 0], [0, 0, t, 0], [0, 0, 0, t]]]

现在，要回到 \(A_0 \otimes \mathbf{1}\)，可以通过切片访问 \(p_{0,m,n}\)：

>>> p1[0,:,:]
[[x, 0, 0, 0], [0, x, 0, 0], [0, 0, x, 0], [0, 0, 0, x]]

张量收缩对指定的轴进行求和，例如收缩位置 \(a\) 和 \(b\) 意味着

\(A_{i_1,\ldots,i_a,\ldots,i_b,\ldots,i_n} \implies \sum_k A_{i_1,\ldots,k,\ldots,k,\ldots,i_n}\)

记住Python的索引是从零开始的，因此要合并第a轴和第b轴，需要指定`a-1`和`b-1`

>>> from sympy import tensorcontraction
>>> C = Array([[x, y], [z, t]])

矩阵的迹等价于一个秩为2的数组的收缩：

\(A_{m,n} \implies \sum_k A_{k,k}\)

>>> tensorcontraction(C, (0, 1))
t + x

要创建一个表示张量收缩的表达式，该表达式不会立即求值，请使用 ArrayContraction，如果其后跟随 .as_explicit()，则等效于 tensorcontraction(...)：

>>> from sympy.tensor.array.expressions import ArrayContraction
>>> ArrayContraction(C, (0, 1))
ArrayContraction([[x, y], [z, t]], (0, 1))
>>> ArrayContraction(C, (0, 1)).as_explicit()
t + x

矩阵乘积等价于两个秩为2的数组的张量积，随后对第2和第3轴（在Python索引中为轴号1, 2）进行收缩。

\(A_{m,n}\cdot B_{i,j} \implies \sum_k A_{m, k}\cdot B_{k, j}\)

>>> D = Array([[2, 1], [0, -1]])
>>> tensorcontraction(tensorproduct(C, D), (1, 2))
[[2*x, x - y], [2*z, -t + z]]

可以验证矩阵乘积是等价的：

>>> from sympy import Matrix
>>> Matrix([[x, y], [z, t]])*Matrix([[2, 1], [0, -1]])
Matrix([
[2*x,  x - y],
[2*z, -t + z]])

或者等效地

>>> C.tomatrix()*D.tomatrix()
Matrix([
[2*x,  x - y],
[2*z, -t + z]])

对角运算符¶

tensordiagonal 函数的作用类似于 tensorcontraction，但连接的索引不会被求和，例如对位置 \(a\) 和 \(b\) 进行对角化意味着

\(A_{i_1,\ldots,i_a,\ldots,i_b,\ldots,i_n} \implies A_{i_1,\ldots,k,\ldots,k,\ldots,i_n} \implies \tilde{A}_{i_1,\ldots,i_{a-1},i_{a+1},\ldots,i_{b-1},i_{b+1},\ldots,i_n,k}\)

其中 \(\tilde{A}\) 是与 \(A\) 的对角线在位置 \(a\) 和 \(b\) 上移动到最后一个索引槽的数组等效项。

比较收缩和对角运算符之间的差异：

>>> from sympy import tensordiagonal
>>> from sympy.abc import a, b, c, d
>>> m = Matrix([[a, b], [c, d]])
>>> tensorcontraction(m, [0, 1])
a + d
>>> tensordiagonal(m, [0, 1])
[a, d]

简而言之，使用 tensordiagonal 不会发生求和。

数组求导¶

通常的导数操作可以扩展以支持对数组的导数，前提是该数组中的所有元素都是符号或适合进行导数的表达式。

数组的导数定义如下：给定数组 \(A_{i_1, \ldots, i_N}\) 和数组 \(X_{j_1, \ldots, j_M}\)，数组的导数将返回一个新数组 \(B\)，其定义为

\(B_{j_1,\ldots,j_M,i_1,\ldots,i_N} := \frac{\partial A_{i_1,\ldots,i_N}}{\partial X_{j_1,\ldots,j_M}}\)

函数 derive_by_array 执行这样的操作：

>>> from sympy import derive_by_array
>>> from sympy.abc import x, y, z, t
>>> from sympy import sin, exp

对于标量，它的行为与普通导数完全相同：

>>> derive_by_array(sin(x*y), x)
y*cos(x*y)

由数组基导出的标量：

>>> derive_by_array(sin(x*y), [x, y, z])
[y*cos(x*y), x*cos(x*y), 0]

数组按数组基推导：\(B^{nm} := \frac{\partial A^m}{\partial x^n}\)

>>> basis = [x, y, z]
>>> ax = derive_by_array([exp(x), sin(y*z), t], basis)
>>> ax
[[exp(x), 0, 0], [0, z*cos(y*z), 0], [0, y*cos(y*z), 0]]

结果数组的收缩：\(\sum_m \frac{\partial A^m}{\partial x^m}\)

>>> tensorcontraction(ax, (0, 1))
z*cos(y*z) + exp(x)

类¶

class sympy.tensor.array.ImmutableDenseNDimArray(

iterable,

shape=None,

**kwargs,

)[源代码][源代码]¶

属性:

args: 返回 ‘self’ 的参数元组。
assumptions0: 返回对象 \(type\) 假设。
canonical_variables: 返回一个字典，将 self.bound_symbols 中定义的任何变量映射到与表达式中任何自由符号不冲突的符号。
expr_free_symbols
free_symbols: 从自身的原子中返回那些自由符号。
func: 表达式中的顶级函数。
is_algebraic
is_antihermitian
is_commutative
is_comparable: 如果 self 可以计算为一个具有精度的实数（或已经是一个实数），则返回 True，否则返回 False。
is_complex
is_composite
is_even
is_extended_negative
is_extended_nonnegative
is_extended_nonpositive
is_extended_nonzero
is_extended_positive
is_extended_real
is_finite
is_hermitian
is_imaginary
is_infinite
is_integer
is_irrational
is_negative
is_noninteger
is_nonnegative
is_nonpositive
is_nonzero
is_odd
is_polar
is_positive
is_prime
is_rational
is_real
is_transcendental
is_zero
种类
shape: 返回数组形状（维度）。

方法

`applyfunc`(f)	对 N 维数组的每个元素应用一个函数。
`as_content_primitive`([radical, clear])	一个存根，允许在计算表达式的内容和基本组件时跳过基本参数（如元组）。
`as_dummy`()	返回表达式，其中任何具有结构绑定符号的对象都被替换为在其出现的对象中唯一的规范符号，并且仅对交换性具有默认假设为True。
`atoms`(*types)	返回构成当前对象的原子。
`class_key`()	类的好顺序。
`compare`(other)	如果对象在规范意义上小于、等于或大于其他对象，则返回 -1、0、1。
`count`(query)	计算匹配的子表达式的数量。
`count_ops`([visual])	用于返回操作计数的 count_ops 的包装器。
`diff`(args, *kwargs)	计算数组中每个元素的导数。
`doit`(**hints)	评估默认情况下不评估的对象，如极限、积分、求和和乘积。
`dummy_eq`(other[, symbol])	比较两个表达式并处理哑符号。
`find`(query[, group])	查找所有匹配查询的子表达式。
`fromiter`(args, **assumptions)	从可迭代对象创建一个新对象。
`has`(*patterns)	测试是否有任何子表达式匹配任何模式。
`has_free`(*patterns)	如果 self 包含对象 `x` 作为自由表达式，则返回 True，否则返回 False。
`has_xfree`(s)	如果 self 有 s 中的任何一个模式作为自由参数，则返回 True，否则返回 False。
`is_same`(b[, approx])	如果 a 和 b 结构相同则返回 True，否则返回 False。
`match`(pattern[, old])	模式匹配。
`matches`(expr[, repl_dict, old])	用于 match() 的辅助方法，用于在 self 中的通配符符号与 expr 中的表达式之间寻找匹配。
`rank`()	返回数组的秩。
`rcall`(*args)	通过表达式树递归应用于参数。
`refine`([assumption])	请参阅 sympy.assumptions 中的 refine 函数。
`replace`(query, value[, map, simultaneous, exact])	将 `self` 中匹配的子表达式替换为 `value`。
`reshape`(*newshape)	返回具有新形状的 MutableDenseNDimArray 实例。
`rewrite`(*args[, deep])	使用定义的规则重写 self。
`simplify`(**kwargs)	请参阅 sympy.simplify 中的 simplify 函数。
`sort_key`([order])	返回一个排序键。
`subs`(args, *kwargs)	在简化参数后，在表达式中用新内容替换旧内容。
`tolist`()	将 MutableDenseNDimArray 转换为一维列表
`tomatrix`()	将 MutableDenseNDimArray 转换为 Matrix。
`xreplace`(rule)	替换表达式中对象的出现。

伴随
as_immutable
as_mutable
共轭
复制
could_extract_minus_sign
is_hypergeometric
转置
零

class sympy.tensor.array.ImmutableSparseNDimArray(

iterable=None,

shape=None,

**kwargs,

)[源代码][源代码]¶

属性:

args: 返回 ‘self’ 的参数元组。
assumptions0: 返回对象 \(type\) 假设。
canonical_variables: 返回一个字典，将 self.bound_symbols 中定义的任何变量映射到与表达式中任何自由符号不冲突的符号。
expr_free_symbols
free_symbols: 从自身的原子中返回那些自由符号。
func: 表达式中的顶级函数。
is_algebraic
is_antihermitian
is_commutative
is_comparable: 如果 self 可以计算为一个具有精度的实数（或已经是一个实数），则返回 True，否则返回 False。
is_complex
is_composite
is_even
is_extended_negative
is_extended_nonnegative
is_extended_nonpositive
is_extended_nonzero
is_extended_positive
is_extended_real
is_finite
is_hermitian
is_imaginary
is_infinite
is_integer
is_irrational
is_negative
is_noninteger
is_nonnegative
is_nonpositive
is_nonzero
is_odd
is_polar
is_positive
is_prime
is_rational
is_real
is_transcendental
is_zero
shape: 返回数组形状（维度）。

方法

`applyfunc`(f)	对 N 维数组的每个元素应用一个函数。
`as_content_primitive`([radical, clear])	一个存根，允许在计算表达式的内容和基本组件时跳过基本参数（如元组）。
`as_dummy`()	返回表达式，其中任何具有结构绑定符号的对象都被替换为在其出现的对象中唯一的规范符号，并且仅对交换性具有默认假设为True。
`atoms`(*types)	返回构成当前对象的原子。
`class_key`()	类的好顺序。
`compare`(other)	如果对象在规范意义上小于、等于或大于其他对象，则返回 -1、0、1。
`count`(query)	计算匹配的子表达式的数量。
`count_ops`([visual])	用于返回操作计数的 count_ops 的包装器。
`diff`(args, *kwargs)	计算数组中每个元素的导数。
`doit`(**hints)	评估默认情况下不评估的对象，如极限、积分、求和和乘积。
`dummy_eq`(other[, symbol])	比较两个表达式并处理哑符号。
`find`(query[, group])	查找所有匹配查询的子表达式。
`fromiter`(args, **assumptions)	从可迭代对象创建一个新对象。
`has`(*patterns)	测试是否有任何子表达式匹配任何模式。
`has_free`(*patterns)	如果 self 包含对象 `x` 作为自由表达式，则返回 True，否则返回 False。
`has_xfree`(s)	如果 self 有 s 中的任何一个模式作为自由参数，则返回 True，否则返回 False。
`is_same`(b[, approx])	如果 a 和 b 结构相同则返回 True，否则返回 False。
`match`(pattern[, old])	模式匹配。
`matches`(expr[, repl_dict, old])	用于 match() 的辅助方法，用于在 self 中的通配符符号与 expr 中的表达式之间寻找匹配。
`rank`()	返回数组的秩。
`rcall`(*args)	通过表达式树递归应用于参数。
`refine`([assumption])	请参阅 sympy.assumptions 中的 refine 函数。
`replace`(query, value[, map, simultaneous, exact])	将 `self` 中匹配的子表达式替换为 `value`。
`rewrite`(*args[, deep])	使用定义的规则重写 self。
`simplify`(**kwargs)	请参阅 sympy.simplify 中的 simplify 函数。
`sort_key`([order])	返回一个排序键。
`subs`(args, *kwargs)	在简化参数后，在表达式中用新内容替换旧内容。
`tolist`()	将 MutableDenseNDimArray 转换为一维列表
`tomatrix`()	将 MutableDenseNDimArray 转换为 Matrix。
`xreplace`(rule)	替换表达式中对象的出现。
`zeros`(*shape)	返回一个全零的稀疏N维数组。

伴随
as_immutable
as_mutable
共轭
复制
could_extract_minus_sign
is_hypergeometric
重塑
转置

class sympy.tensor.array.MutableDenseNDimArray(

iterable=None,

shape=None,

**kwargs,

)[源代码][源代码]¶

属性:

自由符号
种类
shape: 返回数组形状（维度）。

方法

`applyfunc`(f)	对 N 维数组的每个元素应用一个函数。
`diff`(args, *kwargs)	计算数组中每个元素的导数。
`rank`()	返回数组的秩。
`reshape`(*newshape)	返回具有新形状的 MutableDenseNDimArray 实例。
`tolist`()	将 MutableDenseNDimArray 转换为一维列表
`tomatrix`()	将 MutableDenseNDimArray 转换为 Matrix。

伴随
as_immutable
as_mutable
共轭
转置
零

class sympy.tensor.array.MutableSparseNDimArray(

iterable=None,

shape=None,

**kwargs,

)[源代码][源代码]¶

属性:

自由符号
shape: 返回数组形状（维度）。

方法

`applyfunc`(f)	对 N 维数组的每个元素应用一个函数。
`diff`(args, *kwargs)	计算数组中每个元素的导数。
`rank`()	返回数组的秩。
`tolist`()	将 MutableDenseNDimArray 转换为一维列表
`tomatrix`()	将 MutableDenseNDimArray 转换为 Matrix。
`zeros`(*shape)	返回一个全零的稀疏N维数组。

伴随
as_immutable
as_mutable
共轭
重塑
转置

函数¶

sympy.tensor.array.derive_by_array(expr, dx)[源代码][源代码]¶

数组求导。支持数组和标量。

数组表达式的等价运算符是 array_derive。

示例

>>> from sympy import derive_by_array
>>> from sympy.abc import x, y, z, t
>>> from sympy import cos
>>> derive_by_array(cos(x*t), x)
-t*sin(t*x)
>>> derive_by_array(cos(x*t), [x, y, z, t])
[-t*sin(t*x), 0, 0, -x*sin(t*x)]
>>> derive_by_array([x, y**2*z], [[x, y], [z, t]])
[[[1, 0], [0, 2*y*z]], [[0, y**2], [0, 0]]]

sympy.tensor.array.permutedims( expr, perm=None, index_order_old=None, index_order_new=None, )[源代码][源代码]¶

置换数组的索引。

参数指定索引的排列。

数组表达式的等价运算符是 PermuteDims ，它可以用来保持表达式未求值。

参见

sympy.tensor.array.expressions.array_expressions.PermuteDims

示例

>>> from sympy.abc import x, y, z, t
>>> from sympy import sin
>>> from sympy import Array, permutedims
>>> a = Array([[x, y, z], [t, sin(x), 0]])
>>> a
[[x, y, z], [t, sin(x), 0]]
>>> permutedims(a, (1, 0))
[[x, t], [y, sin(x)], [z, 0]]

如果数组是二阶的，可以使用 transpose：

>>> from sympy import transpose
>>> transpose(a)
[[x, t], [y, sin(x)], [z, 0]]

高维度的例子：

>>> b = Array([[[1, 2], [3, 4]], [[5, 6], [7, 8]]])
>>> permutedims(b, (2, 1, 0))
[[[1, 5], [3, 7]], [[2, 6], [4, 8]]]
>>> permutedims(b, (1, 2, 0))
[[[1, 5], [2, 6]], [[3, 7], [4, 8]]]

指定与前几行相同排列的另一种方法涉及传递旧和新索引，可以作为列表或字符串传递：

>>> permutedims(b, index_order_old="cba", index_order_new="abc")
[[[1, 5], [3, 7]], [[2, 6], [4, 8]]]
>>> permutedims(b, index_order_old="cab", index_order_new="abc")
[[[1, 5], [2, 6]], [[3, 7], [4, 8]]]

Permutation 对象也是允许的：

>>> from sympy.combinatorics import Permutation
>>> permutedims(b, Permutation([1, 2, 0]))
[[[1, 5], [2, 6]], [[3, 7], [4, 8]]]

sympy.tensor.array.tensorcontraction(array, *contraction_axes)[源代码][源代码]¶

在指定轴上收缩类似数组的对象。

数组表达式的等价运算符是 ArrayContraction ，它可以用来保持表达式未求值。

参见

sympy.tensor.array.expressions.array_expressions.ArrayContraction

示例

>>> from sympy import Array, tensorcontraction
>>> from sympy import Matrix, eye
>>> tensorcontraction(eye(3), (0, 1))
3
>>> A = Array(range(18), (3, 2, 3))
>>> A
[[[0, 1, 2], [3, 4, 5]], [[6, 7, 8], [9, 10, 11]], [[12, 13, 14], [15, 16, 17]]]
>>> tensorcontraction(A, (0, 2))
[21, 30]

矩阵乘法可以通过 tensorcontraction 和 tensorproduct 的适当组合来模拟。

>>> from sympy import tensorproduct
>>> from sympy.abc import a,b,c,d,e,f,g,h
>>> m1 = Matrix([[a, b], [c, d]])
>>> m2 = Matrix([[e, f], [g, h]])
>>> p = tensorproduct(m1, m2)
>>> p
[[[[a*e, a*f], [a*g, a*h]], [[b*e, b*f], [b*g, b*h]]], [[[c*e, c*f], [c*g, c*h]], [[d*e, d*f], [d*g, d*h]]]]
>>> tensorcontraction(p, (1, 2))
[[a*e + b*g, a*f + b*h], [c*e + d*g, c*f + d*h]]
>>> m1*m2
Matrix([
[a*e + b*g, a*f + b*h],
[c*e + d*g, c*f + d*h]])

sympy.tensor.array.tensorproduct(*args)[源代码][源代码]¶

标量或类数组对象之间的张量积。

数组表达式的等价运算符是 ArrayTensorProduct，它可以用来保持表达式未求值。

参见

sympy.tensor.array.expressions.array_expressions.ArrayTensorProduct

示例

>>> from sympy.tensor.array import tensorproduct, Array
>>> from sympy.abc import x, y, z, t
>>> A = Array([[1, 2], [3, 4]])
>>> B = Array([x, y])
>>> tensorproduct(A, B)
[[[x, y], [2*x, 2*y]], [[3*x, 3*y], [4*x, 4*y]]]
>>> tensorproduct(A, x)
[[x, 2*x], [3*x, 4*x]]
>>> tensorproduct(A, B, B)
[[[[x**2, x*y], [x*y, y**2]], [[2*x**2, 2*x*y], [2*x*y, 2*y**2]]], [[[3*x**2, 3*x*y], [3*x*y, 3*y**2]], [[4*x**2, 4*x*y], [4*x*y, 4*y**2]]]]

将此函数应用于两个矩阵将产生一个秩为4的数组。

>>> from sympy import Matrix, eye
>>> m = Matrix([[x, y], [z, t]])
>>> p = tensorproduct(eye(3), m)
>>> p
[[[[x, y], [z, t]], [[0, 0], [0, 0]], [[0, 0], [0, 0]]], [[[0, 0], [0, 0]], [[x, y], [z, t]], [[0, 0], [0, 0]]], [[[0, 0], [0, 0]], [[0, 0], [0, 0]], [[x, y], [z, t]]]]

sympy.tensor.array.tensordiagonal(array, *diagonal_axes)[源代码][源代码]¶

在指定的轴上对类数组对象进行对角化。

这相当于通过克罗内克δ将轴连接起来乘以表达式。

对角索引被放在轴的末尾。

数组表达式的等价运算符是 ArrayDiagonal，它可以用来保持表达式未求值。

参见

sympy.tensor.array.expressions.array_expressions.ArrayDiagonal

示例

tensordiagonal 通过轴 0 和 1 作用于一个二维数组，等价于矩阵的对角线：

>>> from sympy import Array, tensordiagonal
>>> from sympy import Matrix, eye
>>> tensordiagonal(eye(3), (0, 1))
[1, 1, 1]

>>> from sympy.abc import a,b,c,d
>>> m1 = Matrix([[a, b], [c, d]])
>>> tensordiagonal(m1, [0, 1])
[a, d]

对于更高维的数组，对角化的维度会被移除并在最后作为一个单一维度追加：

>>> A = Array(range(18), (3, 2, 3))
>>> A
[[[0, 1, 2], [3, 4, 5]], [[6, 7, 8], [9, 10, 11]], [[12, 13, 14], [15, 16, 17]]]
>>> tensordiagonal(A, (0, 2))
[[0, 7, 14], [3, 10, 17]]
>>> from sympy import permutedims
>>> tensordiagonal(A, (0, 2)) == permutedims(Array([A[0, :, 0], A[1, :, 1], A[2, :, 2]]), [1, 0])
True