张量的推导与计算十分繁杂，因而使用Mathematica进行张量分析本身便是一件很自然的事，但很无奈，网上几乎没有相应的中文教程，与之相应的是，各种论坛上关于使用Mathematica进行张量分析的求助帖基本没有回应。因此，本帖系统地总结了一下Mathematica中各种张量函数的用法，并辅以部分例题，为大家展示如何使用Mathematica进行张量分析。限于篇幅，希望本帖能起到抛砖引玉的作用。
1. 张量的表示
在Mathematica中，使用多层列表表示张量。例如，二阶对称克罗内克符号为二阶张量，我们可以手动输入，或者借助内置函数产生。里奇-列维塔张量亦是如此。
{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}} // MatrixForm
Array[KroneckerDelta, {3, 3}] // MatrixForm
Table[Mod[(j - i) (k - j) (i - k), 3, -1], {i, 1, 3}, {j, 1, 3}, {k,
1, 3}] // MatrixForm
LeviCivitaTensor[3] // MatrixForm
实际上，Mathematica并不能直接区分协变张量和逆变张量，所以我们需要借助度量张量来实现。例如在极坐标下，我们已知某一矢量在自然协变基矢量下的逆变分量，求其在逆变基矢量下的协变分量。
Gx = CoordinateChartData["Polar", "Metric", {r, \[Theta]}];
Gn = CoordinateChartData["Polar", "InverseMetric", {r, \[Theta]}];
Gx // MatrixForm
Gn // MatrixForm
Pn = {1, r}
Px = Pn.Gn
Px.Gx
张量分析中，有很大一部分是对各种等式的证明，因此抽象表示张量很有必要，Mathematica也提供了这种表示方法，可以看到，这种表示方法可以指定任意维度，数据类型以及对称性等。关于张量等式的证明会在随后章节中给出。
$Assumptions =
A \[Element] Arrays[{4, 4, 4}, Reals, Symmetric[{1, 2, 3}]]
$Assumptions =
M \[Element] Matrices[{3, 3}, Reals, Antisymmetric[{1, 2}]]
$Assumptions = V \[Element] Vectors[3, Reals]
2.张量的代数运算
这里只给出张量分析中较为重要的几个代数运算，并积，缩并，内积，转置。
Clear["Global`*"];
(*并积*)
(*A\[CircleTimes]B=Subscript[A, ij]Subscript[B, mn]*)
A = {a[1], a[2], a[3]}; B = {b[1], b[2], b[3]};
TensorProduct[A, B] // MatrixForm
A\[TensorProduct]B // MatrixForm
Table[x[i, j], {i, 1, 2}, {j, 1, 2}]\[TensorProduct]Table[
y[i, j], {i, 1, 2}, {j, 1, 2}] // MatrixForm
张量的并积用TensorProduct函数来实现，或用"[TensorProduct]"符号来表示，并没有什么特殊的技巧。
Clear["Global`*"];
(*缩并*)
(*Subscript[A, ii]*)
A = Table[a[i, j], {i, 1, 3}, {j, 1, 3}];
A // MatrixForm
TensorContract[A, {{1, 2}}]
Tr[A]
张量的缩并用TensorContract函数来实现，可以具体指定缩并位置。对于矩阵来说，缩并后的结果便是矩阵的迹。
Clear["Global`*"];
(*点积*)
(*A.B=Subscript[A, ij]Subscript[B, jk]*)
A = Table[a[i, j], {i, 1, 2}, {j, 1, 2}];
B = Table[b[i, j], {i, 1, 2}, {j, 1, 2}];
Dot[A, B] // MatrixForm
A.B // MatrixForm
张量的内积用Dot函数来实现，或用"."符号表示，要注意内积是第一个张量最后一层和第二个张量第一层上的运算，因此一般是不可交换的。
Clear["Global`*"];
(*双点积*)
(*A:B=Subscript[A, ij]Subscript[B, ji]*)
A = Table[a[i, j], {i, 1, 2}, {j, 1, 2}];
B = Table[b[i, j], {i, 1, 2}, {j, 1, 2}];
TensorContract[(A\[TensorProduct]B), {{2, 3}, {1, 4}}]
张量的双点积可以使用并积函数与缩并函数一起实现，当然也可以是三点积，四点积等。例如这个例子中，双点积是定义为A矩阵和B矩阵在相邻部分的点积与远离部分的点积。
Clear["Global`*"];
(*转置*)
A = Table[i - j + k, {i, 1, 3}, {j, 1, 3}, {k, 1, 3}];
Transpose[A, {2, 1, 3}] // MatrixForm
Transpose[A] // MatrixForm
张量的转置使用Transpose函数来实现，默认参数下，在张量的前两层内转置，当然，设置参数后，Transpose函数是一种更广义的转置。
3.曲线坐标下的张量分析
本节的主要目的是推导不同坐标系下的张量方程，例如不同坐标系下弹性力学平衡方程，纳维-斯托克斯方程等。通常手动推导这类问题是相当麻烦的，比如教科书中附录会带有不同坐标系下的张量方程，但实际中，我们会遇到更多其他类似的问题却是无法查附录的。故使用Mathematica很有必要。下面先介绍梯度，散度的用法。
Clear["Global`*"];
A = Table[f[i][r, \[Theta]], {i, 1, 2}];
A // MatrixForm
Grad[A, {r, \[Theta]}, "Polar"]
Div[A, {r, \[Theta]}, "Polar"]
B = Table[f[i, j][x, y], {i, 1, 2}, {j, 1, 2}];
B // MatrixForm
Grad[B, {x, y}] // MatrixForm
Div[B, {x, y}] // MatrixForm
要注意两点，其一是这里的梯度散度均是在标准正交曲线坐标下的运算，即基矢量均为单位长度。其二是在Mathematica中，进行梯度散度计算时，增加的维度放在最外层，这意味着直接进行梯度散度运算时，是指右梯度和右散度，需要特别注意，当然，我们可以使用Transpose函数来将其化为左梯度和左散度，毕竟这才是我们最常用的。下面给出一个例子，各位读者可以好好体会一下。（C为二阶张量）
Clear["Global`*"];
c = {{Cxx[x, y], Cxy[x, y]}, {Cxy[x, y], Cyy[x, y]}};
u = {U[x, y], V[x, y]};
u.Transpose[Grad[c, {x, y}], {2, 3, 1}] // MatrixForm
(Grad[u, {x, y}]).c + c.Transpose[(Grad[u, {x, y}])] // MatrixForm
Div[Transpose@c, {x, y}] // MatrixForm
4.例题
1)证明：(a*b)*c=(a.c)b-(b.c)a
Clear["Global`*"];
$Assumptions = {a \[Element] Vectors[3, Reals],
b \[Element] Vectors[3, Reals], c \[Element] Vectors[3, Reals]};
TensorReduce[(a \[Cross]b)\[Cross]c == (a.c) b - (b.c) a]
Clear["Global`*"];
$Assumptions = True;
a = {a1, a2, a3};
b = {b1, b2, b3};
c = {c1, c2, c3};
(a \[Cross]b)\[Cross]c // Simplify // MatrixForm
(a.c) b - (b.c) a // Simplify // MatrixForm
Expand[(a \[Cross]b)\[Cross]c == (a.c) b - (b.c) a]
FullSimplify[(a \[Cross]b)\[Cross]c == (a.c) b - (b.c) a]
从上面可以看出，使用TensorReduce来证明等式本身就没多大意义，相当于验证一下等式，因此推荐写出张量分量，一步一步展示中间过程，最后作判断。这样才能对等式有更直观的了解。
2)推导球坐标系下的 不可压缩Navier-Stokes 方程（动量方程）。即：
Clear["Global`*"];
u = {ur[r, \[Theta], z], u\[Theta][r, \[Theta], z],
uz[r, \[Theta], z]};
ADV = u.Transpose[Grad[u, {r, \[Theta], z}, "Cylindrical"]];
Pre = Grad[p[r, \[Theta], z], {r, \[Theta], z}, "Cylindrical"];
Vis = Laplacian[u, {r, \[Theta], z}, "Cylindrical"];
ADV[[1]]
Pre[[1]]
Vis[[1]] // Simplify
考虑到最终完整的形式很复杂且不是重点，我们只计算了其中的微分算子，其结果跟教材一致。
3)设一点沿着半径为a的球面等速运动，运动起始位置在赤道上，速度v的方向与球的经线夹角[Alpha]保持常数，求点的轨迹（斜驶线）方程，以及到达球极点的时刻[Tau]。
在球面上建立坐标系，考虑到坐标系正交，可直接根据拉梅系数获得速度表达式。
Clear["Global`*"];
r = a;
x = r Sin[\[Theta]] Cos[\[CurlyPhi]];
y = r Sin[\[Theta]] Sin[\[CurlyPhi]];
z = r Cos[\[Theta]];
H\[CurlyPhi] =
D[{x, y, z}, \[CurlyPhi]].D[{x, y, z}, \[CurlyPhi]] // Sqrt //
Assuming[0 <= \[Theta] <= \[Pi] && a > 0, Simplify@#] &
H\[Theta] =
D[{x, y, z}, \[Theta]].D[{x, y, z}, \[Theta]] // Sqrt //
Assuming[0 <= \[Theta] <= \[Pi] && a > 0, Simplify@#] &
速度可以表示为，vφvφ=HφHφ˙φφ˙，vθvθ=HθHθ˙θθ˙，根据速度条件tan[Alpha]=-vθvθ/vφvφ,积分一次可以给出斜驶线方程，以及求出路径长度与时间。
d\[CurlyPhi] = -Tan[\[Alpha]] H\[Theta] d\[Theta] /(H\[CurlyPhi]);
ds = (H\[CurlyPhi]^2 d\[CurlyPhi]^2 + H\[Theta]^2 d\[Theta]^2) //
Sqrt // Assuming[
0 <= \[Theta] <= \[Pi] && a > 0 && d\[Theta] > 0 &&
0 < \[Alpha] < \[Pi]/2, Simplify@#] &;
s = Integrate[-ds/d\[Theta], {\[Theta], \[Pi]/2, 0}]
\[Tau] = s/v
我们可以借助Mathematica来直观感受一下斜驶线。
Clear["Global`*"];
r = 1;
\[Alpha] = \[Pi]/4;
\[Theta] = 2 ArcTan[E^(-Cot[\[Alpha]] \[CurlyPhi])];
x = r Sin[\[Theta]] Cos[\[CurlyPhi]];
y = r Sin[\[Theta]] Sin[\[CurlyPhi]];
z = r Cos[\[Theta]];
p1 = ParametricPlot3D[{ Sin[\[Theta]1] Cos[\[CurlyPhi]1],
Sin[\[Theta]1] Sin[\[CurlyPhi]1], Cos[\[Theta]1]}, {\[CurlyPhi]1,
0, 2 \[Pi]}, {\[Theta]1, 0, \[Pi]}, BoxRatios -> {1, 1, 1},
ColorFunction -> {Green, Green},
Mesh -> None,
AxesOrigin -> {0, 0, 0},
Boxed -> False, Ticks -> {{-1, 1}, {-1, 1}, {-1, 1}},
ViewPoint -> {8, -4, -6}, PlotPoints -> 80,
ImageSize -> 500];
p2 = ParametricPlot3D[{x, y, -z}, {\[CurlyPhi], 0, 10},
PlotStyle -> {Thick, Blue}];
p3 = ParametricPlot3D[{ Cos[\[Theta]1], Sin[\[Theta]1],
0}, {\[Theta]1, 0, 2 \[Pi]}, PlotStyle -> {Thick, Yellow}];
Show[p1, p2, p3]
到此，如果以上函数都比较熟悉了，那么其实就可以尝试进行各种张量分析了，比如定义协变导数，张量微积分。更进一步，微分几何与相对论的东西也可以计算了。限于篇幅，本帖不可能面面俱到。如果读者有什么疑问或者想求解的问题，可以在下方留言，下期文章会整理做出解答。
本站仅提供存储服务，所有内容均由用户发布，如发现有害或侵权内容，请
点击举报。
