要搞 Lorentz 指标结构的投影算符就算是拿 Mathematica 暴力解方程也就那么几秒钟的事儿

• 原文: https://zhuanlan.zhihu.com/p/2043326743751156721
• 发布日期: 2026-05-29
• 分类: 粒子物理 / QFT / 场论计算

其实也没必要写但还是稍微提了两嘴的背景

也是挺奇妙的吧, 一天到晚不务正业我还以为自己早把物理全忘光了, 结果一写又发现自己还是好他妈懂啊. 都给我懂完了所以必须比比赖赖两下, 你可以直接跳过也没关系，，

在场论中我们可以用场算符 XJB 构造流算符, 比如可以用旋量场 $\psi \left(x\right)$ 构建双线性型 $\bar{\psi }\Gamma \psi$, 其中的 $\Gamma$ 一般就是 gamma 矩阵的乘积, 也可能带导数算符之类的. 更常见的操作应该就是拿夸克场去构建具有特定量子数的流算符, 而判断流算符量子数的方法大致可参考下文:

https://zhuanlan.zhihu.com/p/591789284像这样去构建具有确定夸克成分与量子数的流算符通常是会希望它能耦合到拥有相同夸克组分及量子数的粒子上, 所谓的耦合其实就是能对着真空搞出那么点儿来的意思. 比如说一个矢量流算符 $J^{\mu }\left(x\right)$ 如果能耦合到矢量粒子 $X$ 的话就意味着存在关系:

• $⟨\Omega |J^{\mu }\left(0\right)|X\left(\vec{p},\lambda \right)⟩=f_X{\epsilon }^{\mu }\left(\vec{p},\lambda \right)\ne 0.$

内积不为零就表明 $⟨\Omega |J^{\mu }\left(0\right)$ 中存在粒子态 $⟨X\left(\vec{p},\lambda \right)|$ 的分量, 这也就是说 $J^{\mu }\left(0\right)$ 能从真空中激发出 $X$ 来. 其中的 ${\epsilon }^{\mu }\left(\vec{p},\lambda \right)$ 是极化矢量, 它满足极化求和关系:

• \sum\nolimits_{\lambda =\pm ,0}{\varepsilon _{\mu }^{*}\left( \vec{p},\lambda \right){\varepsilon }_{\nu }\left( \vec{p},\lambda \right)}=\frac{p}_{\mu }{p}_{\nu }{m}^{2}-{g}_{\mu \nu }.

这里的来龙去脉既简单又复杂, 在此就不赘述了··· 总之我们要找到了这么个 $J^{\mu }\left(x\right)$, 那么就可以从它的两点关联函数 $⟨\Omega |\text{T}\left[J^{\mu }\left(x\right)J^{\nu †}\left(0\right)\right]|\Omega ⟩$ 中抽出矢量粒子 $X$ 的质量谱来.

如果你的 $J^{\mu }\left(x\right)$ 恰巧是个满足连续性方程 ${\partial }_{\mu }{J}^{\mu }\left(x\right)=0$ 的守恒流, 那么它在动量空间的两点关联函数 ${\Pi }^{\mu \nu }\equiv ⟨\Omega |\text{T}\left[J^{\mu }\left(x\right)J^{\nu †}\left(0\right)\right]|\Omega ⟩$ 自然就满足 Ward 恒等式:

