嵌入式系统在人类生活中发挥着越来越大的作用, 而作为嵌入式系统灵魂的嵌入式软件在其中所占有的比重也越来越大.因此, 嵌入式软件的可靠性将变得更加重要.诸如航空、航天、军事、交通、医疗等关键应用领域都对嵌入式系统的可靠性和安全性要求非常高, 任何错误的发生都可能带来灾难性后果.嵌入式系统具有3个重要属性, 即可达性、终止性、不变式.可达性是指一个系统能否从一个给定状态到达另一个可接受状态.某些混成系统的可达性被证明是能用计算机代数工具来检验的.不变式是系统变量在循环迭代时永远保持的特性.终止性是研究系统中是否会发生死循环.不包括终止性分析的验证被称为程序的部分正确性证明.因此, 程序的终止性分析是确保程序完全正确性的必要基础.

尽管程序的终止性问题早已被证明是不可判定的^[1], 但对其进行研究不仅具有理论意义, 更具有实际意义.当前, 国内外主要通过计算秩函数来进行循环终止性分析.秩函数是验证循环程序终止性的一条重要途径.对给定的一个程序, 若能找到其秩函数, 则表明该程序必然是终止的.例如在2001年, Colon等人^{[2, 3]}就利用多面体理论计算线性程序的线性秩函数(ranking functions).2004年, 针对线性循环程序, Poldelski等人^[4]提出一个完备的线性秩函数生成方法, 并开发出工具RANKFINDER.他们的方法将线性程序的线性秩函数计算问题归结为线性规划(LP)问题.在2005年, 针对线性程序, Manna等人^{[5, 6]}提出了字典秩函数的概念, 并利用Farkas引理呈现了线性字典秩函数的计算方法.最近, 文献[7, 8]研究了现有几种线性秩函数生成算法的时间复杂度.文献[8]中进一步分析了程序变量在整数集合上取值时, 这类线性程序的线性秩函数合成方法, 并研究了这类问题的复杂度.但即便是线性程序, 其秩函数也未必是线性的.因此, 针对非线性程序, Cousot^[9]通过半正定规划方法SDP研究了多项式秩函数的计算问题.该方法能够找到高次数的秩函数, 但由于SDP工具内部采用的是数值计算, 因此由计算误差导致最终得到的函数未必是真正的秩函数.同时, 这个方法也不能回答一个程序是否具有预定义形式的秩函数.

符号计算理论方法及其工具被逐渐应用于程序自动验证.例如, 我国科学家陈应华等人和杨路等人将秩函数和不变式的计算归约为半代数系统的求解, 并运用基于柱形代数分解的工具DISCOVERER, 提出了多项式不等式型不变式和非线性秩函数生成方法^{[2, 10]}.不同于文献[9]中的方法, 他们的方法能够精确回答循环程序是否有给定模板的秩函数或不变式.

除了秩函数法以外, 还出现了其他方法去进行循环程序的终止性研究.如, 文献[11, 12]提出了试探性方法去探索给定循环程序的非终止性.众所周知, 一般的程序终止性问题是不可判定的, 但对某些类特殊的程序而言, 人们总希望能证明其终止性问题是可判定的.由此, 鉴于终止性问题本身的不可判定性所带来的实际困难, 证明循环程序终止性的另一途径就是避开秩函数的合成, 而采用数学方法严格证明某类或某些类循环程序的终止性是可判定的, 并建立相对完备的判定算法.这首先需要对程序进行分类.当限定程序中的赋值语句和条件语句均为线性多项式时, 2004年, Tiwari在文献[13]中首次证明了下列单重无分支线性循环程序在实域上的终止性是可判定的:

$ {\rm{ }}{\mathbf{while}}\;\;{\rm{ }}Bx > 0\;\;{\rm{ }}{\mathbf{do}}{\rm{ }}\;\;\{ x: = Ax\} {\rm{ }}\;\;{\mathbf{endwhile}}{\rm{, }} $

这里, A, B均是实矩阵.相似的结论在2006年被Braveman^[14]推广到整数环上.此外, 为了避免Jordan标准型的浮点计算, 文献[15]中提出了精确的符号计算方法对程序P进行终止性判定.既然一般形式的线性程序终止性是不可判定的^[13], 那么非线性多项式循环程序由于其更为复杂的动力行为使得其终止性分析将变得更加困难.一个程序被称为非线性的, 是指循环中的赋值映射或循环条件中的约束是非线性表达式.2005年, 利用有限差分树理论, 文献[16]提出了试探性算法对一类含多分支语句的多项式程序的终止性问题进行判定.2010年, 针对一类赋值为线性且循环条件由非线性多项式不等式构成的循环程序, 文献[17]分析了其终止性问题, 证明了当这类程序满足给定的NZM条件时, 其终止性是可判定的.2013年, 文献[18]通过分析多项式映射fi的发散区间讨论了一类多项式循环(迭代赋值型如x_i:=f_i(x_i))的终止性问题.针对一类循环条件为多项式等式, 赋值语句为多项式的多分支循环程序, 詹乃军等人在文献[19]中给出了其可终止或不可终止的充分条件.最近, 针对这类多分支程序.他们进一步证明了各条路径形成的理想具有一致上界, 从而证明了这类多分支程序的终止性是可判定的.记R为实数域, Q为有理数域.R[x](Q[x])表示由关于x的实系数(有理系数)多项式构成的多项式环.

本文主要研究具有以下形式的多分支线性循环程序在实数域上的终止性问题.

$\left. \begin{gathered} P{\rm{ }}~{\mathbf{while}}~{\rm{ }}Bx \geqslant b{\rm{ }}~{\mathbf{do}} \\ {\tau _1}:\mathit{\boldsymbol{x}}: = {F_1}(\mathit{\boldsymbol{x}}); \\ {\rm{or}} \\ ... \\ {\rm{or}} \\ {\tau _m}:\mathit{\boldsymbol{x}}: = {F_m}(\mathit{\boldsymbol{x}}); \\ \end{gathered} \right\}$

(1)

这里, B∈Q^n×n为一可逆矩阵, b∈Qⁿ, x∈Rⁿ, F_i(x)=(f_i1(x), …, f_in(x))^T为关于x的n维多项式映射, f_ij(x)∈Q[x], i=1, ..., m, j=1, ..., n.记Ω={x∈Rⁿ:Bx≥b}.为方便, 令P$ \overset{\mathit{\Delta} }{\mathop{=}}\,$P(Ω, F₁, …, F_m).在上述多分支循环程序中, 我们忽略了每个分支的条件判断语句, 这为多分支循环的终止性研究带来了便利.这一处理方式并非新的, 它早已存在于文献[19]中.

对任意给定的v个n维多项式映射g₁(x), …, g_v(x), 令$ \oplus _{j = 1}^v{g_j}(\mathit{\boldsymbol{x}}) = {g_v}(\mathit{\boldsymbol{x}}) \circ ... \circ {g_1}(\mathit{\boldsymbol{x}})$个映射的复合.我们说多分支循环程序P在实域上是不可终止的, 如果存在x₀∈Rⁿ和无穷下标列表$\mathit{\Lambda} = [{i_j}]_{j = 1}^{ + \infty }({i_j} \in \{ 1, 2, ..., m\} )$, 使得对任意的非负整数n, 有:

$ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) = {F_{{i_n}}} \circ ... \circ {F_{{i_2}}}^\circ {F_{{i_1}}}^\circ {F_{{i_0}}}({\mathit{\boldsymbol{x}}_0}) \in \mathit{\Omega}. $

特别地, 记${F_{{i_0}}} = I$为恒等映射.这里, ${F_{{i_j}}}$∈{I, F₁, ..., F_m}, j=1, 2, ...如果那样的x₀不存在, 则表明程序P是可终止的.目前, 对程序P, 往往采用合成秩函数(ranking functions)的方法来进行终止性分析^{[2, 3, 5-9, 20]}.也即当程序P中各个分支的赋值语句都是线性表达式时, 基于Farkas引理, 文献[3, 5, 6, 20]中的方法能够被用于计算P的线性秩函数.一个线性循环程序可能存在秩函数, 但其秩函数未必是线性的.更重要的是, 上述基于Farkas引理的方法并不能计算非线性循环程序的线性(或非线性)秩函数.文献[2, 9]中的方法相同点在于:事先都要设定一个带参数的秩函数模板, 然后通过不同的计算方法去计算得到参数的一个值.具体地说, 给定某多项式循环程序, 文献[9]将秩函数计算问题转化为半正定规划(SDP)问题, 从而利用半正定规划方法和工具去计算该循环程序的多项式秩函数.但由于SDP采用浮点算法会导致计算得到的秩函数有可能并不是真正的秩函数, 因此当SDP求解器计算出一个参数向量的取值时, 还需要最后验证其对应的函数是否为程序的秩函数.不同于文献[9]的方法, 文献[2]的方法采用柱形代数分解(CAD)算法的量词消去方法去计算参数模板的参数空间.他们的方法是纯符号的精确方法, 所有的计算都不会涉及到近似计算.但由于CAD算法在最糟糕情形下的复杂度被证明是双指数的, 故当秩函数模板次数稍高(大于2次)或变元稍多(大于3个)时, 计算时间往往不可接受.本文中, 基于Polya定理, 我们提出了一种新的纯符号方法去探测程序P的多项式秩函数.该方法将程序P的多项式秩函数计算问题归结为线性规划(LP)求解问题.尽管求解LP问题存在多项式时间算法, 如椭球算法和内点算法, 但类似于上述的SDP方法, 这些算法由于采用浮点计算从而导致所得到的结果未必精确可信.因此本文中, 我们将采用线性规划求解算法中的经典算法——单纯形算法, 去计算得到可行域中的一个精确点.同时, 采用单纯形算法, 我们不仅能够探测含次数更高, 变元更多的复杂多项式秩函数, 而且保证了计算所得到的函数就是程序的秩函数, 从而不必像文献[9]中的方法那样需要最后验证所求得的函数是否为秩函数.

考虑程序P.令y=M^-1(x)=Bx-b, Ωy={y:y≥0}.既然程序P中的矩阵B是可逆的, 则M^-1为可逆映射, 记M^-1的逆映射为M(y)=B^-1(y+b).因此, M, M^-1互为Rⁿ$ \to $Rⁿ的可逆映射.显然有, M^-1是Ω到Ω_y上的双射.同样地, M是Ω_y到Ω上的双射.故下面的定理表明P在实数域上的终止性问题等价于下列程序在实数域上的终止性问题.

$ \left. \begin{align} & U~\mathbf{while}~y\ge 0~\mathbf{do} \\ & {{\tau }_{1}}:\mathit{\boldsymbol{y}}:={{M}^{-1}}{}^\circ {{F}_{1}}{}^\circ M(\mathit{\boldsymbol{y}});\text{ } \\ & \text{or } \\ & \cdots \text{ } \\ & \text{or} \\ & \text{ }{{\tau }_{m}}:\mathit{\boldsymbol{y}}:={{M}^{-1}}{}^\circ {{F}_{m}}{}^\circ M(\mathit{\boldsymbol{y}}); \\ \end{align} \right\} $

(2)

与程序P的不可终止定义类似, 我们说程序U在实域上是不可终止的, 如果存在y₀∈Rⁿ和无穷下标列表$ \wedge = [{i_j}]_{j = 1}^{ + \infty }({i_j} \in \{ 1, 2, ..., m\} )$, 使得对任意的非负整数n, 有:

$\mathop \oplus \limits_{j = 0}^n ({M^{ - 1}} \circ {F_{ij}} \circ M({\mathit{\boldsymbol{y}}_0})) = ({M^{ - 1}} \circ {F_{{i_n}}} \circ M) \circ ... \circ ({M^{ - 1}} \circ {F_{{i_1}}} \circ M) \circ ({M^{ - 1}} \circ {F_{{i_0}}} \circ M)({\mathit{\boldsymbol{y}}_0}) \in {\mathit{\Omega} _\mathit{\boldsymbol{y}}}{\rm{.}}$

这里, 对任意的$j = 0, 1, ..., {M^{ - 1}} \circ {F_{{i_j}}} \circ M \in \{ {M^{ - 1}} \circ {F_{{i_0}}} \circ M = I, {M^{ - 1}} \circ {F_1} \circ M, ..., {M^{ - 1}} \circ {F_m} \circ M\} $.

定理1.记号同上.程序P在实数域上不可终止等价于程序U在实数域上不可终止.

证明:若程序P在实数域上不可终止, 则必存在x₀∈Rⁿ和无穷下标列表$ \wedge = [{i_j}]_{j = 1}^{ + \infty }({i_j} \in \{ 1, 2, ..., m\} )$, 使得对任意的n≥0, 有:

$ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) = {F_{{i_n}}} \circ ... \circ {F_{{i_2}}} \circ {F_{{i_1}}} \circ {F_{{i_0}}}({\mathit{\boldsymbol{x}}_0}) \in \mathit{\Omega} $

(3)

这里, ${F_{{i_j}}} \in \{ {F_1}, ..., {F_m}\} , j = 1, 2, ..., $且$ i_{j} \in \mathit{\Lambda}, j=1, 2, \ldots$令y₀=M^-1(x₀)=Bx₀-b, 同时有x₀=M(y₀).因此对任意的非负整数n≥0, 有:

$ \begin{aligned} \mathop \oplus \limits_{j = 0}^n \left(M^{-1} \circ F_{i_{j}} \circ M\right)\left(\boldsymbol{y}_{0}\right) &=\left(M^{-1} \circ F_{i_{n}} \circ M\right) \circ \ldots \circ\left(M^{-1} \circ F_{i_{1}} \circ M\right) \circ\left(M^{-1} \circ F_{i_{0}} \circ M\right) \\ &=M^{-1} \circ F_{i_{n}} \circ F_{i_{n-1}} \circ \ldots \circ F_{i_{0}} \circ M\left(\boldsymbol{y}_{0}\right) \\ &=M^{-1} \circ F_{i_{n}} \circ F_{i_{n-1}} \circ \ldots \circ F_{i_{0}}\left(\boldsymbol{x}_{0}\right) \\ &=M^{-1} \circ F_{i_{j}}\left(\boldsymbol{x}_{0}\right) \end{aligned} $

(4)

这里, $ i_{j} \in \mathit{\Lambda}, j=1, 2, \ldots$, ${F_{{i_0}}} = I{\rm{.}}$根据公式(3)、公式(4)中$ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) \in \mathit{\Omega} {\rm{, }}$既然M^-1是从Ω到Ω_y上的双射, 我们有:

${M^{ - 1}} \circ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) \in \mathit{\Omega}. $

也即对任意的非负整数$n \geqslant 0, \oplus _{j = 0}^n({M^{ - 1}} \circ {F_{{i_j}}} \circ M)({\mathit{\boldsymbol{y}}_0}) \in {\mathit{\Omega} _\mathit{\boldsymbol{y}}}$.因此, 存在一个点y₀和一个无穷下标列表$\mathit{\Lambda} = [{i_j}]_{j = 1}^{ + \infty }({i_j} \in \{ 1, 2, ..., m\} )$, 使得对任意n≥0, 我们有$ \oplus _{j = 0}^n({M^{ - 1}} \circ {F_{{i_j}}} \circ M)({\mathit{\boldsymbol{y}}_0}) \in {\mathit{\Omega} _\mathit{\boldsymbol{y}}}{\rm{.}}$根据程序不可终止的定义得知, 程序在y₀处不可终止.

反之, 若程序U在实数域上不可终止, 那么必存在一点y₀和一个无穷下标列表$\mathit{\Lambda} = [{i_j}]_{j = 1}^{ + \infty }({i_j} \in \{ 1, 2, ..., m\} )$, 使得对任意的n≥0, 我们有$ \oplus _{j = 0}^n({M^{ - 1}} \circ {F_{{i_j}}} \circ M)({\mathit{\boldsymbol{y}}_0}) \in {\mathit{\Omega} _\mathit{\boldsymbol{y}}}{\rm{.}}$令x₀=M(y₀).根据公式(4)和M^-1的可逆性, 有:

$M \circ \oplus _{j = 0}^n({M^{ - 1}} \circ {F_{{i_j}}} \circ M)({\mathit{\boldsymbol{y}}_0}) = M \circ {M^{ - 1}} \circ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) = \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}){\rm{.}}$

既然$ \oplus _{j = 0}^n({M^{ - 1}} \circ {F_{{i_j}}} \circ M)({\mathit{\boldsymbol{y}}_0}) \in {\mathit{\Omega} _\mathit{\boldsymbol{y}}}{\rm{, }}$根据M是从Ω_y到Ω上的双射的性质可知, $ \oplus _{j = 0}^n{F_{{i_j}}}({\mathit{\boldsymbol{x}}_0}) \in \mathit{\Omega} {\rm{.}}$既然这里的n是任意非负整数, 故按其不可终止的定义, 程序U在x₀上不可终止.

根据定理1可知, P的终止性问题总能等价归结为U的终止性问题.为方便分析, 不失一般性, 将原来U中的程序变元y替换为x, M^-1$ \circ $F_i$ \circ $M(y)替换为T_i(x), 则程序U可被重写为以下形式.

$\left. \begin{gathered} U{\rm{ }}~{\mathbf{while}}~{\rm{ }}x \geqslant 0{\rm{ }}~{\mathbf{do}} \\ {\tau _1}:\mathit{\boldsymbol{x}}: = {T_1}(\mathit{\boldsymbol{x}}); \\ {\rm{or}} \\ ... \\ {\rm{or}} \\ {\tau _m}:\mathit{\boldsymbol{x}}: = {T_m}(\mathit{\boldsymbol{x}}); \\ \end{gathered} \right\}$

(5)

这里, x=(x₁, x₂, …, x_n)^T∈Rⁿ, T_i(x)=(t_i1(x), …, t_in(x))^T为关于x的n维多项式映射, t_ij(x)∈Q[x], i=1, …, m, j=1, …, n.记:

$ \mathit{\Omega}_{x}=\left\{\boldsymbol{x} \in \mathbf{R}^{n}: x_{1} \geqslant 0, \dots, x_{n} \geqslant 0\right\}. $

既然P等价于程序U的终止性, 因此为判定P的终止性, 我们可以先分析U的终止性.

下文中, 我们将通过计算多项式秩函数的方法去判定公式(5)中的程序U的终止性.

1 预备知识

本节将介绍有关秩函数、Polya定理的相关基本知识.

定义1.给定如公式(5)定义的循环程序U.函数ρ(x)∈R[x]被称为程序U的秩函数, 如果存在一个正数c > 0, 使得ρ(x)满足:

● (有界): $ \forall \mathit{\boldsymbol{x}}$, (x∈Ω_x $ \Rightarrow $ρ(x)≥0);