• {p}_{\mu }{\Pi }^{\mu \nu }\left( p \right)\equiv {p}_{\mu }\text{i}\int{\text{d}^{4}x}\ {\text{e}^{\text{i}p\cdot x}\left\langle \Omega \right|\text{T}\left[ {J}^{\mu }\left( x \right){J}^{\nu \dagger }\left( 0 \right) \right]\left| \Omega \right\rangle
• \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =\int{\text{d}^{4}x}\ {\partial }_{\mu }{\text{e}^{\text{i}p\cdot x}\left\langle \Omega \right|\text{T}\left[ {J}^{\mu }\left( x \right){J}^{\nu \dagger }\left( 0 \right) \right]\left| \Omega \right\rangle
• \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =-\int{\text{d}^{4}x}\ {\text{e}^{\text{i}p\cdot x}\left\langle \Omega \right|\text{T}\left[ {\partial }_{\mu }{J}^{\mu }\left( x \right){J}^{\nu \dagger }\left( 0 \right) \right]\left| \Omega \right\rangle =0.\[1]

你知道的, 这里能构建 Lorentz 指标结构就只有 $p^{\mu }$ 和 $g^{\mu \nu }$, 然后又考虑到 $p_{\mu }{\Pi }^{\mu \nu }\left(p\right)=0$, 那请问 ${\Pi }^{\mu \nu }\left(p\right)$ 的结构是怎样的呢? 好难猜呀, 会事 \frac{p}^{\mu }{p}^{\nu }{p}^{2}-{g}^{\mu \nu } 吗，，

但要是 ${\partial }_{\mu }{J}^{\mu }\left(x\right)\ne 0$ 的话关联函数肯定就能分成像这样的两部分了:

• {\Pi }^{\mu \nu }\left( p \right)={\Pi }_{1}\left( {p}^{2} \right)\left( \frac{p}^{\mu }{p}^{\nu }{p}^{2}-{g}^{\mu \nu } \right)+{\Pi }_{0}\left( {p}^{2} \right)\frac{p}^{\mu }{p}^{\nu }{p}^{2}.

我们这里要研究矢量粒子, 也就是自旋为 1 的粒子, 那就还是得把 ${\Pi }_1\left(p^2\right)$ 单独提出来, 因为矢量粒子的极化都是四维横向的, 这里的横向体现在 {p}_{\mu }\left( \frac{p}^{\mu }{p}^{\nu }{p}^{2}-{g}^{\mu \nu } \right)=0 上, 所以所谓的四维横向指的就是与动量的方向正交.

这里我就不再详细展开了, 你可以结合前面的 $⟨\Omega |J^{\mu }\left(0\right)|X\left(\vec{p},\lambda \right)⟩=f_X{\epsilon }^{\mu }\left(\vec{p},\lambda \right)\ne 0$ 与下文给出的公式试图去理解它:

https://zhuanlan.zhihu.com/p/396673735那问题就来了, 怎么提取呢?

投影算符呗:

是的我用到了 FeynCalc, 我搞物理还从来没试过不用 FeynCalc.

这里的 % 与 %% 分别是上一个计算结果和上上个计算结果的意思, 所以 %*%% // Calc 的意思就是把前面算出来的两个结果乘在一起算一下.

也就是说直接拿投影算符 \frac{1}{3}\left( \frac{p}_{\mu }{p}_{\nu }{p}^{2}-{g}_{\mu \nu } \right) 怼上去就完了, 但你咋知道投影算符就长这样呢? 呃呃好吧, 其实这个还真就是用屁股都能猜出来，，那我再给你上点儿强度吧: 如果我们要找自旋为 3/2 的粒子的质量谱又会怎样呢?

自旋 1/2 的费米子没啥好说的, 但要是自旋 3/2 的粒子你怎么构建流算符呢? 就像自旋 1/2 与自旋为 1 分别对应到旋量场与矢量场那样, 自旋 3/2 也有它相应的 Rarita–Schwinger 场, 所以此时的插入流就会像 Rarita–Schwinger 场那样既有旋量指标又有矢量指标.

但这问题可就大了, 因为作为 Rarita–Schwinger 场的矢量旋量平面波解的 $u^{\mu }\left(p,\lambda \right)$ 是通过下面这三条关系收缩的自由度:

• \left( {p}\!\!\!/-m \right){u}^{\mu }=0,{\gamma }_{\mu }{u}^{\mu }=0,{p}_{\mu }{u}^{\mu }=0.\

这三条限制将 16 个分量压缩成了 4 个物理自由度对应到自旋 3/2 的部分, 进而还会有如下的极化求和关系:

• ${\sum }_s{u}^{\mu }\left(\vec{p},s\right){\bar{u}}^{\nu }\left(\vec{p},s\right)$
• =\left( {p}\!\!\!/+m \right)\left[ -{g}_{\mu \nu }+\frac{1}{3}{\gamma }_{\mu }{\gamma }_{\nu }+\frac{1}{3m}\left( {\gamma }_{\mu }{p}_{\nu }-{\gamma }_{\nu }{p}_{\mu } \right)+\frac{2}{3{m}^{2}{p}_{\mu }{p}_{\nu } \right].\

于是两点关联函数就可以按照如下形式展开:

• ${\Pi }^{\mu \nu }\left(p^2\right)$
• =\text{i}\int{\text{d}^{4}x}\ {\text{e}^{\text{i}p\cdot x}\langle \Omega |\text{T}\left[ {J}^{\mu }\left( x \right){J}^{\nu \dagger }\left( 0 \right) \right]\left| \Omega \right\rangle
• =\lambda _{3/2,-}^{2}\frac{p}\!\!/+{m}_{3/2,-}{m_{3/2,-}^{2}-{p}^{2}\left[ -{g}^{\mu \nu }+\frac{1}{3}{\gamma }^{\mu }{\gamma }^{\nu }+\frac{1}{3{m}_{3/2,-}\left( {\gamma }^{\mu }{p}^{\nu }-{\gamma }^{\nu }{p}^{\mu } \right)+\frac{2}{3m_{3/2,-}^{2}{p}^{\mu }{p}^{\nu } \right]
• \ \ \ +\lambda _{3/2,+}^{2}\frac{p\!\!/-{m}_{3/2,+}{m_{3/2,+}^{2}-{p}^{2}\left[ -{g}^{\mu \nu }+\frac{1}{3}{\gamma }^{\mu }{\gamma }^{\nu }+\frac{1}{3{m}_{3/2,+}\left( {\gamma }^{\mu }{p}^{\nu }-{\gamma }^{\nu }{p}^{\mu } \right)+\frac{2}{3m_{3/2,+}^{2}{p}^{\mu }{p}^{\nu } \right]
• \ \ \ +\Pi _{1/2}^{\mu \nu }\left( {p}^{2} \right)+\text{radial excitations}\ \text{and}\ \text{continuum}\text{.}\

好消息是这个 ${\Pi }_{1/2}^{\mu \nu }\left(p^2\right)$ 里绝对不含 $g^{\mu \nu }$ 这样的东西, 这是因为旋量矢量流 $J^{\mu }\left(x\right)$ 耦合到自旋 1/2 粒子时必须呈现 $⟨\Omega |J^{\mu }\left(0\right)|X_{1/2}⟩=C_1{p}^{\mu }+C_2{\gamma }^{\mu }$ 这样的结构, 要不然你的指标还能扔哪儿呢? 你的平面波旋量系数又没有 Lorentz 指标.

于是你再结合下文就不难发现这个 ${\Pi }_{1/2}^{\mu \nu }\left(p^2\right)$ 里绝对不会包含 $g^{\mu \nu }$ 这样的东西了, 毕竟已经有个指标让 $p^{\mu }$ 或者 ${\gamma }^{\mu }$ 给占去了嘛:

https://zhuanlan.zhihu.com/p/396673735这意味着什么? 这意味着只要你能拎出所有 $g^{\mu \nu }$ 的系数不就能卡掉所有的 1/2 污染留下我们要的 3/2 的部分了吗? 而且这个操作非常完美, 因为括号里的那一大串儿结构对应的系数都是完全相同的, 所以只要抽出其中一项的系数就够用了.

类比到前面那些内容也是一样的道理, 前面我们面对的体系是这样的:

{\Pi }^{\mu \nu }\left( p \right)={\Pi }_{1}\left( {p}^{2} \right)\left( \frac{p}^{\mu }{p}^{\nu }{p}^{2}-{g}^{\mu \nu } \right)+{\Pi }_{0}\left( {p}^{2} \right)\frac{p}^{\mu }{p}^{\nu }{p}^{2}.

那对这个体系你可以直接用投影算符抽出 \frac{p}_{\mu }{p}_{\nu }{p}^{2}-{g}_{\mu \nu } 部分的系数 ${\Pi }_1\left(p^2\right)$ 来, 也可以干脆就提取 $g_{\mu \nu }$ 的系数嘛, 大不了就补个负号的事儿···

甚至你会发现其实这个投影算符就是在取 $g^{\mu \nu }$ 的系数:

好的讲了那么多屁话, 问题终于来了——你能轻松地抽出所有 $g^{\mu \nu }$ 的系数吗?

直接拿投影算符 \frac{1}{3}\left( \frac{p}_{\mu }{p}_{\nu }{p}^{2}-{g}_{\mu \nu } \right) 抽可事行不通的噢:

很显然 $\Pi \left(1\right)$ 的系数 ${\gamma }^{\mu }.{\gamma }^{\nu }$ 就筛不掉···

那我们究竟该如何去寻找想要的投影算符呢?

答案就是暴力算, 是的就直接用 Mathematica 的 Solve 就能解出来，，

首先写出所有我们可能会碰到的结构, 并依次命名为 $\Gamma \left(n\right)$, 然后再写一个能 general 地包含所有结构的多项式:

• LorBasis={[CapitalGamma][1]=GA[[Mu],[Nu]],[CapitalGamma][2]=MT[[Mu],[Nu]],[CapitalGamma][3]=GA[[Mu]]*FV[p,[Nu]],[CapitalGamma][4]=GA[[Nu]]*FV[p,[Mu]],[CapitalGamma][5]=FV[p,[Mu]]*FV[p,[Nu]]}
• [CapitalPi][All]=Sum[[CapitalGamma][n]*[CapitalPi][n],{n,1,Length[LorBasis]}]

再就是分别用每一个结构去收缩这个包含所有结构的多项式并将结果依次命名为 Eq(n):

• Table[Eq[nEq]==Calc[\[CapitalPi][All].\[CapitalGamma][nEq]],{nEq,1,Length[LorBasis]}]//TableForm

接着我们再把这些式子连立起来用 Solve 反解出所有 $\Pi \left(n\right)$ 关于 Eq(n) 的表达式:

• Solve[Apply[And,Table[Eq[nEq]==Calc[\[CapitalPi][All].\[CapitalGamma][nEq]],{nEq,1,Length[LorBasis]}]],Table[\[CapitalPi][n],{n,1,Length[LorBasis]}]]//Simplify//Flatten[#,1]&//TableForm

然后就只要把所有的 Eq(n) 都替换为相应的结构 $\Gamma \left(n\right)$ 就能得到我们要的投影算符了:

• EqDot/:EqDot[AnyShit_]Coeff_:=EqDot[AnyShit].Coeff;
• ReplaceAll[{([CapitalPi]n_->Oper_):>([CapitalPi]n->DotSimplify[DiracSimplify[Oper]])}][ReplaceAll[{EqDot->[CapitalGamma]}][ReplaceAll[{Eq->EqDot}][Solve[Apply[And,Table[Eq[nEq]==Calc[[CapitalPi][All].[CapitalGamma][nEq]],{nEq,1,Length[LorBasis]}]],Table[[CapitalPi][n],{n,1,Length[LorBasis]}]]]]]//Flatten[#,1]&//TableForm

上面定义了个 EqDot 是为了能正确地引入 gamma 矩阵之间的乘积, 而这个函数的效果就是无论你给它乘个啥都会变成矩阵点乘:

然后再通过两次替换——从 Eq 换成 EqDot 再从 EqDot 转到 $\Gamma$ ——就能将所有的 $\Gamma \left(n\right)$ 都正确地以矩阵乘法的形式加入式子并进行运算得到所有我们想要的那些投影算符了.

代码我反正是给你讲清楚了, 至于说这个方程··· 那这要看不懂就是你自己的问题了, 这边儿的建议事退学呢，，

总之打包好的代码如下, 你就算看不懂也能用:

• LorBasis={[CapitalGamma][1]=GA[[Mu],[Nu]],[CapitalGamma][2]=MT[[Mu],[Nu]],[CapitalGamma][3]=GA[[Mu]]*FV[p,[Nu]],[CapitalGamma][4]=GA[[Nu]]*FV[p,[Mu]],[CapitalGamma][5]=FV[p,[Mu]]*FV[p,[Nu]]};
• [CapitalPi][All]=Sum[[CapitalGamma][n]*[CapitalPi][n],{n,1,Length[LorBasis]}];
• EqDot/:EqDot[AnyShit_]Coeff_:=EqDot[AnyShit].Coeff;
• [CapitalPi]sAndOperators=ReplaceAll[{([CapitalPi]n_->Oper_):>([CapitalPi]n->DotSimplify[DiracSimplify[Oper]])}][ReplaceAll[{EqDot->[CapitalGamma]}][ReplaceAll[{Eq->EqDot}][Solve[Apply[And,Table[Eq[nEq]==Calc[[CapitalPi][All].[CapitalGamma][nEq]],{nEq,1,Length[LorBasis]}]],Table[[CapitalPi][n],{n,1,Length[LorBasis]}]]]]];
• FCEId[Term_]:=FCE[Term];
• ProjectionOperator[n[CapitalPi]]:=FCEId[Part[#,Apply[Sequence,Flatten[Position[#,[CapitalPi][n[CapitalPi]]->]]],2]&[[CapitalPi]sAndOperators]];

其中那个叫 ProjectionOperator 的函数就是用来给出任意投影算符的, 其变量对应各个结构的编号. 比如说你要取出 $\Gamma \left(2\right)=g^{\mu \nu }$ 来的话, 那就输入 ProjectionOperator[2] 就可以得到你要的投影算符:

要算其它体系的投影算符只需改改第一行的指标结构就行, 比如更简单的:

不过要记得这些投影算符必须得放右边, 有时候能放左边仅仅就只是碰巧罢了.

你可能会说那我几把直接 Coefficient 提系数不完了?

嗯, 你说得对, 但要有人把其中几个 ${\gamma }^{\mu }{\gamma }^{\nu }$ 写成 $2g^{\mu \nu }-{\gamma }^{\nu }{\gamma }^{\mu }$ 你不炸了:

不过其实还是可以这么做的, 只要你提前 DiracOrder 一下就行:

当然你用投影算符肯定无论如何都不会出错吧:

但其实我也会选择 DiracOrder + Coefficient 而不是乘个投影算符, 因为我平时工作动不动就要算你吗上万项··· 那乘上万个投影算符花的时间我都不敢想，，

你也可以写个或许看起来更无脑些的函数?

• FCEId[Term_]:=FCE[Term];
• ProjectionOperator[Basis_,n[CapitalPi]Tab_]:=Module[{EqDot,[CapitalGamma],[CapitalPi],LorBasis,[CapitalPi]sAndOperators,Eq},EqDot/:EqDot[AnyShit_]Coeff_:=EqDot[AnyShit].Coeff;
• LorBasis=Table[[CapitalGamma][n[CapitalGamma]]=Part[#,n[CapitalGamma]],{n[CapitalGamma],Length[#]}]&[Basis];
• [CapitalPi][All]=Sum[[CapitalGamma][n[CapitalGamma]]*[CapitalPi][n[CapitalGamma]],{n[CapitalGamma],1,Length[LorBasis]}];
• [CapitalPi]sAndOperators=ReplaceAll[{([CapitalPi]n_->Oper_):>([CapitalPi]n->DotSimplify[DiracSimplify[Oper]])}][ReplaceAll[{EqDot->[CapitalGamma]}][ReplaceAll[{Eq->EqDot}][Solve[Apply[And,Table[Eq[nEq]==Calc[[CapitalPi][All].[CapitalGamma][nEq]],{nEq,1,Length[LorBasis]}]],Table[[CapitalPi][n[CapitalPi]],{n[CapitalPi],1,Length[LorBasis]}]]]]];
• If[Head[n[CapitalPi]Tab]===Integer,FCEId[Part[#,Apply[Sequence,Flatten[Position[#,[CapitalPi][n[CapitalPi]Tab]->]]],2]]&[[CapitalPi]sAndOperators],TableForm[Table[FCEId[ReplaceAll[{[CapitalPi]->[CapitalGamma]}][Part[#,Apply[Sequence,Flatten[Position[#,[CapitalPi][n[CapitalPi]]->]]]]]]&[[CapitalPi]sAndOperators],{n[CapitalPi],n[CapitalPi]Tab}]]]];

上面这款的话就什么都不用修改了, 输入涉及的结构及序号就能输出投影算符:

如上图所示, 你就只需分别填写结构及序号, 而序号可以只填个数字也可以写个数表.

当然这里生成也都是些右投影算符, 也就是要放在右边.

Mathematica 是真够方便的, 我操，，

呃呃我想了想还是随手补个数学层面的证明过程吧

首先摆出所有涉及到的 Lorentz 指标结构, 即 \left\{ \Gamma _{a}^{\mu \nu \cdots \rho }|a\in {\mathbb{N}_{+} \right\}.

比如说 \left\{ \Gamma _{1}^{\mu \nu }={\gamma }^{\mu }{\gamma }^{\nu },\Gamma _{2}^{\mu \nu }={g}^{\mu \nu },\Gamma _{3}^{\mu \nu }={\gamma }^{\mu }{p}^{\nu },\Gamma _{4}^{\mu \nu }={\gamma }^{\nu }{p}^{\mu },\Gamma _{5}^{\mu \nu }={p}^{\mu }{p}^{\nu } \right\}.\

于是一个通项就可以被表达为 ${\Pi }^{\mu \nu \cdots \rho }\equiv {\sum }_a{\Pi }^a{\Gamma }_a^{\mu \nu \cdots \rho }\equiv {\Pi }^a{\Gamma }_a^{\mu \nu \cdots \rho }$, 而所谓的投影算符 $P_{\mu \nu \cdots \rho }^a$ 就是可以办到 $P_{\mu \nu \cdots \rho }^a{\Pi }^{\mu \nu \cdots \rho }={\Pi }^b{P}_{\mu \nu \cdots \rho }^a{\Gamma }_b^{\mu \nu \cdots \rho }={\Pi }^a$ 的算符.

这就意味着 $P_{\mu \nu \cdots \rho }^a{\Gamma }_b^{\mu \nu \cdots \rho }={\delta }_b^a$, 你会发现这个方程右边是没有 Lorentz 指标的.

这自然就会想到要定义矩阵元 $M_b^a\equiv {\Gamma }_{\mu \nu \cdots \rho }^a{\Gamma }_b^{\mu \nu \cdots \rho }$ , 接着就显然有:

• $P_{\mu \nu \cdots \rho }^a{\Gamma }_b^{\mu \nu \cdots \rho }={\delta }_b^a={\left(M^{-1}\right)}_c^a{M}_b^c={\left(M^{-1}\right)}_c^a{\Gamma }_{\mu \nu \cdots \rho }^c{\Gamma }_b^{\mu \nu \cdots \rho }$
• \Rightarrow P_{\mu \nu \cdots \rho }^{a}={\left( {M}^{-1} \right)}^{a}_{b}\Gamma _{\mu \nu \cdots \rho }^{b}.\

所以从头到尾要求的仅仅就只是个逆矩阵 $M^{-1}$ 罢了, 其中 $M_b^a\equiv {\Gamma }_{\mu \nu \cdots \rho }^a{\Gamma }_b^{\mu \nu \cdots \rho }.$

感觉这就是线性代数的入门级内容吧, 不过估计多数人都解不出来所以我就还是补了个过程.

配合这个证明过程就很好理解我前面的代码中的操作了吧:

我就是先定义了 \text{E}{\text{q}_{a}\equiv {\Pi }_{\mu \nu }\Gamma _{a}^{\mu \nu }={\Pi }_{b}\Gamma _{\mu \nu }^{b}\Gamma _{a}^{\mu \nu }\equiv {\Pi }_{b}{M}^{b}_{a},a\in {\mathbb{N}_{+} 嘛, 然后再解方程组试图将 ${\Pi }_a$ 展开为 \text{E}{\text{q}_{a} 的线性组合, 这其实就是求逆得到 {\Pi }_{a}=\text{E}{\text{q}_{b}{\left( {M}^{-1} \right)}^{b}_{a}.

所以我的下一步就是将 \text{E}{\text{q}_{b} 都换成 ${\Gamma }_b^{\mu \nu }$ 于是就得到了 ${\Gamma }_b^{\mu \nu }{\left(M^{-1}\right)}_a^b\equiv P_a^{\mu \nu }.$

被你发现了, ${\Gamma }_b^{\mu \nu }{\left(M^{-1}\right)}_a^b\ne {\left(M^{-1}\right)}_b^a{\Gamma }_{\mu \nu }^b$ 对吧? 其实一般来说最后搞出的投影算符都是要扣帽子分左右的, 所以前面写的那些才不能放左边, 那能不能写一款左投影算符的呢?

可以但···

你知道的, 左的都不是什么好东西, 包括投影算符.

所以理论上直接构造矩阵求逆也可以写一套代码

• FCEId[Term_]:=FCE[Term];
• ProjectionOperator[LR_,Basis_,nPO_]:=Module[{NumOfBas,InverBB},NumOfBas=Length[Basis];InverBB=Calc[Inverse[Table[Calc[Part[Basis,aBas].Part[Basis,bBas]],{aBas,1,NumOfBas},{bBas,1,NumOfBas}]]];If[Head[nPO]===Integer,#[nPO],TableForm[Table[Part[Basis,nBas]->#[nBas],{nBas,nPO}]]]&[FCEId[Calc[Sum[If[LR===L,Part[InverBB,#,aBas].Part[Basis,aBas],Part[Basis,aBas].Part[InverBB,aBas,#]],{aBas,1,NumOfBas}]]]&]];

这个版本就非常终极了, 按理说这么牛逼的东西不应该免费的呀:

如图所示, ProjectionOperator 函数里就先扣帽子分左右, 再填齐结构及序号即可, 而其中的序号部分还是跟以前一样数字、数表皆可填的.

参考

• ^倒数第二步是有点儿微妙但这里就先别管了吧.