● (下降): $ \mathop \wedge \limits_{i = 1}^m (\forall \mathit{\boldsymbol{x}}\forall \mathit{\boldsymbol{x}}', (\mathit{\boldsymbol{x}} \in {\mathit{\Omega} _\mathit{\boldsymbol{x}}} \wedge \mathit{\boldsymbol{x}}' = {T_i}(\mathit{\boldsymbol{x}}) \Rightarrow \rho (\mathit{\boldsymbol{x}}) - \rho (\mathit{\boldsymbol{x}}') \geqslant c)).$

注:根据定义1及下面的定义3, 计算程序U的秩函数等价于在Ω_x上探测m+1个半正定多项式:

$ \rho(\boldsymbol{x})-\rho\left(T_{i}(\boldsymbol{x})\right)-c. $

众所周知, 如果一个循环程序具有秩函数, 那么这个程序必定是终止的.对程序U, 文献[3, 5, 6, 20]中的方法均不能计算其秩函数.这是因为U中的表达式可以是非线性多项式表达式, 而上述文中基于Farkas引理的方法仅能计算线性循环程序的线性秩函数.针对程序U, 目前有两种方法可以计算它的非线性多项式秩函数.

●一种是由杨路等人提出的基于符号计算的多项式秩函数计算方法, 他们的方法是纯符号的方法, 并不涉及近似计算.这一方法的特点是能够完备判定给定的程序是否具有预设模板形式的秩函数, 但该方法是基于复杂度较高的柱形代数分解算法的, 故可用以计算含变元较少、次数较低的秩函数合成问题.

●另一种是由Cousot提出的基于SDP的秩函数计算方法.利用SDP的求解是多项式时间复杂度, 该方法能够探测到含有较高次数、较多变元的秩函数.但不足之处在于, 由于SDP采用浮点计算, 从而导致最终得到的秩函数可能是不精确的, 故还需要进一步验证.

在下一节中, 我们将给出一种新的多项式秩函数计算方法.我们的方法是基于下面的Polya定理.首先给出一些多项式相关的定义.

定义2.给定多元齐次多项式f(x)∈R[x].我们说f是关于x的齐次多项式, 如果其中的所有非零系数项的次数都相同.

例如, $f({x_1}, {x_2}) = 2x_1^2 - 6{x_1}{x_2}$就是关于x₁, x₂的齐次多项式, 其次数为2.特别地, 如果f是非齐次多项式, 那么f中次数最高的单项式的次数被称为f的次数.例如, 非齐次多项式f(x₁, x₂)=x₁x₂+3x₁-x₂+1中, 单项式(忽略系数)有x₁x₂, x₁, x₂, 1, 其关于x₁, x₂的次数分别为2, 1, 1, 0, 因此, f的次数为2=max{2, 1, 1, 0}.${x^\alpha } = x_1^{{\alpha _1}}...x_n^{{\alpha _n}}, d = |\alpha | = \sum\nolimits_{i = 1}^n {{\alpha _i}} {\rm{, }}$记$c(\alpha ) = \frac{{d!}}{{{a_1}!...{a_n}!}}$.根据上述记号, 一个齐d次n元多项式能够被写成以下形式.

$ f(x)=\sum\limits_{|\alpha|=d} a_{\alpha} x^{\alpha}=\sum\limits_{|\alpha|=d} c(\alpha) b_{\alpha} x^{\alpha} $

(6)

这里, a_α=c(α)·b_α, 其中, c(α)为α的函数.例如, 考虑2元齐2次多项式:

$ f\left(x_{1}, x_{2}\right)=2 x_{1}^{2}-6 x_{1} x_{2}=2 x^{(2, 0)}-6 x^{(1, 1)}=c((2, 0)) \cdot b_{(2, 0)} x^{(2, 0)}+c((1, 1)) \cdot b_{(1, 1)} x^{(1, 1)}, $

其中, c((2, 0))=1, b_{(2, 0)}=2, c((1, 1))=2, b_{(1, 1)}=-3.

定义3.给定n元齐次(或非齐次)多项式f(x₁, …, x_n)∈R[x]和S $ \subseteq $Rⁿ.f在S上是正定的, 如果对任意的x∈S\{0}, 都有f(x) > 0;同时, f在S上是半正定的, 如果对任意的x∈S, 有f(x)≥0.

显然, 若f在S上是正定的, 那么f在S上必然是半正定的.记单形${\mathit{\Delta} _n} = \left\{ {({x_1}, ..., {x_n}):{x_i} \geqslant 0, \sum\nolimits_{i = 1}^n {{x_i}} = 1} \right\}{\rm{.}}$

定理2(Polya定理).给定n元齐次多项式f(x)∈R[x], 如果f在Δ_n上是正定的, 那么当N充分大时, 表达式(x₁+…+x_n)^Nf(x₁, …, x_n)展开后的所有系数均为正.

例如, 不难验证$f({x_1}, {x_2}) = x_1^2 - {x_1}{x_2} + x_2^2$在Δ₂={(x₁, x₂):x₁≥0, x₂≥0, x₁+x₂=1}上是正定的.这是因为f(x₁, x₂)= $x_1^2 - {x_1}{x_2} + x_2^2 = \left( {\frac{{x_1^2}}{2} + \frac{{x_2^2}}{2}} \right) + {\left( {\frac{{{x_1}}}{{\sqrt 2 }} - \frac{{{x_2}}}{{\sqrt 2 }}} \right)^2}$.因为$\frac{{x_1^2}}{2} + \frac{{x_2^2}}{2}$只能在原点O(0, 0)点处的值为0, 且o(0, 0)$ \notin $Δ₂, 故对任意的$({x_1}, {x_2}) \in {\mathit{\Delta} _2}, \frac{{x_1^2}}{2} + \frac{{x_2^2}}{2} > 0$.又因为$\left( {\frac{{{x_1}}}{{\sqrt 2 }} - \frac{{{x_2}}}{{\sqrt 2 }}} \right) \geqslant 0$, 因此f在Δ₂上正定.尽管f中的系数不是全部为正, 但根据Polya定理, 存在N, 使得(x₁+x₂)^N·f的系数全为正.为此, 令N=1, 展开${({x_1} + {x_2})^1} \cdot f = x_1^3 + x_2^3$, 其系数均为正.

注:Polya定理表明, 如果给定的齐次多项式f在单形上是正定的, 那么必然存在N, 使得${\left( {\sum\nolimits_{i = 1}^n {{x_i}} } \right)^N}f$展开后各项系数为正; 反之, 若存在N, 使得${\left( {\sum\nolimits_{i = 1}^n {{x_i}} } \right)^N}f$展开后各项系数为正, 则f在单形是半正定的(正定多项式一定是半正定的).因此, Polya定理给出了f在单形上是正定的必要条件, 同时也给出了f是单形上的半正定多项式的充分条件.但Polya定理并没告诉我们N的界.在文献[21]中, Powers等人构造了Polya定理中N的界, 并证明了那样的界仅仅依赖于多项式f的次数、系数和f在单形上的最小值.

定理3^[21].给定在单形Δ_n上正定的n元齐d次多项式f(x).如果

$ N>\frac{d(d-1)}{2} \frac{L}{\lambda}-d $

(7)

那么, (x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正.其中, L=L(f)=max{|b_α|:|α|=d}, λ=λ(f)=min{f(x):x∈Δ_n}.

注:关于定理3有几点说明.

(A) 若存在正整数N, 使得(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正, 那么对任意的非负整数j, 有(x₁+…+ x_n)^N+jf(x₁, …, x_n)的所有系数为正.这是因为(x₁+…+x_n)中各项前的系数均为正数, 故(x₁+…+x_n)^j中各项前的系数必为正数.所以当(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正时, (x₁+…+x_n)^N+jf(x₁, …, x_n)=(x₁+…+ x_n)^j(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数仍然为正.

(B) N的下界公式, 即公式(7)中, N的下界值与L的值成正比而与最小值成反比.在下一节中, 我们将通过对公式(7)中L或λ的值进行放缩来调整N的下界公式, 并在新的下界公式下建立相应算法.

(C) 若存在N, 使得(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正, 那么f(x₁, …, x_n)至少是单形上的半正定多项式.这是因为所有变量取值于单形, 故均为非负变量; 同时, 每项系数都为正.

2 程序U的秩函数计算

本节中, 我们将程序U的秩函数计算问题归结为单形上的正定多项式的探测问题, 并基于Polya和Powers等人的结果, 将单形上正定多项式探测的问题归结为线性规划问题.

定义1中的秩函数定义是传统秩函数定义, 呈现在许多文献中.在定义1中, 秩函数定义的有界和下降两个条件是由有限个蕴涵式刻画, 且蕴涵式后件为非严格不等式.但为了便于利用正定多项式的性质, 我们将定义1中程序U的秩函数定义做了稍微修改, 即将定义1中蕴涵式后件中的非严格不等式全部替换为严格不等式.

定义4.给定如公式(5)定义的循环程序U, 函数ρ(x)是程序U的秩函数, 若存在一个正数c > 0, 使得ρ(x)满足:

● (有界):$ \forall \mathit{\boldsymbol{x}}$, (x∈Ω_x$ \Rightarrow $ρ(x) > 0);

● (下降): $\mathop \wedge \limits_{i = 1}^m (\forall x\forall x', (X \in {\mathit{\Omega} _x} \wedge x' = {T_i}(x) \Rightarrow \rho (x) - \rho (x') > c)). $

注:修改后的秩函数定义与传统秩函数定义是等价的.这可被简单证明如下.

首先, 若ρ(x)是满足定义1中两个条件的秩函数, 即有$ \forall x$, (x∈Ω_x$ \Rightarrow $ρ(x)≥0)且$\mathop \wedge \limits_{i = 1}^m (\forall x\forall x', x \in {\mathit{\Omega} _x} \wedge x' = {T_i}(x) \Rightarrow $ ρ(x)ρ'(x)≥c$ \wedge $c > 0), 进而有$ \forall \mathit{\boldsymbol{x}}$, (x∈Ω_x$ \Rightarrow $ρ(x)≥0 > 1)且$\mathop \wedge \limits_{i = 1}^m (\forall x\forall x', x \in {\mathit{\Omega} _x} \wedge x' = {T_i}(x) \Rightarrow \rho (x) - \rho (x') \geqslant c > c/2 \wedge $ c > 0).那么令$\rho '(x) = \rho (x) + 1, c' = \frac{c}{2}$, 则有:存在一个正数c' > 0, 使得ρ'(x)满足定义4中的有界和下降两个条件.反之, 若ρ(x)是满足定义4中两个条件的秩函数, 那么ρ(x)显然也满足定义1中的两个条件.因此, 根据定义4, 计算U的秩函数等价于在Ω_x上探测m+1个正定多项式.既然定义1与定义4等价, 故计算U的秩函数$ \Leftrightarrow $在Ω_x上探测m+1个半正定多项式$ \Leftrightarrow $在Ω_x上探测m+1个正定多项式.根据定义4, 要计算程序U的秩函数, 仅需探寻函数ρ(x), 使其满足定义4中的有界和下降条件.

● (有界):

$ \forall \mathit{\boldsymbol{x}}, \left( {\mathit{\boldsymbol{x}} \in {\mathit{\Omega} _\mathit{\boldsymbol{x}}} \Rightarrow \rho (\mathit{\boldsymbol{x}}) > 0} \right) $

(8)

● (下降):

$ \mathop \wedge \limits_{i = 1}^m \left( {\forall x, \left( {x \in {\mathit{\Omega} _x} \Rightarrow {H_i}(x) > c} \right)} \right) $

(9)

这里, H_i(x)=ρ(x)-ρ(T_i(x)), c > 0.一般地, 公式(8)和公式9)中蕴含式后件的多项式表达式未必是关于x的齐次多项式, 但我们总可以引入新的变元进行齐次化, 并能保证齐次化后的蕴含式与原蕴含式是等价的.令$\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} = (\mathit{\boldsymbol{x}}, z){\rm{, }}$记$ \mathit{\Omega}_{\mathit{\boldsymbol{\hat{x}}}}=\left\{\hat{\boldsymbol{x}} \in \mathbf{R}^{n+1}: x_{1} \geqslant 0, \dots, x_{n} \geqslant 0, z>0\right\}$, 记${0_{\mathit{\boldsymbol{\hat x}}}}$为Rⁿ⁺¹中的原点.显然, ${0_{\mathit{\boldsymbol{\hat x}}}} \notin {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}}.$

命题1.给定非齐次多项式G(x)∈R[x], 记G的次数为d, 引进新变元z, 对G进行齐次化得到G_H(x, z).记w₁= $ \forall \mathit{\boldsymbol{x}}$, (x∈Ω_x$ \Rightarrow $G(x) > 0)和${w_2} \triangleq \forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \Rightarrow {G_H}({\mathit{\boldsymbol{\hat x}}}) > 0)$, 则有w₁$ \Leftrightarrow $w₂.

证明:w₂$ \Rightarrow $w₁是显然成立的.这是因为当蕴涵式w₂成立时, 令w₂中的z=1即可得到w₁也成立.下证w₁$ \Rightarrow $w₂.

若蕴涵式w₁成立, 即$ \forall \mathit{\boldsymbol{x}}$, (x∈Ω_x$ \Rightarrow $G(x) > 0).又${G_H}({\hat x}) = {G_H}(x, z) = {z^d} \cdot G\left( {\frac{{{x_1}}}{z}, ..., \frac{{{x_n}}}{z}} \right)$, 令${t_1} = \frac{{{x_1}}}{z}, ..., {t_n} = \frac{{{x_n}}}{z}$.若z > 0, 则(x₁, …, x_n)∈Ω_x等价于t₁≥0, …, t_n≥0.既然已知w₁成立, 因此有t₁≥$0 \wedge $…$ \wedge $t_n≥0$ \Rightarrow $G(t₁, …, t_n) > 0.又因为z > 0, 故

$ {t_1} \ge {0 \wedge }{ \ldots \wedge }{t_n} \ge 0 \wedge z > 0 \Rightarrow {z^d} \cdot G\left( {{t_1}, \ldots , {t_n}} \right) > 0 $

(10)

将${t_i} = \frac{{{x_i}}}{z}, i = 1, ..., n$代入到公式(10), 则有:

$ \frac{x_{1}}{z} \geqslant 0 \wedge \ldots _\wedge \frac{x_{n}}{z} \geqslant 0 \wedge z>0 \Rightarrow z^{d} \cdot G\left(\frac{x_{1}}{z}, \ldots, \frac{x_{n}}{z}\right)=G_{H}(x, z)>0 $

(11)

又因为$\left\{ {(x, z):\frac{{{x_1}}}{z} \geqslant 0 \wedge ... \wedge \frac{{{x_n}}}{z} \geqslant 0 \wedge z > 0} \right\} = \{ (x, z):{x_1} \geqslant 0 \wedge ... \wedge {x_n} \geqslant 0 \wedge z > 0\} $, 将公式(11)中的前件替换为x₁≥0$ \wedge $…$ \wedge $x_n≥0$ \wedge $z > 0, 即可得w₂成立.

假设公式(8)和公式(9)中的多项式ρ(x), H_i(x)-c的次数分别为$ d_{\rho}, d_{H_{i}}$, 然后引入新变元z, 对ρ(x), H_i(x)-c, i=1, …, m进行齐次化.记齐次化所得的表达式分别为

$ \hat{\rho}(\hat{\boldsymbol{x}})=\hat{\rho}(\boldsymbol{x}, z), M_{i}(\hat{\boldsymbol{x}})=M_{i}(\boldsymbol{x}, z)=\hat{H}_{i}(\hat{\boldsymbol{x}})-c \cdot z^{d_{H_{i}}}. $

显然, ${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$的次数分别为${{d}_{\rho }},{{d}_{{{H}_{i}}}} $, 其中, ${\hat H_i}({\mathit{\boldsymbol{\hat x}}})$是将H_i(x)进行齐次化后所得到的表达式.根据命题1, 定义4中刻画秩函数的两个条件:公式(8)和公式(9)分别等价于以下公式.

$ \left. \begin{array}{l} {\forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \Rightarrow {\hat \rho } (\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) > 0} \right)}\\ {\mathop \wedge \limits_{i = 1}^m \left( {\forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \Rightarrow {M_i}(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) > 0} \right)} \right.} \end{array} \right\} $

(12)

令${\mathit{\bar \Omega} _{\mathit{\boldsymbol{\hat x}}}} = \{ {\mathit{\boldsymbol{\hat x}}} = (x, z) \in {{\mathbf{R}}^{n + 1}}:{x_1} \geqslant 0, ..., {x_n} \geqslant 0, z \geqslant 0\} {\rm{, }}$既然${0_{\mathit{\boldsymbol{\hat x}}}} \notin {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}}$且${\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \subseteq {\mathit{\bar \Omega} _{\mathit{\boldsymbol{\hat x}}}}$, 显然有${\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \subseteq {\mathit{\bar \Omega} _{\mathit{\boldsymbol{\hat x}}}}\backslash \{ {0_{\mathit{\boldsymbol{\hat x}}}}\} $.

令$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o = {\mathit{\bar \Omega} _{\mathit{\boldsymbol{\hat x}}}}\backslash \{ {0_{\mathit{\boldsymbol{\hat x}}}}\} $, 即$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中不包含原点, 因此, 将公式(12)中的${\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}}$替换为$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 则得到:

$\left. \begin{gathered} \forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) > 0) \\ \mathop \wedge \limits_{i = 1}^m (\forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow {M_i}({\mathit{\boldsymbol{\hat x}}}) > 0)) \\ \end{gathered} \right\}$

(13)

既然${\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}} \subseteq \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 因此, 若公式(13)成立, 则公式(12)也成立.因此, 计算满足公式(8)和公式(9)的函数ρ(x)可归结为判定公式(13)是否成立.

令Δ_n+1={(x₁, …, x_n, z)∈Rⁿ⁺¹:x₁≥0, …, x_n≥0, z≥0, x₁+…+x_n+z=1}.

显然, ${0_{\mathit{\boldsymbol{\hat x}}}} \notin {\mathit{\Delta} _{n + 1}}$且${\mathit{\Delta} _{n + 1}} \subseteq {\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}}$.因此, ${\mathit{\Delta} _{n + 1}} \subseteq \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$.令

$\left. \begin{gathered} \forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}} \Rightarrow {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) > 0) \\ \mathop \wedge \limits_{i = 1}^m (\forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}} \Rightarrow {M_i}({\mathit{\boldsymbol{\hat x}}}) > 0)) \\ \end{gathered} \right\}$

(14)

给定两个非零点${\hat x}, \hat y \in {{\mathbf{R}}^{n + 1}}$, 我们记$\hat y \sim {{\hat x}}$.如果存在正数λ > 0, 使得$ {\mathit{\boldsymbol{\hat y}}} = \lambda \cdot {\mathit{\boldsymbol{\hat x}}}$, 令$[{\mathit{\boldsymbol{\hat x}}}] = \{ {\mathit{\boldsymbol{\hat y}}} \in {{\mathbf{R}}^{n + 1}}:{\mathit{\boldsymbol{\hat y}}} \sim {\mathit{\boldsymbol{\hat x}}}\} {\rm{, }}$即$[{\mathit{\boldsymbol{\hat x}}}]$代表了一条从原点出发且通过点${\mathit{\boldsymbol{\hat x}}}$的射线, 而${\mathit{\boldsymbol{\hat x}}}$被称为该射线的代表元.显然, 射线上任一点均可被称为该射线的代表元, 记$\sum {({\mathit{\boldsymbol{\hat x}}})} = {x_1} + ... + {x_n} + z{\rm{.}}$下面的结论表明, 公式(13)与公式(14)是等价的.

命题2.记号同上.公式(13) $ \Leftrightarrow $公式(14).

证明:既然${\mathit{\Delta} _{n + 1}} \subseteq \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 那么公式(13) $ \Rightarrow $公式(14)是平凡成立的.下面证公式(14)$ \Rightarrow $公式(13)也成立.根据$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$的定义可知, $\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$是由从原点出发但不包含原点的射线构成, 因此, 给定${\mathit{\boldsymbol{\hat x}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 对任意的λ > 0, 都有:$\lambda \cdot {\mathit{\boldsymbol{\hat x}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \cdot {\mathit{\Delta} _{n + 1}}$是$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$与平面(x₁+…+x_n+z=1的交集.下面证明$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中的每一条射线均与x₁+…+x_n+z=1定义的平面相交.

任取一点${{\mathit{\boldsymbol{\hat x}}}_0} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 根据$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$的定义, ${{\mathit{\boldsymbol{\hat x}}}_0} \ne 0$.取${\lambda _0} = \frac{1}{{\sum {{{\mathit{\boldsymbol{\hat x}}}_0}} }}$.既然$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中的任一非零点的所有分量均非负且不全为0, 故λ₀ > 0.因此, ${\lambda _0}{{{\mathit{\boldsymbol{\hat x}}}_0}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$.同时有$\sum {({\lambda _0}{{{\mathit{\boldsymbol{\hat x}}}}_0})} = {\lambda _0}\sum {({{\mathit{\boldsymbol{\hat x}}}_0})} = 1$, 故${\lambda _0}{{{\mathit{\boldsymbol{\hat x}}}_0}}$落在x₁+…+x_n+z+…=1定义的平面上.由此可得, ${\lambda _0}{{{\mathit{\boldsymbol{\hat x}}}_0}} \in {\mathit{\Delta} _{n + 1}}$即为交点.这就证明$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中的每一条射线均与x₁+…+x_n+z=1定义的平面相交.其交点与从原点出发且过该交点的射线一一对应.因此, 这样的交点实为与之对映的射线的代表元.为方便, 称它们为标准代表元.因此, $\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中的任何一点${\mathit{\boldsymbol{\hat x}}}$都在Δ_n+1上存在唯一一点$\frac{{{\mathit{\boldsymbol{\hat x}}}}}{{\sum {{\mathit{\boldsymbol{\hat x}}}} }}$与之对应.故, $\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$中所有射线对映的标准代表元构成了Δ_n+1.给定齐d次多项式$U({\mathit{\boldsymbol{\hat x}}})$, 它具有如下简单性质:若对${\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}}$有$U({\mathit{\boldsymbol{\hat x}}}) > 0{\rm{, }}$那么对任意λ > 0, 有$U(\lambda {\mathit{\boldsymbol{\hat x}}}) = $ ${\lambda ^d}U({\mathit{\boldsymbol{\hat x}}}) > 0$.也即若齐次多项式$U({\mathit{\boldsymbol{\hat x}}})$在单形Δ_n+1上某点的取值为正, 那么该齐次多项式也必在过该点且从原点出发但不包含原点的射线上的任一点处取值为正.进而可得:如果齐次多项式$U({\mathit{\boldsymbol{\hat x}}})$在单形Δ_n+1上的每一点处取值为正, 则根据上述分析得知, $U({\mathit{\boldsymbol{\hat x}}})$在$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$的每一点处为正.因此, 若公式(14)成立, 即对$\forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}} \Rightarrow {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) > 0)$且$ \wedge _{i = 1}^m(\forall {\mathit{\boldsymbol{\hat x}}}, ({\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}} \Rightarrow {M_i}({\mathit{\boldsymbol{\hat x}}}) > 0))$, 那么因为$\rho ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$为齐次多项式且其次数分别为$d_{\rho}, d_{H_{i}} $, 故有:

$ \forall \lambda \forall \hat{\boldsymbol{x}}, \left(\lambda>0 \wedge \hat{\boldsymbol{x}} \in \mathit{\Delta}_{n+1} \Rightarrow \lambda^{d_{p}} \hat{\rho}(\hat{\boldsymbol{x}})=\hat{\rho}(\lambda \hat{\boldsymbol{x}})>0\right)且\wedge_{i=0}^{m}\left(\forall \lambda \forall \hat{\boldsymbol{x}}, \left(\lambda>0 \wedge \hat{\boldsymbol{x}} \in \mathit{\Delta}_{n+1} \Rightarrow \lambda^{d_{M_{i}}} M_{i}(\hat{\boldsymbol{x}})=M_{i}(\lambda \hat{\boldsymbol{x}})>0\right).\right. $

又因为集合$\{ \lambda {\mathit{\boldsymbol{\hat x}}}:\lambda > 0 \wedge {\mathit{\boldsymbol{\hat x}}} \in {\mathit{\Delta} _{n + 1}}\} $正好就是$\mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o$, 故公式(13)成立.

注:若将公式(13)和公式(14)中的所有不等式均改为非严格不等式, 则命题2仍然成立.

要判定公式(14)是否成立, 等价于判定齐次多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$在单形Δ_n+1上是否正定的.根据上述分析, 为了计算程序U的秩函数, 需要以下5个步骤.

S1.给出秩函数模板——带参数的多项式表达式.不失一般性, 可假设为ρ(x)=a^T·Γ, 其中, a为参系数向量, Γ是所有单项式构成的向量.如ρ(x)=a₁x²y+a₂x+a₃y+a₄, 则a=(a₁, …, a₄)^T, Γ=(x²y, x, y, 1)^T.

S2.将ρ(x)代入到定义4中给出的秩函数条件中, 从而构造公式(8)和公式(9).其中, 可记公式(9)中的H_i(x)=$\rho (x) - \rho ({T_i}(x)) = \mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^T \cdot {\mathit{\Gamma '}}$.这里, b_a为一向量, 其每个分量均是关于a的齐次线性表达式; Γ'为单项式构成的向量.如, 给定循环程序的第i个赋值语句T_i(x)(x:=3x-4y², y:=y)^T, 则H_i(x)=ρ(x)-ρ(T_i(x))=-8a₁x²y-2a₂x+24a₁y³x+4a₂y²-16a₁y⁵.故, a=(-8a₁, -2a₂, 24a₁, 4a₂, -16a₁)^T, Γ'=(x²y, x, y³x, y², y⁵)^T.

S3.对公式(8)和公式(9)中的ρ(x), H_i(x)-c进行齐次化, 得到:

$ {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) = {\mathit{\boldsymbol{a}}^T} \cdot {\mathit{\Gamma} _0}, {M_i}({\mathit{\boldsymbol{\hat x}}}) = {\hat H_i}({\mathit{\boldsymbol{\hat x}}}) - c \cdot {z^{{d_{{H_i}}}}} = \mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^T \cdot {{\mathit{\Gamma '}}_0} - c \cdot {z^{{d_{{H_i}}}}} = (\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^T, - c) \cdot {({\mathit{\Gamma '}}_0^{\rm T}, {z^{{d_{{H_i}}}}})^T}. $

从而构造公式(12).这里, ${d_{{H_i}}}$为H_i(x)的次数, ${\mathit{\Gamma} _o}, {{\mathit{\Gamma '}}_o}$分别为Γ, Γ'齐次化后得到的向量.

S4.根据上述分析, 因为公式(13)$ \Rightarrow $公式(12), 所以要使得公式(12)成立, 只需保证公式(13)成立即可.再根据命题2可知, 公式(13)成立等价于公式(14)成立.进而将问题归结为保证公式(14)成立.

S5.探测是否存在参系数$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的一组取值$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$, 使得${{\hat{\rho }}^{*}}({\mathit{\boldsymbol{\hat x}}})={{a}^{*T}}\cdot {{\tau }_{0}},M_{i}^{*}({\mathit{\boldsymbol{\hat x}}})=(b_{a}^{*T},-{{c}^{*}})\cdot {{({\tau }_{0}^{'T},{{z}^{{{d}_{{{H}_{i}}}}}})}^{T}}$满足公式(14).若那样的一组取值存在, 则表明程序U是终止的.

从上面5个步骤可以看出, 计算程序U秩函数的问题被归约为探测单形上的正定多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$的问题.也即若单形上的正定多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$被探测到, 则表明程序U具有秩函数, 从而证明程序U是终止的.

根据下面的命题可知, 在上述S5中, 存在参系数$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的一组取值$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$, 使得${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}) = {a^{*T}} \cdot {\mathit{\Gamma} _0}, M_i^*({\mathit{\boldsymbol{\hat x}}}) = $ $(b_a^{*T}, - {c^*}) \cdot {({\mathit{\Gamma '}}_0^T, {z^{{d_{{H_i}}}}})^T}$满足公式(14), 等价于存在参系数$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $在${[ - 1, 1]^\ell }$中的一组取值$({\mathit{\boldsymbol{a}}^{**}},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^{**},{c^{**}}){\rm{, }}$使得${{\hat \rho } ^{**}}({\mathit{\boldsymbol{\hat x}}}) = {a^{**T}} \cdot {\mathit{\Gamma} _0}, M_i^{**}({\mathit{\boldsymbol{\hat x}}}) = (b_a^{**T}, - {c^{**}}) \cdot {({\mathit{\Gamma '}}_0^T, {z^{{d_{{H_i}}}}})^T}$满足式(14).这里, $\ell $=|$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $|为向量$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $中分量的个数. $(\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c) \in {[ - 1, 1]^\ell }$表明, 向量$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $中每个分量都在区间[-1, 1]中取值.

命题3.记号同上.存在$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c_*}){\rm{ }} \in {{\mathbf{R}}^\ell }{\rm{, }}$使得${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}) = {a^{*T}} \cdot {\tau _0}, M_i^*({\mathit{\boldsymbol{\hat x}}}) = (b_a^{*T}, - {c^*}) \cdot {(\tau _0^{'T}, {z^{{d_{{H_i}}}}})^T}$满足公式(14), 等价于存在$({\mathit{\boldsymbol{a}}^{**}},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^{**},{c^{**}}) \in {[ - 1, 1]^\ell }$, 使得${{\hat \rho } ^{**}}({\mathit{\boldsymbol{\hat x}}}) = {a^{**T}} \cdot {\tau _0}, M_i^{**}({\mathit{\boldsymbol{\hat x}}}) = (b_a^{**T}, - {c^{**}}) \cdot {(\tau _0^{'T}, {z^{{d_{{H_i}}}}})^T}$满足公式(14).

证明:该证明非常简单.“ $ \Leftarrow $”是平凡成立的, 既然${[ - 1, 1]^\ell } \subset {{{R}}^\ell }$.故我们仅证“$ \Rightarrow $”.

若存在$\nu = ({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*}) \in {{{R}}^\ell }$, 使得${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), M_i^*({\mathit{\boldsymbol{\hat x}}})$满足公式(14)中的两个蕴涵式, 那么我们可以令

$\nu ' = ({\mathit{\boldsymbol{a}}^{**}},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^{**},{c^{**}}) = \left( {\frac{{{\mathit{\boldsymbol{a}}^*}}}{{abs(v)}}, \frac{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*}{{abs(v)}}, \frac{{{c^*}}}{{abs(v)}}} \right).$

这里, $abs(v) = \sum {i|{v_i}|} $为向量v的分量的绝对值之和.显然有, $v' \in {\rm{ }}{[ - 1, 1]^l}$且${{\hat \rho } ^{**}}({\mathit{\boldsymbol{\hat x}}}) = \frac{1}{{abs(v)}}{{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}) = \frac{1}{{abs(v)}} \cdot {a^{*T}} \cdot $ ${\tau _o}, M_i^{**}({\mathit{\boldsymbol{\hat x}}}) = \frac{1}{{abs(v)}}M_i^*({\mathit{\boldsymbol{\hat x}}}) = \frac{1}{{abs(v)}} \cdot (b_a^{*T}, - {c^*}) \cdot {(\tau _o^{'T}, {z^{{d_{{H_i}}}}})^T}$满足公式(14).

根据定义4、命题1~命题3, 我们可以建立下面的结论.

定理4.给定程序U, 若存在$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$的一组取值$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*}) \in {[ - 1, 1]^\ell }$, 使得${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}) = {a^{*T}} \cdot {\tau _0}, M_i^*({\mathit{\boldsymbol{\hat x}}}) = (b_a^{*T}, - {c^*}) \cdot $ ${(\tau _0^{'T}, {z^{{d_{{H_i}}}}})^T}$满足公式(14), 那么程序U必然终止.

在S5中, 我们需要判断是否存在$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$的一组取值$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$, 使得${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), M_i^*({\mathit{\boldsymbol{\hat x}}})$满足公式(14), 这等价于在单形Δ_n+1上探测是否存在正定多项式${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), M_i^*({\mathit{\boldsymbol{\hat x}}})$.我们将利用上述的定理2(Polya定理)和定理3去探测单形上的正定多项式.定理2表明, 若齐次多项式f在单形Δ_n上是正定的, 则必然存在充分大的正整数N, 使得(x₁+…+x_n)^Nf(x₁, …, x_n)展开后的所有系数均为正.但Polya并没有给出N的上界.在定理3中, Powers等人根据f的次数、系数以及f在单形上的最小值, 构造了N的界.

$N > \frac{{d(d - 1)}}{2}\frac{L}{\lambda } - d$

(15)

这里, L=L(f)=max{|b_α|:|α|=d}, λ=λ(f)=min{f(x):x∈Δ_n}.根据下面的结论, 我们可以将公式(15)中的L替换为1.

命题4.记号同上.给定齐d次n元多项式$f(\mathit{\boldsymbol{x}}){\rm{ }} = \sum\nolimits_{|\alpha | = d} {{a_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} = \sum\nolimits_{|\alpha | = d} {c(\alpha ){b_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} , $这里, 非负整数向量α= (α₁, …, α_n)∈Nⁿ, 则(α₁+…+α_n)!≥α₁!…α_n!.

证明:考虑n=2.

$ \left(\alpha_{1}+\alpha_{2}\right) !=\left(\alpha_{1}-0+\alpha_{2}\right)\left(\alpha_{1}-1+\alpha_{2}\right)\left(\alpha_{1}-2+\alpha_{2}\right) \ldots\left(\alpha_{1}-\left(\alpha_{1}-1\right)+\alpha_{2}\right) \cdot\left(\alpha_{2}\right) ! $

显然, (α₁-0+α₂)≥α₁, (α₁-1+α₂)≥α₁-1, (α₁-2+α₂)≥α₁-2, …, (α₁-(α₁-1)+α₂)≥1.故

$ \left(\alpha_{1}+\alpha_{2}\right) ! \geqslant \alpha_{1} ! \alpha_{2} ! $

(16)

考虑n=3.

根据公式(16), 有(α₁+α₂+α₃)!=(α₃+(α₁+α₂))!≥(α₃)!(α₁+α₂)!≥α₃!α₁!α₂!.

以此类推, 可得(α₁+…+α_n)!≥α₁!…α_n!.

注:给定齐d次n元多项式$f(\mathit{\boldsymbol{x}}){\rm{ }} = \sum\nolimits_{|\alpha | = d} {{a_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} = \sum\nolimits_{|\alpha | = d} {c(\alpha ){b_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} $, 在公式(15)中,

$L = L(f) = \max \left\{ {|{b_\alpha }| = \frac{{|{a_\alpha }|}}{{c(\alpha )}}:|\alpha | = d} \right\}.$

根据命题4可知, $c(\alpha ) = \frac{{d!}}{{{\alpha _1}!...{\alpha _n}!}} = \frac{{({\alpha _1} + ... + {\alpha _n})!}}{{{\alpha _1}!...{\alpha _n}!}} \geqslant 1.$假如f的所有系数a_α均在区间[-η, η]中, 那么$\frac{{|{a_\alpha }|}}{{c(\alpha )}} \leqslant \eta $.因此, 倘若f的所有系数a_α均在区间[-η, η]中, 那么$L = L(f) = \max \left\{ {|{b_\alpha }| = \frac{{|{a_\alpha }|}}{{c(\alpha )}}:|\alpha | = d} \right\} \leqslant \eta $.根据上面的分析, 若系数a_α在区间[−1, 1], 可以得到下列结果:

定理5.给定齐d次n元多项式$f(\boldsymbol{x})=\sum_{|\alpha|=d} a_{\alpha} \boldsymbol{x}^{\alpha} $, 且其所有系数a_α均在区间[−1, 1]中.如果f在Δ_n上是正定的, 则当:

$N > \frac{{d(d - 1)}}{2}\frac{1}{\lambda } - d$

(17)

有(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正.其中, λ=λ(f)=min{f(x):x∈Δ_n.

根据定理4得知, 可以在系数空间${[ - 1, 1]^\ell }$中探测满足公式(14)的齐次多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) = {\mathit{\boldsymbol{a}}^T} \cdot {\tau _0}, {M_i}({\mathit{\boldsymbol{\hat x}}}) = (\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^T, - c) \cdot $ ${(\tau _0^{'T}, {z^{{d_{{H_i}}}}})^T}$.这里, $(\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c) \in {[ - 1, 1]^\ell }$.同时, 满足公式(14)的${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$恰好是单形Δ_n+1上的正定多项式.因此, 既然${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$在单形Δ_n+1上正定且其所有系数均在${[ - 1, 1]^\ell }$中取值, 那么根据定理5, 必然存在一个N(其值仅仅依赖于最小值λ), 使得${({x_1} + ... + {x_n} + z)^N}{\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {({x_1} + ... + {x_n} + z)^N}{M_i}({\mathit{\boldsymbol{\hat x}}})$的所有系数为正.因此, 为了得到$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的一组取值, 我们可以从${({x_1} + ... + {x_n} + z)^N}{\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {({x_1} + ... + {x_n} + z)^N}{M_i}({\mathit{\boldsymbol{\hat x}}})$中抽取出所有系数(它们均是关于$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的齐次线性表达式), 并令所有系数大于0, 且加上约束$c > 0 \wedge (\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c) \in {[ - 1, 1]^\ell }$, 从而构造出关于$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的不等式组Sys.若Sys有实解, 则它的任一组解$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$就正好给出了单形Δ_n+1上的(半)正定多项式${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), M_i^*({\mathit{\boldsymbol{\hat x}}})$(这是因为根据定理2或定理3的注解可知, 若存在N, 使得${\left( {\sum\nolimits_{i = 1}^n {{x_i}} } \right)^N}f$展开后各项系数为正, 则f在单形是半正定的).下面的定理表明, 若系统Sys有实解, 则程序U可终止.

定理6.记号同上.给定程序U, 若系统Sys有实解, 则程序U是终止的.

证明:不妨令$({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*})$为系统Sys的一组解, 则有${c^*} > 0 \wedge ({\mathit{\boldsymbol{a}}^*},\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}^*,{c^*}) \in {[ - 1, 1]^\ell };$且存在N, 使得${({x_1} + ... + {x_n} + z)^N} \cdot $ ${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), {({x_1} + ... + {x_n} + z)^N}M_i^*({\mathit{\boldsymbol{\hat x}}})$的各项系数为正.这里, $M_i^*({\mathit{\boldsymbol{\hat x}}}) = (b_a^{*T}, - {c^*}) \cdot {({\mathit{\Gamma '}}_0^T, {z^{{d_{{H_i}}}}})^T}, {{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}) = {a^{*T}} \cdot {\mathit{\Gamma} _0}.$因此, 根据定理2或定理3的注解可知, ${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}}), M_i^*({\mathit{\boldsymbol{\hat x}}})$是单形Δ_n+1上的半正定多项式, 即

$ \forall \hat{\boldsymbol{x}}, \left(\hat{\boldsymbol{x}} \in \mathit{\Delta}_{n+1} \Rightarrow \hat{\rho}^{*}(\hat{\boldsymbol{x}}) \geqslant 0\right)且\mathop \wedge \limits_{i = 1}^m \left(\forall \hat{\boldsymbol{x}}, \left(\hat{\boldsymbol{x}} \in \mathit{\Delta}_{n+1} \Rightarrow M_{i}^{*}(\hat{\boldsymbol{x}}) \geqslant 0\right)\right) $

根据命题2及其注解, 上式等价于

$ \forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow {{{\hat \rho } }^*}(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) \ge 0} \right)且\mathop \wedge \limits_{i = 1}^m \left( {\forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow M_i^*(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) \ge 0} \right)} \right). $

也即$\forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow {{{\hat \rho } }^*}(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) \ge 0} \right) $且$\mathop \wedge \limits_{i = 1}^m \left( {\forall \mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}, \left( {\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}} \in \mathit{\Omega} _{\mathit{\boldsymbol{\hat x}}}^o \Rightarrow M_i^*(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) = M_i^*(\mathit{\boldsymbol{x}}, z) = \hat H_i^*(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) - {c^*} \cdot {z^{{d_{{H_i}}}}} \ge 0} \right)} \right) $

在上式中令z=1, 得到:

$ \forall {\rm{ }}\mathit{\boldsymbol{x}}, \left( {\mathit{\boldsymbol{x}} \in {\mathit{\Omega} _\mathit{\boldsymbol{x}}} \Rightarrow {{{\hat \rho } }^*}(\mathit{\boldsymbol{{\mathit{\boldsymbol{\hat x}}}}}) \ge 0} \right)且\mathop \wedge \limits_{i = 1}^m \left( {\forall \mathit{\boldsymbol{x}}, \left( {\mathit{\boldsymbol{x}} \in {\mathit{\Omega} _\mathit{\boldsymbol{x}}} \Rightarrow {\rho ^*}(\mathit{\boldsymbol{x}}) - {\rho ^*}\left( {{T_i}(\mathit{\boldsymbol{x}})} \right) - {c^*} \ge 0} \right)} \right). $

根据定义1可知, 满足上式的ρ^*(x)是程序U的秩函数.故程序U可终止.

注:根据定理6, 若系统Sys有解, 则能够得到一个最小值为预设定值λ^*的具有预定形式的秩函数; 反之, 若系统Sys无解, 则表明在单形上没有最小值为预设定值λ^*且具有预定形式的正定多项式.

根据定理4, ${{\hat \rho } ^*}({\mathit{\boldsymbol{\hat x}}})$的存在, 表明了程序U是终止的.综上所述, N的计算是关键的.在公式(17)中, N的值仅依赖于多项式的次数d和齐次多项式在单形上的最小值λ.但由于${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$的所有系数均是关于参数$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的线性表达式, 故它们在单形Δ_n+1上的最小值是关于$\left( {\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c} \right) $的函数, 即${\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}} = \lambda (\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c), {\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}} = \lambda (\mathit{\boldsymbol{a}}{\rm{,}}\mathit{\boldsymbol{ }}{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}c)$.因此在实际的计算中, 我们需事先给定${\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}}, {\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}}$的值, 然后根据公式(17)计算得到N的值并抽取${({x_1} + ... + {x_n} + z)^N}{\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), $ ${({x_1} + ... + {x_n} + z)^N}{M_i}({\mathit{\boldsymbol{\hat x}}})$的所有系数构造上述的不等式组Sys.再根据定理6, 若Sys有实解, 则获得一个满足预定模板形式的秩函数; 但若该不等式组没有实解, 则需要再次设定新的最小值${\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}}, {\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}}$, 并使其值小于上一次设定的最小值(这是因为根据公式(17)可知, 当最小值越小时, N所需的下界值越大).因此, 上述秩函数的探测过程是试探性(heuristic)的.当然, 为了保证探测过程终止, 可人为固定探测的深度depth.根据上述过程, 我们建立下列算法1.

算法1.

输入:一个程序U, 探测深度depth, ${\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}} = {\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}} = \lambda > 0, $变元个数n, 次数d.

输出:一个具有预定形式的n元d次多项式秩函数.

Step 1:设定多项式秩函数模板$ \rho (\mathit{\boldsymbol{x}}) = {\mathit{\boldsymbol{a}}^T} \cdot \mathit{\Gamma} $

Step 2:构造齐次多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) = {\mathit{\boldsymbol{a}}^T} \cdot {\mathit{\Gamma} _0}, {M_i}({\mathit{\boldsymbol{\hat x}}}) = {\hat H_i}({\mathit{\boldsymbol{\hat x}}}) - c \cdot {z^{{d_{{H_i}}}}}$

Step 3: For i=1 to depth

  Step 3.1:根据${\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}}, {\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}}$的值以及公式(17)分别计算${N_{{\hat \rho } }}, {N_{{M_i}}}$

  Step 3.2:抽取${({x_1} + ... + {x_n} + z)^{{N_{{\hat \rho } }}}}{\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {({x_1} + ... + {x_n} + z)^{{N_{{M_i}}}}}{M_i}({\mathit{\boldsymbol{\hat x}}})$的所有系数构造不等式组Sys

  Step 3.3:用线性规划工具Simplex计算Sys:若Sys有解, 则输出多项式秩函数ρ^*(x); 否则, 转Step 3.4

  Step 3.4:令$\lambda : = \frac{\lambda }{2};{\lambda _{{\hat \rho } ({\mathit{\boldsymbol{\hat x}}})}}: = \lambda ;{\lambda _{{M_i}({\mathit{\boldsymbol{\hat x}}})}}: = \lambda $(串行赋值).

下面通过一个例子来阐述本文的方法.

例1:考虑循环程序:

$ \left. {{\begin{array}{*{20}{l}} {{U_1}{\rm{ }}\;{\mathbf{while}}\;{\rm{ }}x \ge 0 \wedge y \ge 0\;{\mathbf{do}}}\\ {{\tau _1}:x: = 1 - 3x - 4y - {x^3};}\\ {y: = - x - y;}\\ {{\rm{ or }}}\\ {{\tau _2}:x: = - y - 1;}\\ {y: = - 7{x^2} + y;} \end{array}}} \right\} $

(18)

记T₁=(1-3x-4y-x³, -x-y), T₂=(-y-1, -7x²+y).

(a) 设定秩函数模板ρ(x)=a₁x³+a₂x²y+a₃xy²+a₄x²+a₅xy+a₆y²+a₇x+a₈y+a₉, 令

$ H_{1}(\boldsymbol{x})=\rho(\boldsymbol{x})-\rho\left(T_{1}(\boldsymbol{x})\right), H_{2}(\boldsymbol{x})=\rho(\boldsymbol{x})-\rho\left(T_{2}(\boldsymbol{x})\right). $

(b) 对H_i(x)-c进行齐次化, 得到${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_1}({\mathit{\boldsymbol{\hat x}}}), {M_2}({\mathit{\boldsymbol{\hat x}}})$, 其中,

${\hat \rho } (x, y, z) = {a_1}{x^3} + {a_2}{x^2}y + {a_3}x{y^2} + {a_4}{x^2}z + {a_5}xyz + {a_6}{y^2}z + {a_7}x{z^2} + {a_8}y{z^2} + {a_9}{z^3}.$

由于${M_1}({\mathit{\boldsymbol{\hat x}}}), {M_2}({\mathit{\boldsymbol{\hat x}}})$表达式过长, 在此省略.

(c) 判定是否存在(a₁, ..., a₉)的一组取值$(a_1^*, ..., a_9^*)$, 使得${\hat \rho } (x, y, z)$和M₁(x, y, z), M₂(x, y, z)在单形Δ₂₊₁上都同时正定.要计算得到那样的一组(a₁, ..., a₉)的取值, 根据定理5, 需首先计算出N的值.根据公式(17), N依赖于${\hat \rho } , {M_1}, {M_2}$各自次数以及各自在Δ₃上的最小值λ.但这里, ${\hat \rho } (x, y, z)$和M₁(x, y, z), M₂(x, y, z)均是参系数齐次多项式, 故其在单形上的最小值是随着参数取值的变动而变动的.因此, 在N的实际计算中, 我们需要事先固定λ的值, 比如可令${\lambda _{{\hat \rho } }} = {\lambda _{{M_1}}} = {\lambda _{{M_2}}} = \lambda = \frac{1}{3}, $然后, 既然该例中${\hat \rho } (x, y, z)$和M₁(x, y, z), M₂(x, y, z)的次数分别为3, 9, 5, 那么根据公式(17), 分别计算对应的N的取值范围为

${N_{{\hat \rho } }} > N\left( {{\hat \rho } , \lambda = \frac{1}{3}} \right) = 6, {\rm{ }}{N_{{M_1}}} > N\left( {{M_1}, \lambda = \frac{1}{3}} \right) = 99, {\rm{ }}{N_{{M_2}}} > N\left( {{M_2}, \lambda = \frac{1}{3}} \right) = 25.$

故取${N_{{\hat \rho } }} = 7, {N_{{M_1}}} = 100, {N_{{M_2}}} = 26$(当然, 我们也可以令${N_{{\hat \rho } }} = {N_{{M_1}}} = {N_{{M_2}}} = 100$, 即取公共界).根据定理5, 分别展开多项式${(x + y + z)^{{N_{{\hat \rho } }}}}{\hat \rho } (x, y, z), {(x + y + z)^{{N_{{M_1}}}}}{M_1}(x, y, z), {(x + y + z)^{{N_{{M_2}}}}}{M_2}(x, y, z), $合并同类项后, 各自抽取出关于x, y, z项的所有系数(均为关于a₁, ..., a₉, c的线性表达式), 然后再令所有系数大于0, 构造出关于a₁, ..., a₉, c的严格不等式组(为参系数的第1个约束系统), 分别记为$co{e_{{\hat \rho } }}, co{e_{{M_1}}}, co{e_{{M_2}}}.$同时, 注意定理5成立的前提条件是, ${\hat \rho } , {M_1}, {M_2}$关于x, y, z的所有系数都被限定在[-1, 1]时, 公式(17)才成立——即N的值才仅与其次数d和其在单形上的最小值λ相关.因此, 再次分别从多项式${\hat \rho } , {M_1}, {M_2}$中提取其关于x, y, z的所有系数(均为关于a₁, ..., a₉, c的线性表达式), 然后分别令每个系数≥-1且≤1, 由此构造出${\hat \rho } , {M_1}, {M_2}$关于x, y, z的系数的第2个约束系统, 分别记为${\mathit{\Theta} _{{\hat \rho } }}, {\mathit{\Theta} _{{M_1}}}, {\mathit{\Theta} _{{M_2}}}.$比如, 在该例中, 我们有:

$ \begin{gathered} {\mathit{\Theta} _{{\hat \rho } }}: = {a_1} \geqslant - 1 \wedge {a_1} \leqslant 1 \wedge {a_2} \geqslant - 1 \wedge {a_2} \leqslant 1 \wedge {a_3} \geqslant - 1 \wedge {a_3} \leqslant 1 \wedge {a_4} \geqslant - 1 \wedge {a_4} \leqslant 1 \wedge {a_5} \geqslant - 1 \wedge {a_5} \leqslant 1 \wedge \\ {\rm{ }}{a_6} \geqslant - 1 \wedge {a_6} \leqslant 1 \wedge {a_7} \geqslant - 1 \wedge {a_7} \leqslant 1 \wedge {a_8} \geqslant - 1 \wedge {a_8} \leqslant 1 \wedge {a_9} \geqslant - 1 \wedge {a_9} \leqslant 1. \\ \end{gathered} $

同时, 根据秩函数的定义4可知, c > 0.最后, 求解不等式组$Sys: = co{e_{{\hat \rho } }} \wedge co{e_{{M_1}}} \wedge co{e_{{M_2}}} \wedge c > 0 \wedge {\mathit{\Theta} _{{\hat \rho } }} \wedge {\mathit{\Theta} _{{M_1}}} \wedge {\mathit{\Theta} _{{M_2}}}$.若Sys有解, 则表明该程序具有预定形式的3次秩函数.根据定理4可知, 该程序是终止的.但是因为Maple自带的线性规划工具Simplex仅能够求解非严格的不等式组——不等式组中的所有不等式都是非严格的, 所以为了使用工具Simplex, 在实际计算中, 我们让${(x + y + z)^{{N_{{\hat \rho } }}}}{\hat \rho } (x, y, z), {(x + y + z)^{{N_{{M_1}}}}}{M_1}(x, y, z), {(x + y + z)^{{N_{{M_2}}}}}{M_2}(x, y, z), $合并同类项后, 各自抽取出关于x, y, z项的所有系数, 然后再令所有系数≥某个正数δ.这相当于对$co{e_{{\hat \rho } }}, co{e_{{M_1}}}, $ $co{e_{{M_2}}}$中的所有不等式做了一个小小的扰动, 即将 > 0替换为≥δ(δ > 0).记扰动后的系数集为${\widetilde {coe}_{{\hat \rho } }}, {\widetilde {coe}_{{M_1}}}, {\widetilde {coe}_{{M_2}}}$.同时, 将c > 0扰动为c≥δ.显然, 扰动后的系统$\widetilde {Sys}: = {\widetilde {coe}_{{\hat \rho } }} \wedge {\widetilde {coe}_{{M_1}}} \wedge {\widetilde {coe}_{{M_2}}} \wedge c \geqslant \delta \wedge {\mathit{\Theta} _{{\hat \rho } }} \wedge {\mathit{\Theta} _{{M_1}}} \wedge {\mathit{\Theta} _{{M_2}}}$为非严格不等式组, 如果它有解, 则原系统Sys也必有解.为方便, 我们让Sys, $\widetilde {Sys}$不仅表示不等式系统, 而且还表示不等式系统所对应的解集.在本例中, 我们取$\delta = \frac{1}{{1000}}$.利用Maple中的线性规划工具包Simplex, 求解不等式组$\widetilde {Sys}$, 得到:

$\left\{ \begin{gathered} {a_1} = \frac{{29}}{{7000}}, {a_2} = \frac{{51}}{{12250}}, {a_3} = \frac{1}{{980}}, {a_4} = - \frac{{43}}{{7000}}, {a_5} = 0, {a_6} = \frac{1}{{1000}}, {a_7} = 0, \\ {a_8} = \frac{{955875294522846847882517}}{{154330444367487258828048000}}, {a_9} = \frac{9}{{875}}, c = \frac{1}{{1000}} \\ \end{gathered} \right\}$

(19)

由于采用单纯形算法, 故该点是满足$\widetilde {Sys}$的一个精确点而非浮点向量.因此, 我们可以得到程序U₁的一个3次秩函数$\rho (x, y) = \frac{{29{x^3}}}{{7000}} + \frac{{51{x^2}y}}{{12250}} + \frac{{x{y^2}}}{{980}} - \frac{{43{x^2}}}{{7000}} + \frac{{{y^2}}}{{1000}} + \frac{{955875294522846847882517}}{{154330444367487258828048000}}y + \frac{9}{{875}}.$由于单纯形算法给出的是满足系统$\widetilde {Sys}( \subseteq Sys)$的精确解而非数值解, 故求得的函数必然是秩函数, 因而不必再对其验证是否为秩函数.该秩函数的存在表明, 该循环程序是终止的.

注:对于该程序, 文献[3, 5, 6, 20]中的方法已不能计算其预定形式的3次秩函数.因为这些方法仅能计算线性循环程序的线性秩函数, 而本例中的循环是非线性的, 且所要计算的秩函数也是3次的.此外, 对该程序, 我们也尝试利用一些基于量词消去技术^[22]的工具, 如Redlog、RegularChains, 去计算该程序的预定形式的3次秩函数(其中, RegularChains集成了当前功能强大的实解分类工具DISCOVERER^[23], 作者此前的诸多工作也是在该工具的支持下予以开展).这等价于利用这些量词消去工具从秩函数定义4中的两个蕴含公式(8)和公式(9)中消去变元x, y即可.但由于量词消去算法的双指数复杂度以及所处理的系统都是非线性的, 从而使得这两个工具均无法计算得到该程序的3次秩函数.同时, 我们也尝试用Cousot在文献[9]中的方法——半正定规划(SDP), 并借助半正定规划求解器YALMIP^[24]计算得到一个函数:

$ 42.87498608+4.1519 x+25.9747 y+1.5329 x^{3}-0.8212 y^{2}-8.7046 x^{2}+3.9678 x y-1.8358 x^{2} y+1.8393 x y^{2}. $

但经过检验发现, 该函数并不满足秩函数定义中的有界条件, 因此它不是程序的秩函数.

例2:考虑循环程序:

$\begin{array}{*{20}{l}} {{U_2}{\rm{ }}~{\mathbf{while}}~{\rm{ }}x \geqslant 0 \wedge y \geqslant 0{\rm{ }}~{\mathbf{do}}} \\ {{\tau _1}:x: = {\rm{ }}1 + x - y;} \\ {y: = 2x + y + 3;} \\ \begin{gathered} {\rm{or}} \\ {\tau _2}:x: = {\rm{ }}x - y; \\ y: = y + 2; \\ \end{gathered} \end{array}$

通过本文方法, 我们得到该循环的一个3次秩函数:

$\begin{gathered} \rho (x, y) = \frac{{56978573}}{{514034040000}}{x^3} - \frac{{10604929}}{{257017020000}}{x^2}y + \frac{{1155787}}{{79082160000}}x{y^2} - \frac{{16927733}}{{171344680000}}xy + \\ {\rm{ }}\frac{1}{{52000}}{y^2} + \frac{{7408419}}{{26360720000}}x - \frac{7}{{13000}}y + \frac{{203}}{{26000}}{\rm{.}} \\ \end{gathered} $

例3:考虑循环程序:

$\left. {\begin{array}{*{20}{l}} {{U_3}{\rm{ }}~{\mathbf{while}}~{\rm{ }}x \geqslant 0 \wedge y \geqslant 0{\rm{ }}~{\mathbf{do}}} \\ {{\tau _1}:x: = - x - 3y + 2;} \\ {y: = - 4{x^2} + y - 1;} \\ \begin{gathered} {\rm{or}} \\ {\tau _2}:x: = x - 1; \\ y: = - {y^2} + x - 3; \\ \end{gathered} \end{array}} \right\}$

(20)

通过本文方法, 我们得到该循环的一个3次秩函数$\rho (x, y) = \frac{1}{{4600}}{x^3} + \frac{{63}}{{23000}}y.$

基于Polya定理, 算法1提供了一个试探性方法去探测程序U的多项式秩函数.下面我们将证明, 给定齐d次n元多项式$f(x){\rm{ }} = \sum\nolimits_{|\alpha | = d} {{a_\alpha }{x^\alpha }} = \sum\nolimits_{|\alpha | = d} {c(\alpha ){b_\alpha }{x^\alpha }} , $若

(A) 所有系数均为整数, 即a_α∈Z,

(B) 其系数的绝对值被界定, 即存在正数η, 使得|a_α|≤η,

则Polya定理中关于N的下界公式可以被改写为

$N > \frac{{d(d - 1)}}{2}\frac{\eta }{{\sigma (n, d, \eta )}} - d$

(21)

显然, 上述N的下界公式仅依赖于多项式f(x)的次数、变元个数以及系数的绝对值上界.根据公式(21), 我们将建立不同于算法1的新算法.在证明上述结论之前, 我们首先引入文献[25]中建立的关于整系数正定多项式在单形上的最小值下界的一个重要结果.

定理7^[25].给定d次n元齐次整系数多项式$f(\mathit{\boldsymbol{x}}){\rm{ }} = \sum\nolimits_{|\alpha | = d} {{a_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} = \sum\nolimits_{|\alpha | = d} {c(\alpha ){b_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} \in Z[\mathit{\boldsymbol{x}}]$, 其系数均被正数η界定, 即|a_α|≤η.记${\mathit{\Delta} _n} = \{ ({x_1}, ..., {x_n}):{x_i} \geqslant 0, \sum\nolimits_{i = 1}^n {{x_i} = 1} \} $为单形.如果f(x)在单形Δ_n上是正定的, 那么有:

$\mathop {\min }\limits_{{\mathit{\Delta} _n}} f(\mathit{\boldsymbol{x}}) \geqslant {(2\eta )^{ - {d^n}}}{n^{ - {d^{n + 1}} - d}}{d^{ - n{d^n}}}$

(22)

注:因本文仅考虑单形上的正定多项式, 故定理7只引述了文献[25]中引理3.3中的前部分结论, 其后部分所涉及到单形上的负定多项式的内容在此被舍去.

根据定理2、定理3及其定理7, 我们建立下列定理8.

定理8.给定d次n元齐次整系数多项式$f(\mathit{\boldsymbol{x}}){\rm{ }} = \sum\nolimits_{|\alpha | = d} {{a_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} = \sum\nolimits_{|\alpha | = d} {c(\alpha ){b_\alpha }{\mathit{\boldsymbol{x}}^\alpha }} \in Z[x]$, 且其系数均被正数η界定, 即|a_α|≤η.如果f在Δ_n上是正定的, 那么当

$N > \frac{{d(d - 1)}}{2}\frac{\eta }{{\sigma (n, d, \eta )}} - d$

(23)

时, 有(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正.这里, $\sigma (n, d, \eta ) = {(2\eta )^{ - {d^n}}}{n^{ - {d^{n + 1}} - d}}{d^{ - n{d^n}}}.$

证明:首先, 根据定理3可知, 如果n元齐d次多项式f(x)在单形上正定, 那么存在正整数N, 当其满足:

$N > \frac{{d(d - 1)}}{2}\frac{L}{\lambda } - d$

(24)

时, 有(x₁+…+x_n)^Nf(x₁, …, x_n)的所有系数为正.其中, L=L(f)=max{|b_α|:|α|=d}, λ=λ(f)=min{f(x):x∈Δ_n}.再根据命题4及其注解可知, 倘若f的所有系数a_α均在区间[-η, η]中, 那么,

$L = L(f) = \max \left\{ {|{b_\alpha }| = \frac{{|{a_\alpha }|}}{{c(\alpha )}}:|\alpha | = d} \right\} \leqslant \eta $

(25)

同时, 根据已知题设, 既然f在单形上正定且其所有系数均为整数并均被正数η界定, 那么由定理7可得, f在单形Δ_n上的最小值:

${\lambda _f} = \mathop {\min }\limits_{{\mathit{\Delta} _n}} f(x) \geqslant {(2\eta )^{ - {d^n}}}{n^{ - {d^{n + 1}} - d}}{d^{ - n{d^n}}}$

(26)

根据公式(25)和公式(26)将公式(24)进行放缩, 即得到公式(23).

注:定理8中, 公式(23)中N的下界公式仅仅依赖多项式的次数、变元个数以及系数绝对值上界.

给定带参数的多项式$G({y_1}, ..., {y_n}) = \sum\nolimits_{|\alpha | = d} {{g_\alpha }(\mathit{\boldsymbol{v}}){\mathit{\boldsymbol{y}}^\alpha }} $ (这里, v为参数向量, g_α(v)为参数v的齐次整系数线性函数), 倘若我们限定其系数g_α(v)在某个区间取值, 即-η≤g_α(v)≤η, 那么根据定理8, 完成下列两个步骤便可以构造出使得多项式G在单形Δ_n上正定且其系数满足|g_α(v)|≤η的参数v应满足的必要条件$Sys_G^{_{\mathbf{z}}}$(与之前的Sys以示区别).

(ⅰ)根据公式(23)计算N的值, 记为N_G;

(ⅱ)抽取${({y_1} + ... + {y_n})^{{N_G}}}G({y_1}, ..., {y_n})$的所有系数g_α(v), 分别令g_α(v) > 0, 且|g_α(v)|≤η, 构造不等式组$Sys_G^{_Z}.$

这里, g_α(v)为齐次整系数线性函数这一限定具有一般性.这是因为, 若多项式G的系数均是有理数, 那么总可以通过乘上所有有理数的分母将其系数变为整数.根据公式(5), 程序U中的所有系数均为有理数, 故参数模板多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) = {\mathit{\boldsymbol{a}}^T} \cdot {\mathit{\Gamma} _0}, {M_i}({\mathit{\boldsymbol{\hat x}}}) = \sum\nolimits_{|\alpha | = {d_{{M_i}}}} {{g_\alpha }({{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}\mathit{\boldsymbol{c}}}){{{\mathit{\boldsymbol{\hat x}}}}^\alpha }} $的所有系数均为有理数.因此, 可以将所有有理数系数的分母(若有负号, 则将负号放到分子上)都乘在一起, 记为β.则在公式(14)中的不等式两端同时乘上β并不会改变不等式的符号, 故有$\beta \cdot {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), \beta \cdot {M_i}({\mathit{\boldsymbol{\hat x}}})$的系数均为整数.也即β·g_α(${{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}\mathit{\boldsymbol{c}}} $)是关于参数${{\mathit{\boldsymbol{b}}_\mathit{\boldsymbol{a}}}{\rm{, }}\mathit{\boldsymbol{c}}} $的齐次整系数线性函数.因此, 若${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$含有有理系数, 则可以乘上一个正数β, 使得其所有系数为整数.根据公式(23), 分别计算${N_{{\hat \rho } }}, {N_{{M_i}}}$并抽取$\left(x_{1}+\ldots+x_{n}+z\right)^{N} \hat{\rho}({\mathit{\boldsymbol{\hat x}}}),\left(x_{1}+\ldots+x_{n}+z\right)^{N_{M_{i}}} M_{i}({\mathit{\boldsymbol{\hat x}}}) $的所有系数分别构造不等式组$Sys_{{\hat \rho } }^{_{\mathbf{z}}}, Sys_{{M_i}}^{\mathbf{z}}.$令$Sy{s_{_{\mathbf{z}}}} = $ $Sys_{{\hat \rho } }^{_{\mathbf{z}}} \wedge ( \wedge _{i = 1}^mSys_{{M_i}}^{\mathbf{z}})$.既然所有参数均被设定为整数, 故需要在整数环上求解系统$Sy{s_Z}.$类似定理6, 下面的定理9表明, 若$Sy{s_\mathbb{Z}}$有整数解, 则程序U是终止的.

定理9.记号同上.给定程序U.若系统$ S y s_{\mathbf{Z}}$有整数解, 那么程序U有预定形式的且系数绝对值上界为η的秩函数.

证明:该证明完全类似于定理6的证明.在此省略.

注:定理9表明, 若系统$S y s_{\mathbf{Z}}$有整数解, 则程序U有一个整系数多项式秩函数; 反之, 如果系统$S y s_{\mathbf{Z}}$没有整数解, 则表明在单形上没有预定形式的且系数在[-η, η]的正定多项式.此时需要增大上界η的值继续构造新的系统并求解.因此, 整个过程仍是试探性的.

算法2.

输入:一个程序U, 变元个数n, 次数d, 系数绝对值上界η.

输出:一个具有预定形式的n元d次整系数多项式秩函数; 不确定(unknown).

Step 1:设定多项式秩函数模板$ \rho (\mathit{\boldsymbol{x}}) = {\mathit{\boldsymbol{a}}^T} \cdot \mathit{\Gamma} $

Step 2:构造齐次多项式${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}) = {\mathit{\boldsymbol{a}}^T} \cdot {\mathit{\Gamma} _0}, {M_i}({\mathit{\boldsymbol{\hat x}}}) = {\hat H_i}({\mathit{\boldsymbol{\hat x}}}) - c \cdot {z^{{d_{{H_i}}}}}$

Step 3:如果${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$含有有理系数, 则将所有有理数系数的分母都乘在一起, 记为β.令${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}): = \beta \cdot {\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), $ ${M_i}({\mathit{\boldsymbol{\hat x}}}): = \beta \cdot {M_i}({\mathit{\boldsymbol{\hat x}}}) $, 转Step 4;

Step 4:根据公式(23)分别计算${N_{{\hat \rho } }}, {N_{{M_i}}}$

Step 5:抽取$ \left(x_{1}+\ldots+x_{n}+z\right)^{N_{\hat{\rho}}} \hat{\rho}(\hat{\boldsymbol{x}}),\left(x_{1}+\ldots+x_{n}+z\right)^{N_{M_{i}}} M_{i}(\hat{\boldsymbol{x}})$的所有系数构造分别构造不等式组:

$Sys_{{\hat \rho } }^{_\mathbf{Z}}, Sys_{{M_i}}^{\mathbf{z}}$

Step 6:利用整数线性规划算法计算$Sy{s_Z} = Sys_{{\hat \rho } }^Z \wedge ( \wedge _{i = 1}^mSys_{{M_i}}^Z)$中的一个整数解.若$S y s_{\mathbf{Z}}$有整数解, 则输出整系数多项式秩函数ρ^*(x); 否则, 输出unknown

注:算法2中, 我们仅是简单地使${\hat \rho } ({\mathit{\boldsymbol{\hat x}}}), {M_i}({\mathit{\boldsymbol{\hat x}}})$的所有系数绝对值具有相同的上界η.实际上, 可以对它们设置不同的系数绝对值上界.在判定$S y s_{\mathbf{Z}}$是否有整数解时, 可以利用Pugh在文献[26]中提出的Omega test方法.

3 总结

针对一类多项式循环程序, 本文给出了一种新的方法去计算这类程序的多项式秩函数.该方法将这类程序的秩函数计算归结为单形上的正定多项式的探测问题; 然后, 利用Polya定理, 将单形上的正定多项式探测问题归结为线性不等式约束系统的可行问题.从这一角度看, 我们的方法是一个“线性化”的方法.而线性不等式系统的可行问题则可以利用(整数)线性规划工具进行求解.相对于现有诸如Redlog、RegularChains等基于柱形代数分解的量词消去工具, 本文的算法1可以在可接受时间内进行复杂秩函数的计算.同时, 也不同于基于SDP的方法, 通过本文方法计算得到的函数是精确的秩函数, 因此不必再次验证计算所得的函数是否为秩函数.

致谢 感谢匿名审稿人对本文工作提出的宝贵意见.同时, 也感谢中国科学院成都计算机应用研究所的杨路先生对本文的修改提供的好建议.