琼州大学-谱分析：

最大熵谱估计 清华大学《现代信号处理》讲义 -张贤达http://www.docin.com/p-391781202.html

著名物理学家马大猷院士的著作《现代声学理论基础》中的阐述

一个系统中，熵的增加与热量增加成正比，

dS=dQ/T

绝热过程dQ=0，也可称等熵过程，dS=0或

dS/dt=0

用全微熵，S不因任何运动而改变。全微熵可写做

dS/dt+vÑS=0

欧拉质量守恒方程是

dr/dt+vÑr=0

二式结合，可得熵的连续性方程是

drS/dt+ Ñ×·(rS v)=0

可以证明，在流动（或波动）中，任何能量损失（如黏滞性，热传导等）都要使系统的熵增加。在热机学中，熵的应用更多。在所有过程中，熵不能减少，

dS/dt³0

这是热力学第二定律。（气体的内能增加等于内能增加和对气体所做的功（增加的动能），dE=dQ-PdV/T，这是热力学第一定律。对气体加压力时P时，其体积要缩小，所以这一项用负号，体积V减小，一定质量M的气体，其体积与密度成反比，rV =M。岁以上面的内能式子也可写成dE=dQ+(P/r²)dr ），热力学第一、二定律都普遍适用于一切热力学系统。

谱分析的基本步骤：

1.模型的建立

*非线性行为无效；

*一定要定义ex和dens

2.获得模态解

*只能用subspace，和black lanczos法；

*所提取的模态数应足以表征在感兴趣的频率范围内结构所具有的响应；

*不要同时进行扩展模态计算，以便在没有扩展的情况下进行有选择的扩展；

*有阻尼必须在模态分析中定义；

*必须在施加激励谱的位置添加自由度约束；

*求解结束退出solution处理器。

3.获得谱解

*指定分析类型为single-ptresp（单点响应谱）；

*no.of models for solu(模态扩展数);

*响应谱的类型：位移，速度，加速度，力等（命令：sytyp）；

*激励方向（命令：sed）

*谱值与谱线的关系曲线（freq和sv）；

*阻尼选项；

*开始求解。

4.扩展模态

*只选择有明显意义的模态进行扩展；

*扩展后才能合并；

*选择应力计算；

5.模态合并

*用srss法；

*指定输出结果类型：disp，velo，acel

6.观察结果

*读入jobname.MCOM文件；

*显示结果

谱是谱值和频率的关系曲线，反映了时间-历程载荷的强度和频率之间的关系。

响应谱代表系统对一个时间-历程载荷函数的响应，是一个响应和频率的关系曲线。

谱分析是一种将模态分析结果和已知谱联系起来的计算结构响应的分析方法，主要用于确定结构对随机载荷或随时间变化载荷的动力响应。谱分析可分为时间-历程分析和频域的谱分析。时间-历程谱分析主要应用瞬态动力学分析。谱分析可以代替费时的时间-历程分析，主要用于确定结构对随机载荷或时间变化载荷（地震、风载、海洋波浪、喷气发动机推力、火箭发动机振动等）的动力响应情况。 谱分析的主要应用包括核电站（建筑和部件），机载电子设备（飞机/导弹），宇宙飞船部件、飞机构件，任何承受地震或其他不规则载荷的结构或构件，建筑框架和桥梁等。

功率谱密度（Power Spectrum Density）：是结构在随机动态载荷激励下响应的统计结果，是一条功率谱密度值-频率值的关系曲线，其中PSD可以是位移PSD、速度PSD、加速度PSD、力PSD等形式。数学上，PSD-频率关系曲线下面的面积就是方差，即响应标准偏差的平方值。

ANSYS谱分析分为3种类型：

响应谱分析（SPRS OR MPRS）

ANSYS响应谱分为单点响应谱和多点响应谱，前者指在模型的一个点集（不局限于一个点）定义一条响应谱；后者指在模型的多个点集定义多条响应谱。

动力设计分析（DDAM）

动力分析设计是一种用于分析船舶装备抗震性的技术

随机振动分析（PSD）

随机振动分析主要用于确定结构在具有随机性质的载荷作用下的响应。

与响应谱分析类似，随机振动分析也可以是单点的或多点的。 在单点随机振动分析时，要求在结构的一个点集上指定一个PSD；在多点随机振动分析时，则要求在模型的不同点集上指定不同的PSD。

结构静力分析，动力学分析，显式动力学，屈曲分析，接触非线性分析，疲劳寿命分析，模态频率分析，谐振响应分析，跌落碰撞分析，热结构耦合分析，管道流体分析，流固耦合分析，弹塑性材料CAE分析，材料非线性力学分析，岩土材料非线性分析，橡胶材料超弹性分析，压力容器分析，货架分析，应力应变分析，受力分析，机械设计分析，汽车结构分析，电子电器分析，工程机械分析，瞬态分析，谱分析，汽车有限元分析

§5－10、最大熵谱估计

1967年Burg提出，对未知延迟点上的自相关函数值，按最大熵的原则进行外推，克服了经典谱估计中窗函数法，相当于令未用到的自相关或观测数据均为零的假设。故提高了分辨度，尤其适用于短数据的情况。

1、 最大熵准则：

熵－不确定度，最不确定的事件－熵最大。

定义：若随机向量具有概率密度函数则其熵为

熵－方差－相关函数

☆ 一维高斯分布

其中

代入熵的定义公式，并注意到得：

☆ N维高斯分布

，表示行列式；

故可得N维高斯分布的熵为：

为方便起见改写熵的定义为：

C－常数，可用来确定量度熵的参考基准或熵的绝对值

N－正整数。

因此令时，N维高斯分布的熵可写为：

显然，要使熵最大，就要使最大。

研究熵与功率谱之间的关系，即已知自相关函数的N+1个值,计算下一个最大延迟之外的自相关函数值。选择的准则就是要使熵最大。由自相关矩阵的正定性，行列式是非负的。要使熵最大，就是使最大，于是有：

由于N+1阶的自相关矩阵可写为：

式中

根据分块矩阵行列式的恒等式：（证明见王宏禹书p.297）

于是有：

显然(5-106)式，即将上式对的导数为零，得：

此式为的一次方程，从而求解此式，可得到合适的。同理将此代入，继而可求出。这就是在最大熵的原则下，由已知的N+1个值外推得到,,自相关序列。故可得：B为常数

外推的方法很多，但应该是与已知点上的自相关相符的功率谱中最任意的，即具有最大熵的。这相当于扩大了自相关的信息，故得到的谱估计比传统方法的分辨率高。

当序列无限长时，上述熵的定义会发散，则应采用熵率的定义：

相关图法中有的假设，即为其主要缺陷。

其功率谱估计为：

这里不采用此假设，而是按最大熵的原则即取值最不确定，来选择的值。其数学表达式为：

可以证明最大熵谱估计的表达式为：

式中的系数就是Yule-Walker方程的解。

最大熵谱估计的问题就变成估计这些系数的问题，可以证明最大熵谱估计与AR模型谱估计的是等价的。（有关证明可见陈炳和“随机信号处理”p.415）

§5－11、最大似然谱估计

Capon空间阵列接收数据，线性最小方差，地震、水声信号等。

当信号形式为时，第k 个传感器的接收信号为：

传播延迟，－白噪声过程。该传感器的输出延迟了，再乘以权，于是可得如下形式的输出：

第m 个传感器的输出为：

,且有：

若为互不相关的零均值，方差为的随机过程，则均方功率为：(为信号功率)

若定义,则其采样输出为

其方差为：

式中采样间隔。要使方差在下列约束条件下为最小：

式中

此估计器对x(t)是在高斯环境下进行的，所以此时最大似然估计与最小方差估计一致，即为：

由于为常数，故应使,即为：

易见，估计器参数的选择，是在一定的约束条件(5-119)下，使上式成立。故引入拉格朗日乘子，构造目标函数：

令则可求出：

代入方差表达式得：

当向量中的频率ω等于信号频率时，上式就代表信号功率

当频率ω的值遍及定义域时，上式就可以看作是过程x(t)的功率谱密度的估值。这样就得到了最大似然谱估计的值，为：

可见，自相关函数矩阵--->。

可以证明：最大似然谱估计和最大熵谱估计之间有密切联系，即：最大似然谱估计的倒数等于m=1—p阶的最大熵谱估计的倒数和。

式中为m阶AR模型（最大熵）谱估计。通常，先用自回归AR模型法计算 再用上式(5-125)计算最大似然谱估计。

一般说来，最大似然谱估计的分辨率要比最大熵谱估计的分辨率低。

§5－12、ARMA模型谱估计

自回归滑动平均ARMA信号模型，实际上是自回归AR模型与滑动平均MA模型综合应用的一种信号模型。从数字滤波器设计的角度来看，全极点的无限冲击响应的滤波器(AR)要比全零点的有限冲击响应的滤波器(MA)的阶次低，就可以实现同样特性指标的频率响应。而在无限冲击响应的滤波器中，要实现同样指标的特性，采用零极点型的滤波器(ARMA)又要比全极点型滤波器(AR)阶次更低。所以零极点型滤波器的特性优化逼近程度会优于全极点型的，阶次低，待定系数少，运算速度更快。显然、ARMA信号模型谱估计法的分辨度比最大熵谱估计的还要高。但其模型复杂、算法较繁，有许多问题仍在继续研究。

一般的ARMA信号模型：

对于一般的时域序列的信号模型，往往包含有确定的和非确定的两部分。确定的部分－AR模型，非确定的部分－MA模型。

从而组合成下列形式：

其中表示具有零均值和方差为的白噪声激励。将上式两边取Z变换得：

式中,

于是ARMA模型的系统函数为：

由此可得功率谱为（Z变换的形式）：

ARMA模型系数的确定：

用代入上式，可得功率谱的估值：

若用表示的反变换序列，对照上式，并注意到，则有：

为方便起见令，将上式两边同乘以,并令Z的等幂次项相等，则可得到下列关系式：

式中是ARMA信号的自相关函数。值此，我们可以先用Yule-Walker方程估计自回归AR信号模型的参数，再用上式(5-132)求出，于是可由上述(5-131)式，求出ARMA的功率谱。易见，ARMA模型是自回归AR模型（最大熵）谱估计的一种推广。实际上，求ARMA模型的方法很多。关于ARMA模型阶次的确定方法比较复杂，属于高度非线性优化的问题，目前已有一些结果，可见有关参考书。

例5-4、用ARMA模型估计白噪声与两个正弦波组合的序列。

两正弦波：

第一组：f₁-SNR₁=5dB, f₂-SNR₂=10dB

第二组：f₁- f₂-SNR₁=SNR₂=0dB；取N=100,p=q=10,结果如图。

第三节 相关分析及其应用

一、相关系数与相关函数

对于两变量、之间的相关程度可以采用相关系数表示，相关系数定义为

式中，为数学期望；，分别为随机变量和的均值，，；，分别为随机变量和的标准差，且

根据柯西-许瓦兹不等式：

故知。

二、自相关及其应用

1、自相关函数的定义

2、自相关函数的性质

(1) R_x()为实偶函数，即R_x()＝R_x()。由于

即R_x()＝R_x()，又因为x(t)为实函数，所以自相关函数R_x()为实偶函数。

(2) 时延值不同，R_x()不同。当=0时，R_x()的值最大，并等于信号的均方值。

则

这说明变量x(t)本身在同一时刻的记录样本完全成线性关系，是完全相关的，其自相关系数为1。

(3) R_x()值的范围为

同时，由式得

(4) 当→时，x(t)和x(t+)之间不存在内在联系，彼此无关，即 ，如果均值=0，则R_x() →0。

根据以上性质，自相关函数R_x()的可能图形如图5-16所示。

图5-16 自相关函数的性质

(5) 当信号x(t)为周期函数时，自相关函数R_x()也是同频率的周期函数。

若周期函数为x(t)= x(t+nT)，则其自相关函数为

例 求正弦函数的自相关函数。

解：此处初始相角是一个随机变量，由于存在周期性，所以各种平均值可以用一个周期内的平均值计算。

根据自相关函数的定义

式中，T₀为正弦函数的周期，T₀=，即

可见正弦函数的自相关函数是一个余弦函数，在=0时具有最大值。它保留了变量x(t)的幅值信息x₀和频率信息，但丢掉了初始相位信息。

例 如图5-17所示，用轮廓仪对一机械加工表面的粗糙度检测信号a(t)进行自相关分析，得到了其相关函数Ra()。试根据Ra()分析造成机械加工表面的粗糙度的原因。

图5-17表面粗糙度的相关检测法

解：观察a(t)的自相关函数Ra(t)，发现Ra(t)呈周期性，这说明造成粗糙度的原因之一是某种周期因素。从自相关函数图可以确定周期因素的频率为。

根据加工该工件的机械设备中的各个运动部件的运动频率(如电动机的转速，拖板的往复运动次数，液压系统的油脉动频率等)，通过测算和对比分析，运动频率与所测结果接近的部件的振动，就是造成该粗糙度的主要原因。

三、互相关及其应用

1、互相关函数的定义

若、为两个不同的信号和，则把称为函数x(t)与y(t)的互相关函数，即

2. 互相关函数的性质

(1) 互相关函数是可正、可负的实函数。因为x(t)和y(t)均为实函数，R_xy()也应当为实函数。在=0时，由于x(t)和y(t)可正、可负，故R_xy()的值可正、可负。

(2) 互相关函数是非奇函数、非偶函数，而且R_xy()= R_yx()。对于平稳随机过程，在t时刻从样本采样计算的互相关函数应与t时刻从样本采样计算的互相关函数一致，即

互相关函数不是偶函数，也不是奇函数，R_xy()与R_yx()在图形上对称于纵坐标轴，如图3.20所示。

(3) R_xy()的峰值不在=0处。R_xy()的峰值偏离原点的位置反映了两信号时移的大小，相关程度最高，如图3.21所示。在时，R_xy()出现最大值，它反映x(t)、y(t)之间主传输通道的滞后时间。

图5-18互相关函数的对称性

图5-19 互相关函数的性质

(4) 互相关函数的取值范围：

(5) 两个统计独立的随机信号，当均值为零时，则R_xy()=0。将随机信号x(t)和y(t)表示为其均值和波动分量之和的形式，即，

则

因为信号x(t)与y(t)是统计独立的随机信号，所以。所以。当时，。

(6) 两个不同频率的周期信号的互相关函数为零。由于周期信号可以用谐波函数合成，故取两个周期信号中的两个不同频率的谐波成分

进行相关分析，则

(7) 两个不同频率正余弦函数不相关。

(8) 周期信号与随机信号的互相关函数为零。由于随机信号y(t+)在时间t→t+内并无确定的关系，它的取值显然与任何周期函数x(t)无关，因此，R_xy()=0。

例 求的互相关函数R_xy()。

解：

由此可见，与自相关函数不同，两个同频率的谐波信号的互相关函数不仅保留了两个信号的幅值x₀、y₀信息、频率信息，而且还保留了两信号的相位信息。

3、互相关函数的应用

(1) 在混有周期成分的信号中提取特定的频率成分（在噪声背景下提取有用信息）

图5-20 利用互相关分析仪消除噪声的机床主轴振动测试系统框图

（2）用相关分析法分析复杂信号的频谱。

图5-21利用相关分析法分析信号频谱的工作原理框图

(3) 线性定位和相关测速。

图5-22利用相关分析进行线性定位实例

第四节 功率谱分析及应用

信号经过频谱分析，可以求得信号的频率成分和结构，并进而分析系统的传输特性，通过频谱分析，还可以对被测对象进行振动监测和故障诊断。

信号频域分析是采用傅立叶变换将时域信号x(t)变换为频域信号X(f)，从而帮助人们从另一个角度来了解信号的特征。

一、自功率谱密度函数

1、定义

假定:

(1) x(t)是零均值的随机过程

(2) 没有周期分量

相关函数满足傅立叶变换的条件

是信号平均功率。

2、物理意义

图5-23单边谱与双边谱

自功率谱密度函数曲线下和频率轴所包围的面积就是信号的平均功率。

3、巴赛伐尔定理

在时域中计算的信号总能量，等于在频域中计算的总能量，这就是巴塞伐尔定理。

——能谱，它是沿频率轴的能量分布密度。

利用这一关系式，可以直接对时域信号进行傅立叶变换来计算功率谱。

4、功率谱的估计

在实际中无法按照定义的公式求随机过程的功率谱，只能用有限长度的样本记录来计算样本功率谱，作为信号功率谱的初步估计。

模拟信号：

数字信号：

计算方法：

5、应用：

（1）可分析系统幅频特性

通过对输入、输出自谱的分析，可以得出系统的幅频特性，但丢失相角信息。

（2）判断信号中是否有周期成分

（3）分析信号频域结构

二、互谱密度函数

1、定义

如果互相关函数满足傅立叶变换的条件 ，则定义

称为信号和的互谱密度函数，简称互谱。

2、功率谱的估计

模拟信号：

数字信号：

初步估计：有偏估计、非一致估计

，随机误差大，估计值不能用。

分段平均：分成段，每段时长。

当各段周期图不相关时，

为进一步增大平滑效果，可使相邻各段之间重叠，以便在同样之下增加段数。实践表明，相邻两端之间重叠50%效果最佳。

现代谱分析方法：最大熵谱估计、ARMA谱估计

3、应用：（1）用于计算系统的频率响应函数。

不仅包含幅频信息，同时包含相频信息，提供了一种用实验方法得到系统频率响应函数的方法。

（2）可在强噪声背景下分析系统的传输特性

根据相关特性，同频相关，不同频不相关。

由此可见，通过互相关可以去除系统中的噪音，但无法排除输入端测量噪音的影响。

（3）在线测试。
（4）评价系统输入——输出因果性。

相干函数。

第五章熵原理

第一、二两章我们重点是引出分布函数这个概念和它在气象领域的重要事例，这些都可以在不理会熵概念的条件下进行。第三、四章开始引入熵概念和它在气象上的应用，但这并没有讲明自然界关于熵到底有那些客观规律。

很显然，如果仅使用人们较为陌生的“熵”作概念游戏而不引用新的原理去揭示气象现象的某些内在机理，那么我们整个工作的必要性也就值得怀疑了。诚然在一定意义上讲，正是为了引用熵原理去解决气象问题，我们才费了那么多力量谈分布函数和熵概念。现在就正面介绍熵原理以为下一步用它说明气象问题做准备。

什么是熵原理?一些人会说这也就是指热力学第二定律。但是在笔者想来，熵原理似应包含比热力学第二定律更多一些的内容。我们有这种认识的基本依据是一个世纪以来熵概念的含义有了不少扩展，因而熵的原理也应当适用于后来认识的这些“熵”。

在我看来申农把他定义的不肯定性也称为熵可能是熵认识史上继玻耳兹曼以后的最重大的贡献，但至今也应当遗憾地指出搞物理的仍主要是讲热力熵，而申农的信息熵则主要用于通迅理论和数学中。申农通过熵概念搭起来的桥似乎过往的人并不多，这就造成了横跨热力学与信息论的熵概念究竟有没有普适的熵原理问题长期没有受到应有的理论重视。

在笔者看来，理论界不仅应当给人们一个统一的对熵的认识，而且也应当给出一个统一的可用于一切领域的熵原理。显然，说热力学第二定律适用于天体、地理、生物体在内的任何热力过程还说得过去，可是说热力学第二定律也适用于非热力学过程就显得自相矛盾了。

所以对应用于各领域的熵概念，阐明其普遍适用的原理，在我看来是件重要而尚待完成的工作。在这个普适的熵原理中，过去我们找到的热力学第二定律或者信息熵的规律都应当仅只是它的特例。遗憾的是人们至今还没有完成这种综合工作。我们在这一章中则把几个领域中似应看做熵原理的内容都做些讨论。

§1 熵增加原理

上世纪中叶，继热力学第一定律之后又发现了热力学第二定律，克劳修斯对它的表述是^[1]：让热量自发地从低温传向高温而不引起任何其他影响是不可能的。这句话听来并没有什么有悖常理的惊奇之处，然而由此(也可以从开尔文，普朗克等人的其他等价表述)作为推论的一个前提却可以(逻辑的、数学的、物理的)推出一系列十分重要的论断、定律、函数来。

熵增加原理是说一个与外界隔离(无物质、能量交换)的物质系统，它的熵只能自动的增加或保持不变(不减少)。它是热力学第二定律的一个重要推论，而有的文献上直接把它视为热力学第二定律。

在第三章我们曾指出在可逆的卡诺循环过程中，热机的效率η=（T₁-T₂)／T₁，对于工作于冷源为T₂，热源为T₁的其他的不可逆的循环过程，也可以把一部分热变成为功。但是从热力学第二定律不难推论出，任何其他工作于T₁，T₂之间的热机，其热量转化为功的效率至多是与卡诺循环的效奉相等或更低，而决不会高于它。

在第三章第二节，已经通过(3.8)式表明熵是介质的热力学状态的函数，其变化量dS与可逆过程中的Q/T成正比。把这个结论与不可逆热机效率只会低于卡诺热机效率联系起来就可以推出另一个结论，即不可逆过程中一定的介质如从状态A到达状态B，那么它的熵的变化量(由可逆过程中对δQ/T积分求得)要比此不可逆过程δQ/T的积分值大。这写成式子就是

dS>δQ/T (不可逆) (5.1)

把这个式子与可逆过程的dS=δQ/T[即(3．6)式]合并就得出

dS≥δQ/T (5.2)

上式等号用于可逆过程，不等号用于不可逆过程，它表明不可逆过程中实际的熵增量比可逆过程的熵增量要多。

从热力学第二定律的原始表述可以推出上述结论，其推导过程在一般热力学书上都作了介绍，因而我们不再重达了。

如果一个工作介质或物质系统与外界隔开(孤立系统)，它就不会与外界交换热量，因而δQ=0。此时工作介质(说成物质系统可能妥当)内部可能进行着各种可逆的或不可逆的过程，从(5.1)式看，δQ=0导致dS≥0，即孤立系统的熵仅能加大或不变，这就是熵增加原理。所以一般把下式称为熵增加原理

dS≥0 (5.3)

一只绝热的容器内盛着冷水，当放入一个热的鸡蛋后，鸡蛋就把热量从它自身(高温)传向冷水，这是一个不可逆过程。把容器内的水与鸡蛋合起来看成一个物质系统，则它与外界交换的热量为零，可是物质系统内的热力状态在改变着(水变热，蛋变凉)，其物质系统内的总熵则依(5.3)式就变大。当水与蛋的温度相等后，内部状态不变化了，此时dS=0，即状态不变，其熵也不变(见图5.1)。

图5．1 孤立系统内熵增大到极大值后就保持在熵极大的状态下

请注意，上述物质系统的内部进行的不可逆过程中先是系统的总熵(水和蛋)不断加大，而在内部的温度达到均一之后其熵也就达到了条件允许下的极大值。此后这一孤立系统的热力状态就总是稳定地处于温度均一的不变状态下，此时系统的熵就仃滞在极大值这个数值而不会自动变小。

所以，所谓熵增加原理，一则是说孤立系统熵仅会加大不会减少，另一则是说熵是有上限的(极大值)，它会自动加大到它力所能及的极大值，并会长时间地稳定在熵极大状态下。

　　一个比重大于1的物体落入水中，它会自动下沉，最后它会长时间的停滞在位能达到极小值的状态下。这与熵增加原理是很类似的。

　　人们对于热力学第二定律的原始表述如此浅显易懂而感到亲切，可人们也对由上述原理引出来的熵增加原理感到陌生。一般的原理、定律总是给出一个准确的只含等号的关系式而不含不等式，可熵增加原理竟是个不等式，可以在一定意义上讲它仅只告诉了“+”、“—”号而没有给出计算值。确实，140年来关于熵原理的非议一直不断。

这里我们不去评论这种非议，而是从另一个侧面谈一下科学界的先驱人物如何从这个看上去古怪的原理进而导出十分有用的定律的。

有了熵概念后，亥姆霍斯和吉布斯又提出两个热力学函数，即亥姆霍斯自由能和吉布斯自由焓。自由能F是由下式定义的

F=E-TS (5.4)

E、T，S分别为该物质系统的内能、温度和熵。物质的自由焓G是由下式给出

G=E+PV-TS (5.5)

这里的P、V是该物质的压力与体积。

有的化学教科书上称这是化学中最重要的公式。确实，这些公式对于分析各种化学变化和物态相变过程是十分有用的，其关键点在于从熵增加原理可以进而导出自由能和自由焓的减少原理。在化学反应过程或相变过程中，在等温等容或者等温等压条件下自由能或自由焙分别达到极小值体现了熵增加原理，又进而引出化学或物相变化达到平衡的条件，从这里可以引出一系列十分有用的公式来。

水的饱和蒸汽压随着温度而变化是十分显著的，其变化关系就是气象学中熟知的克拉珀龙(Clapeyron)方程，这个著名公式之所以能得出来可以追朔到熵增加原理。

在化学中如果一个化学反应的生成物与反应物都是气体，反应进行到各种成分各是多少时就达到平衡态呢?在等温、等压条件下这个平衡态对应于自由始最小(熵最大)，由此引出了著名的“质量作用定律”。

这些公式在工程界和科技界都是广为引用的，遗憾的是不少人并不知道导出这些公式时都利用了熵原理。

普朗克的突出贡献当然是引进了能量仅能取离散值的量子化思想，可是他由此假设参与推导而得出黑体辐射与频率的著名公式时，也同时利用了熵原理。

有人说熵原理好象是用来生产“公式”的培养基，确实就在一些人对熵原理感到别扭、一些人议论宇宙走向热寂的同时，经一批著名学者的努力，熵原理早巳为现代科学作出了突出贡献。

为便于看出这几个重要事例与熵原理的逻辑关系，我们会了一个逻辑框图(图5.2)，它反映了熵原理在得出这几个著名公式时所起的作用。

图5．2 在导出几个著名公式时熵原理的地位

§2 最大熵原理

近几十年来统计数学在社会科学、医学、生物、气象等领域广为应用。它的理论基础是概率论。在概率论中，介绍各种各样的概率分布函数是其重要内容，此外还介绍这些分布(如正态分布、二项分布…)常常适用于那些实际问题。

但是在概率或统计的书籍中一般地说，它们并不细究为什么某些自然现象恰恰遵守某种分布函数。而不少科技人员则对于自己收集来的数据能证实符合某种分布(如正态)已经很满意了，他们似更无心多问一句为什么。

自然(包含社会)现象A符合甲分布，现象B符合乙分布……，难道没有一个更深一层的原理来统一地说明这究竟是为什么吗?

医生发观某些疾病符合正态分布律，他会认为这已经升华到了理论殿堂了，可是在理论家看来这仅是证实，而还没有构成一个理论来证明该现象必然符合正态分布。赤池弘次把话说的更难听，他说^[3]：“统计学的课本往往被看做应用上有用方法的大杂烩。”

看来今天的统计理论水平比起动力学来就是要低一个档次，或者说，统计现象的理论框架尚待构筑。

在笔者看来最大熵原理恰好可以填补这一理论空白。从这个原理出发可以对众多外形不同的分布做出统一的理论说明，这不仅使大量的统计现象得到深一层的理论说明(又多回答了一个为什么)，也使统计理论水平提高一步。统计数学不再仅是有用方法的大杂烩了，而最大熵原理与热力学第二定律的融合会更加显示出熵原理的重要意义。看来，爱因斯坦把熵原理作为整个科学的第一法则的提法就更令人感到亲切了。

E．T．Jaynes(杰尼斯)是信息论中的最大熵原理的创导者，他在此领域的工作从50年代跨到80年代^[4，5，6]。他的这一思路逐步为社会所公认^[7，8，9]，我们认为把他的这套思路也归入熵原理是公正的。

2．1 数学思路

离散变量x仅能取分立的若干个值。并有相应的概率值p与之对应，这可构成下表

表5．1 分立的变量值和与之对应的概率值

x1	x2	x3	…	xn
p1	p2	p3		pn

依表5．1的记法，离散变量的熵公式(3.15)在C=1时可简写成

(5.6)

这里要求

p_i≥0

(5.6)式中把H写成H_n。是想说明H的值是n个变量p₁，p₂,…,p_n的函数。

现在把问题反过来考虑，如果不是已知概率p_i值而是已知熵H_n达到了极大值，那么这对P，的要求又会是什么呢?

这个问题是十分有价值的，而其具体的答案则与具体的约束条件有关，即约束条件不同答案就不同(以下论证主要参考了文献[10-12])。

首先我们看到H_n的值至少要受到∑p_i=1的约束，即H极大仅能在各概率合计值恰为单位值1的条件下取得。

一般地说，其他的约束经常以变量x的某种函数f(x)的平均值为已知值的形式出现。如尚有m个约束条件(m<n)， 它是说m个x的已知函数f_l(x)，f₂(x)，…，f_m(x)都有事先确定的平均值F₁，F₂，…，F_m，即

k=1,2,…,m m<n (5.8)

现在问在(5.7)、(5.8)式的约束下各p_i取什么值恰好使(5.6)式的熵H_n达到极大值?

为此参照Lagrange乘子法[13]去构造一个新函数，它是H_n与常数α，β₁，β₂，…和F_k的如下线性关系

而依(5.6)、(5.7)和(5.8)式有

因对数具有x=lne^x的性质，在上式中的对数以e为底时(这不失一般性)可改写成

`

在本章附录中证明1nx≤x-l，将此不等式代入上式有

此式可变成

由于F_k是事先给定的系数，α，β是待定常数，所以上式表明H_n的值是p_i(i=1，2，…，n)的函数，而H_n的极大值应是让上式取等号时达到。依本章附录易于看出这对应于要求每个p_i满足下式

(5.10)

即如令各p_i恰满足(5．10)式，则熵H_n达极大值，这样就初步找出熵极大对p_i的要求。余下的是求出待定常数α和β_k (k=1,2,…,m)。

利用各p_i之和为1的(5.7)式可将(5.10)式变成

即

或

如令

Z=e^α (5.11)

上式可以写成

（5.12）

用(5.11)、(5.12)代入(5.10)就消去了α，使p_i的式子变

　　

顺便指出函数Z在统计物理中称为配分函数。为求得β_k (k=1,2,…,m)的值，把上式代回约束方程(5.8)有

k=1,2,…,m (5.14)

此式也可以简记为

k=1,2,…,m (5.15)

式子中F_k，f_k(x_i)都是已知值，因而真正的未知数是m个β_k值(k=1,2,…,m) 。m个方程应当能解出m个值， 这样，从原理上讲，我们就可以求出熵极大时的各个p_i值了。

把这样求的各p_i值代人熵表达式[即将(5.13)式代入(5.6)式]就求出了熵在对应约束下的极大值是多少。这可以。表示成

利用(5.7)、(5.14)可将上式整理成

这就是熵极大时计算其熵值的公式。

为了把以上推导过程的基本关系揭示清楚，这里给出一个框图，即图5.3。

图5．3 熵最大的条件与结果

2.2 若干分布实例

这里给出一些实例用以揭示如何在不同的约束条件下导出常用的若干概率分布函数。

2.2.1 等概率分布

(5.7)和(5.8)是求熵极大的约束方程，如果约束条件仅有(5.7)式，那么熵极大对应的概率分布是什么呢?

此时可以理解为在(5.8)式中f_k(x_i)≡0，从而使(5.10)式变成

p_i=e^-α

这表明熵极大时p_i应当是常数(各p_i相同)。换言之，如不附加进一步约束，那么等概率的分布也就是熵最大的分布。

如果一个随机变量x有n个可分辨的状态，那么每个状态的出现概率都相等时，这个离散变量的信息熵恰好达到极大值。依(5.7)式可得

p=l/n (5.17)

这就是各p_i相等时的概率值(见图5.4)，而其熵则为H=lnn。

图5．4等概率分布

等概率分布是经常遇到的一种分布，它对应于约束最少情况下熵最大时对概率的要求。统计物理中微观粒子处于各能级状态并无多少约束，因而各微观态出现概率相等就成了它的基本假设。这种分布有时称为正则分布.

2.2.2 均匀分布

前面导出的熵最大时的概率分布是在离散场合下得到的，如果把各公式中的求和()都改由积分代替，而所有的概率（p_i）该由概率密度p(x)代替，那么离散场合的结果也可以用于求解连续变量的场合。现把主要关系汇集于下

与（5.6）对应的连续变量熵公式H变成

（5.18）

b，a是x的上下限。与（5.7）对应的公式是

， p(x)≥0 （5.19）

与(5.8)对应的公式是

k=1,2,…,m （5.20）

与(5.10)对应的是

（5.21）

而(5.12)变成

（5.22）

(5.13)对应为

（5.23）

而（5.14)变成了

k=1,2,…,m （5.24）

而(5.15)、(5.16)两个式子并不改动。

这样我们就把连续变量场合求熵极大时对应的分布函数有关公式都引伸了出来。

如果有一个连续型的随机变量x仅知其出现于区间[a,b]之间，而没有其他约束，则(5.22)式变成

将此代入(5.23)式即可求出概率密度函数p(x)：

p(x)=1/(b-a) a≤x≥b ( 5.25)

这正是均匀分布的概率密度函数(见图5.5)。故依(5.16)式可求出信息熵的极大值

Hmax= ln(b-a) (5．26)

不难看出这个结果与离散变量的等概率分布很类似。

图5.5 均匀分布

2.2.3 指数分布

如连续变量x仅能出现于(0-∞)之间，而且知道其数学期望值(平均值)为有限值u(不是无限大)，此外如不附加进一步的条件，那么熵最大时对概率密度分布p(x)将有什么要求呢?

显然，此时(5.19)式变成了，而(5.20)式仅为一个关系(m=1)：

(5.27)

由此得配分函数Z=[exp(-βa)]/ β 。而

p(x)= βexp[-β(x-a)]，利用(5.27)可得β=1/(u-a)。故

(5.28)

如下限。a =0则

(5.29)

这就是典型的负指数分布律，其图形见于图5.6中。它表明一个变量如有下限(如0)，且有给定的平均值(平均值对变量的概率分布起了约束作用)，那么熵最大对应(要求)的概率分布恰好是指数分布(5.28)式。在下限为零时它简化为(5.29)式。

由于此时Z=[exp(-βa)]/ β，不难从(5.16)式求出负指数分布时的熵值H=ln(u-a)e。a=0时H=lnue。

图5．6 负指数分布

2.2.4 正态分布

如知变量x的平均值为u，标准差为σ，问x遵守什么分布才使熵最大?

此时不仅有，而且(5.20)式实为两个， 即f₁(x)=x, f₂(x)=(x-u)²故有

(5.30)

和 (5.31)

这使配分函数Z变成

经分部积分和代入u的定义式得β₁=0，对上式作积分可得，再利用标准差的公式最后得β₂=1/(2σ²)。

由此代回(5.23)式最后得

(5.32)

这正是统计学中著名的高斯(正态)分布(见图5.7)。这表明给定了标准差σ。(其平方称为方差)和平均值u，在熵最大的要求下概率应遵守高斯分布。这样，正态分布的原因就追朔到熵最大这个一般原理上去了。

图5．7 正态分布

　　不难进一步求得正态分布下的熵值与平均值u的大小无关。其值。

以上几个示例已经看出在统—的熵最大的要求下，约束条件不同，得到的分布形状也可以大不相同，改变约束还可以导出多种分布来。在本书的附录C中， 我们把这些分布的约束与分布函数统一列表，以便读者查用。关于这类结果如何用于气象上的分布问题中将在下一章讨论。

2．3 物理原理

在以上的论述中，抽象的数学推理与数学假设是主体，人们会问，这种与物质世界分开的假设与原理，难道有资格称为物理学(物质科学)的规律吗?

确实，从图5.3看，熵最大仅只是解某些数学问题时的一个要求或假设，前面的推导计算和分布实例仅只是数学演算步骤，它们都不是最大熵原理的证明，都未涉足物理。

那么最大熵原理又是什么呢?文献[11]中对此有一段评述，我们把它引于下面。

“有一些随机事件我们不能直接计算其概率， 还有一些随机事件也不容我们计算出相应的出现频率。通常，我们所能掌握的仅只是关于该随机事件的一个或几个平均值。例如，对微观系统所测得的宏观尺度的数据就是某些随机变量的这种平均数值。然而，对于给定的随机变量的平均值来说 (指观测到此平均值之后——笔者)，会存在着很多种(甚至无限)分布模式都与这个平均值兼容，问题在于如何挑选出最为合适的一个分布来。显然，为此我们就得介绍评选的判据，而最大信息原理就可以看成是其判据。

“依照最大信息原理(指信息熵最大——笔者注)，我们要挑选在一定约束下(常常是给定某些随机变量的平均值)使得熵(或条件熵)能极大化的那种分布作为选定的分布。这个原理是独立地为R．S．Ingarden (1963)^[14]，E．T．Jaynes^[4]，(1957)和S．Kullback与R．A，Leibler^[15]建立的，这可以看成是Laplace(拉普拉斯，1749—1827笔者注)的著名的“充足理由律”的新的和重要的发展。这个定律认为如果我们处于对随机事件一无所知的场合下，应假设它是均匀分布。

“当然，最大信息原理有其主观性，不过在构造一个随机分布时，应当把它看成最客观的主观准则。实际上，只要我们承认熵是计量不肯定性的最合适的标尺，我们就完全有权在给定约束下选择不肯定性最大的那种分布作为随机变量的分布。从其含有的不肯定性的角度看，这种随机分布是最为随机的分布。依照最大信息原理，我们就要在给定的那些约束下挑选不肯定性最大的那种概率分布。最大信息原理在经典和量子统计力学中的成就使我们又扩大了它的应用领域。”

作者在文中指的“最大信息”实际上就是最大信息熵，我们这里统一归入“最大熵”这一名词之下。

作者的上述评述突出了最大熵原理是一种判别标准。利用这一准则在给定约束下求得的概率分布，是主观成分最少，把不确定的东西(以熵计量)作最大估计下的分布。从认识论和方法论的角度看，这么作都是较为客观、较为合理的，这一套办法在实用上是成功的。Burg在60年代创立的频谱分析新方法——最大熵谱得到广泛应用就是其例。

在上述引文中的最后一句话提到了这一原理在统计力学中的成功应用。在笔者看来，这种成功不能仅从认识论和方法论的角度去理解。应当看到它具有本体论的意义。

在不少场合最大熵原理本来就是客观事物自身的客观规律，它并不依赖人类是否存在，人是否进行某种观测。最大熵原理应当与热力学第二定律融合在一起共同构成一个客观规律，它是不以人的意志为转移的客观存在。

依照现代观点，任何物质系统除了都受到或多或少的外部约束外，其内部总是具有一定的自由，这种自由会导致系统内的各个元素处于不同的状态下。而状态的多样性、状态的丰富程度(混乱程度、复杂程度)的定量计量标尺就是熵。

熵最大就是事物状态的丰富程度自动达到最大值，或说事物总是在约束下争取(呈现)最大的自由权(自发达到最混乱、最复杂)。换句话说，当我们利用最大熵原理这一数学方法时，实质上是我们承认物质系统内的熵(这可能是热力学熵，也可能不是，而是某一种或几种信息熵)自动地应当处于约束条件所允许的最大取值状态下。

而上面这种提法就与解热力学问题时说我们利用了“热力学第二定律”几乎是等价了。人们从不把热力学第二定律看成是数学中的“估计”用的判据，而认为它是自然界的根本原理，它是物理原理而不是数学原理。

依此看来，在各个科学，技术领域中为什么有那么多自然现象遵守正态分布、指数分布、F分布……也就易于理解了。实际上，经过本节的数学说明，正是这些自然现象中都在受制于不同的约束的同时，还共同受制于最大熵原理。换言之

最大熵原理是很多自然现象经常遵守几种概率分布的根本原因

从信息论中发展起来的最大熵原理以及它应用于统计力学时取得的成功，几乎使我们应当把统计力学看成信息论的特例。信息论与统计力学的融合使我们看到熵概念的强大生命力，也使我们看到熵原理的重要性。现在看来把最大熵原理作为重要内容的熵原理， 再仅仅称作热力学的第二定律，或者仅仅看成是“统计估计”用的一种数学判据都是不够全面的。我们应当看到信息论和统计力学共同把熵概念和熵原理推向了新的高峰，对此我们不仅应当有相应的哲学认识，也应充分利用它把众多学科推向前进。

§3统计力学思路

科学应当回答支配各种现象的原因。当一个力学现象通过牛顿第二定律和相应的条件(初始的和边界的)而解释清楚时，似乎“科学”已经完成了它的任务，人们不再进一步追问下一个为什么了。

前面讨论了很多分布现象，并且指出通过最大熵原理和相应的条件也可以对此作出解释。

由于最大熵原理与热力学第二定律的亲缘关系，我们有权说，甲现象用热力学第二定律作了说明与乙现象用牛顿第二定律作了说明，两者的“科学水平”是一样的。换句话说，用最大熵原理阐明的事理与用力学原理阐明的事理至少在离终极真理的距离是一样的近。实际上，当今的科技人员中能用力学原理说明自然现象的人数远比能用热力学第二定律(或称最大熵原理)去说明自然现象的人数多的多。造成这种差别的原因如果不是由于“熵原理难以掌握”就是由于熵概念过分神秘，使人望而生畏了。

为克服对熵的神秘感，这里将沿另一思路导出负指数型的统计分布，它与熵原理得出的结果有异曲同功之妙，这显然对认识统计力学原理和认识熵都有益。

下面进行的分析，其思路源于上世纪玻耳兹曼的工作，现常被称为玻耳兹曼能量分布律。不过为节约篇幅，在论述中我们直接用于分析降水量(取代玻耳兹曼的分子动能)，这样作也加强它的气象色彩。

从云物理和动力学角度分析降水现象实在是相当复杂的(空气湿度、凝结核、云滴谱、上升运动、碰并过程等)。现在我们从统计力学思路分析这内中的一个问题：在一场降水中，不同深度的降水所占的面积有什么内在关系。

假设大气中发生了一个天气过程，它要把体积为V的雨水降落在面积为A的地面上。如果大气象一辆洒水车那样把体积为V的水洒向A区域，那么这几乎是非常规律、非常有序，因而也非常均匀的把V洒下来。均匀两字意味着每个点得到的雨水深度是一样的，其值应当是平均雨深。并有

(5.33)

可是气象观测证实从来就没有一次降水形成的雨量会象上述模型那么均匀。降水分布图上的等雨量线经常很复杂，如图2.17所示的，它经常是少数地域雨量很大，而多数地域的雨较少，因而“天气过程”把V体积的雨水抛向地面时是有很大的任意性的。

“任意性”一词在动力学中几乎使人束手无策，可是统计力学中微观尺度的任意性(随机性)常可导至宏观尺度的必然性。

对于A这么大的降水区，一个小小的雨量筒的横截面积应当说是太小了。如把落在A上的V连同A都视宏观量，雨量筒的截面积ΔA和落入其中的降水ΔV都可视为微观量(如设雨区止为20km×20km，雨量筒截面为20cm×20cm，则A/ΔA为10¹⁰，即100亿倍。

现假设在地面上(指降水区内)布满了雨量筒(样子象密蜂窝了)，其个数为N个(上例中N=100亿)。而雨量分布的不均匀就意味着仅有n₁个雨量筒收到了体积为v₁：这么多的降水，同时还有n₂，n₃，．．．．个雨量筒分别收到v₂，v₃，．．．．，这么多的降水。依此分析，每场降水都会得到下表(表5.2)。

表5.2 一场雨中收到不同的雨量的雨量筒的个数也不同

雨量（体积）	v1	v2	…	vi	…	vq
对应的雨量桶数量	n1	n2	…	ni	…	nq

在表5.2中，如v₁，v_q分别代表雨区内收到的雨量的最小值和最大值，那么n₁，n_q就是对应的雨量筒个数。

由于遍地布满了雨量筒，我们易于得出

(5.34)

(5.35)‘

即雨量筒个数总和为N，而各雨量筒雨量之和为V。

现在回过头来再分析降水的任意(随机)性因素。先考虑“从N个雨量筒中取出n₁ 个筒，；使其雨量恰为v₁ 的办法共有多少种”这样一个组合问题。

如果把N个雨量筒都编上号码，那么会有N个办法仅使一个筒里装有v₁的降水。如让两个筒装了v₁，则在装好第一筒后，第二只筒仅有(N-1)个办法了。故装两筒v₁的办法为N(N—1)个。如让n₁个筒装入v₁ ，依上述分析应有N(N—1)…·(N- n₁+1)个办法。

我们把N个筒中有n₁个装v₁的降水视为一个宏观状态，而实现这种结局的每个具体办法视为一个微观状态，这样就看到一个(指N中有n₁个为v₁)宏观态会有多种微观态与之对应。

能说N中取n₁的宏观态与N(N-1)…·(N- n₁+1)个微观态对应吗?细—想还不能。因为上述计算中是计较了谁是第—个装v₁谁又是第二、三个装入的，而这种顺序不同并没有形成有差别的微观态，这就要把由顺序不同(排列不同)引出的“个数”再扣除。而n_i个筒的顺序不同的排列办法有n_i个。所以把以上结果除以n_i!就是N个中取n_i，使其雨量为v_i：的方法个数了。依此对于N个中n₁的方法数自然应当是

此式也可写成

（5.36）

这正是从N个中取 n₁个组合时的方法个数。

从N中仅取n₁个为v₁事情并未了结。我们应再问，此后再令n₂个筒有v₂的雨量会有多少个方法?仿前，这应有

而N个中有n₁个为v₁，有n₂，个为v₂，的方法数W₂应为

仿此，最后会针对表5.2查问N个中有n₁，n₂，…分别为v₁， v₂，…的这种宏观态对应于多少个微观态个数(W)。不难从上面的式子中最后得出

W=N!/ n₁! n₂!... n_q!

(5.37)

这一番排列组合问题的分析有什么用呢?我们说上式表明每一宏观态(对应一组n₁， n₂，...， n_q)对应的微观态的个数W是N和n₁， n₂，...的函数，这可写成一般形式

W=f(N, n₁, n₂,...,n_q) (5.38)

换言之， N, n₁, n₂,...,n_q的值不同W的值也不同。统计力学中的一个基本假设是各个微观态出现的概率(机会)相等，上式给出了各宏观态对应的微观态个数W并不相等，因而易于想到微观态个数最多的那种宏观态是最易于出现的。这个思路十分重要，它引导我们进一步问可否从W达到极大值反求 n₁, n₂,...,n_q 各是多少(设N已知)。

现在我们不急于回答上述问题，而是先求几种n₁， n₂，...值下的W作个对比，以增加对总问题的感性认识。

设某区域A仅用100个雨量筒就布满了地面，且发生了总量为1000的一场降水，此时依(5．33)式， 平均雨量就是1000/100，即10。

这总量为1000的降水现仅选（假设）按如下4种方式降在A区内的100个雨量筒中，我们依(5.37)式求其实现这种宏观分配方式的微观办法的个数W。

方式1：每个雨量筒收集的降水恰好都与平均降水相等。

方式2：有99个雨量筒收集的降水恰好是9。

方式3：有50个筒得到的降水都是5，另50筒都是15。

方式4：使W达最大的那种雨量分配方式。

这4种降水分配方式(对应4种雨量分布图)的代表图形见于图5.8中。

10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10
10	10	10	10	10	10	10	10	10

图5.8a 方式1

此类分布仅此一种实现办法

9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	109	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9
9	9	9	9	9	9	9	9	9	9

图5.8b 方式2

此类分布有100实现办法

5	5	5	5	5	5	5	5	5	5
5	5	5	15	15	15	15	5	5	5
5	5	15	15	15	15	15	15	5	5
5	5	15	15	15	15	15	15	5	5
5	15	15	15	15	15	15	15	15	5
5	15	15	15	15	15	15	15	15	5
5	15	15	15	15	15	15	15	15	5
5	5	15	15	15	15	15	15	5	5
5	5	5	15	15	15	15	5	5	5
5	5	5	5	5	5	5	5	5	5

图5.8c 方式3

此类分布有1.26×10³⁰种实现办法

1	1	1	1	3	3	1	1	1	1
1	3	5	5	7	7	9	5	3	1
3	5	7	13	15	19	11	11	5	3
3	7	11	17	21	27	23	17	9	3
5	9	13	25	31	42	25	19	11	7
7	9	15	27	35	50	29	15	11	7
5	7	17	19	23	29	21	13	9	3
3	5	11	13	13	17	15	11	5	3
1	3	5	7	9	9	9	5	3	1
1	1	1	3	5	7	3	1	1	1

图5.8d 方式4

此类分布有9.47×10⁹⁷种实现办法

图5.8 4种(a,b,c,d)宏观降水分布对应的微观分布的示例 (图中每两种不同的降水量——交换位置就形成一种新的微观态，但是其对应的宏观分布却不变)。

在各筒降水都相等的方式1中n₁=100，而 n₂, n₃,...,n_q皆为零。N=100时代入n₁=100到(5.37)式得W=100！/100！=1。即此宏观雨量分布仅与一种微观分布对应（此时各地的雨两皆为10）。

方式2中n₁=99，v₁=9，这使n₂=1，v₂=109，代入公式得W=100，即99个点的雨量为9，仅一个点为109的这种宏观雨量分布与100种微观分布对应。 这100种就与雨量为109的那个点分别位于100个点中各个点相对应。

方式3中n₁=50，v₁=5而n₂=50，v₂=15。依公式求得W=1.26×10²⁰，即实现这种宏观降水分布的微观分布方式竟多达1.26×10²⁰这么多，而图5.8.c中仅绘了这么多方案中的一个方案。与前两者比，实现本宏观结局的微观办法实在太多了。

方式4的各个n₁，n₂…的值列于表5.3中，它对应于后边介绍的W达极大值时的那种宏观降水分布状态。W极大是说实现此种宏观分布的微观办法的个数实在太多，而且不能有其他的宏观态会对应与比它还多的方法数，此时W= 9.47×10⁹⁷。图5.8.d给出的仅是9.47×10⁹⁷种微观态中的—种。

		n1	n2	n3	n4	n5	n6	n7	n8	n9	n10	n11	n12
		v1	v2	v3	v4	v5	v6	v7	v8	v9	v10	v11	v12
a	n	100	0	0	0	0	0	0	0	0	0	0	0
	v	10	1
b	n	99	109	0	0	0	0	0
	v	9
c	n	50	50	0	0	0	0	0
	v	5	15
d	n	18	15	12	10	8	7	5	4	4	3	10	4
	v	1	3	5	7	9	11	13	15	17	19	21-29	≥31

表5.3 a,b,c,d 4 种分布方式

(实现a,b,c,d分布的微观方法w分别是1，100，1.26×10²⁰，9.47×10⁹⁷)

如果一个袋子中有1万粒好的葡萄干而有1粒是坏的，当你从袋子中任取一粒时抽中好的葡萄干的概率就非常高 (10000／10001)。与此类似，如果天气系统把1000的降水洒在100个雨量筒内时，每种微观态出现的机会相等，那么W值最大的那种宏观态必然出现机会极高，而W值太小的几乎不会出现^[16]。正是这种思路把我们引向求W达极大时n₁，n₂…应当各是多少这样一个数学问题。而每种微观态出现的机会相等是统计物理中的基本假设，在本问题中它对应任一降水元量落在那个雨量筒中的概率是一样的。

如何求W的极大值?W极大时n₁，n₂…各是多少?

由于lnW与W互为单值关系，W达极大时lnW也达极大，我们现在转而研究如何使lnW达极大。这样(5.37)式变成了

在前面的极简单的例子中我们取N=100，而实际构成这一思想实验时，可以让N充分大(如本节开始时的分析中N= 10¹⁰)，N充分大时n₁，n₂…也应相当大。这时注意到斯特林公式

lnn!≈nlnn-n

会使^[17]

(5.41)

　可是，根据古典概率定义N个事件中有n_i个有利于v_i，（N个筒中有n_i个的雨量为v_i)，则降水为v_i的事件从A区域内的N个筒中被抽中的概率p应当是p=n/N。所以依信息熵H定义的（3.13）式或（5.6）式有（取C=1）

lnW=NH (5.42)

这样我们看到求lnW的极大值的问题与求对应的信息熵H的极值问题等价起来了(N为固定值)。

换句话说，在本节中我们沿着较为细致的物理思路去寻找什么样的n₁，n₂…分布是最易于(最可能)出现的分布，可是这个问题最后变成了从N个样本中抽取n₁个的降水v₁, n₂为v₂ 的概率组合成的信息熵在n₁，n₂…为什么值时熵最大的问题。

这里我们看到不同的思路最后归结到同一个学数问题上去了，这种殊途同归一则使我们看到一个清楚的统计力学思路，还看到熵极大提法的简明深刻。

至此可以说已经讲明本节的主题，即统计力学思路导致的要求与熵最大导致的要求是一致的。下面再作些数学处理以得出具体结果。

由于这里(n_i/N)=p_i，故(5.41)式变成了(5.42)式，即lnW=NH。分析此问题时N并不变化，这使lnW达极大时H--信息熵也达极大。而信息熵H达极大是与一定的约束条件有关的。在这个降水问题中(5．34)、(5．35)两个式子分别被N去除，并注意p_i =n_i/N就变成了

(5.43)

(5.44)

这两个式子正是H极大的约束条件。因V/N为常数，所以(5.44)式实际是一般约束式(5.8)的特例。利用(5.6)、(5.43)、(5.44)式就构成了在约束条件(5.43)、(5.44)式下的降水信息熵的极值要求。

仿第二节补入拉哥朗日乘子时，因仅有(5.43)、(5.44)两个约束，故仅有αβ 两个待定常数[参见(5.9)式]。而依(5.10)式求得的解应为如下形式

p_i=exp[-α-βv_i] (5.45)

这是因为(5.10)式的f(x)在本场合下仅为x_i(即v_i)。而(5.12)式的Z应为

如各v_i相隔的比较密，则上式可用积分代替，从而有

依(5.11)，Z=e^α ，故(5.45)式变成了

(5.46)

为求β，把它与改写成积分形式的(5.44)联立，可得

(5.46)

这里用平均雨量代替V/N，v₁仍为雨区A内的最小降水量。

不难看出上式表示的是指数分布形式，它实际上与上一节求得的(5.28)式是一回事，它们都表示降水大的出现的机会少(对应于n_i 少，即占的面积少)，降水少的(v_i小)出现的机会大(n_i大)。

上式是把离散变量改成连续变量处理而得出之结果，因而p值具有概率密度的含义，即它代表了降水出现于v_i →v_i+1这单位增量间时它对应的概率。或说(5.46)式右侧对应于 dp/dv)。将此关系代入上式并改成差分形式，则有

(5.47)

把这个式子用于实际计算时应把Δp_i,再改写成Δp_i=n_i/N，如果说N对应于总的雨量筒个数，则n_i代表雨量在v_i →v_i+1这个范围内的雨量筒个数(它与相应的雨区面积成正比)。

现在用此公式去计算表5.3，最后一行中给出的数据。为此取Δv=2，v₁：=0[1]、N=100,在v=1、3、5，…、19等场合计算了n值，此后把Δv变成10和30以上又求出n分别为10和4。这样依(5.47)式就求出了表5.3中最后一行数据（n_i)。而W=9.47×10⁹⁷是利用这些n₁，n₂…值代入(5.37)式算出来的。

以上的讨论使我们看到统计力学思路的精神以及它与熵极大完全相同的具体结果。改变约束条件，也可沿统计力学思路导出其他的对应结果，这些就不一一举例了。读者处理问题时可自由选择走那条路，但结果都是一致的。

统计力学处理中实际有不少细节在这里没有全讲透，这主要是由于篇幅所限。有兴趣的读者可参阅一般统计力学书籍或者1985年新疆气象台油印的一册“降水统计力学初探”。（2007注：此书电子版已经在2006年公布在张学文个人网站上）

§4 信息的不增殖原理

有些科学问题从古远的年代就提出来了，但仅只近数百年才有公认的结论。“巧妇难成无米之炊”这句中国谚语说明物质的“质量”不能无中生有，可这作为一个科学结论则是化学家拉瓦锡提出的质量不灭定律；这大约是200年前的事。

永动机是否存在?或者问，能否制造一部机器，使其输出的能量多于输入?150年前被证实的能量不灭和转换定律对比也下了科学的，但也是否定性的结论。

这说明，在自然界，使质量和能量增殖的办法、机构、机器是不存在的。占卜术流传至今不只二千年了，占卜家的使命与永动机工程师十分相似：要找到一种办法、机构、机器，使输出的信息多于输入。

相信科学的人都认为占卜术的骗人，可他们也拿不出一个“定律”来对占卜术判死刑。

是信息论为此问题找到一个正确的提法：对变量作变换能否使信息增殖。

又是信息论经过并不冗长的数学处理^[18]使这个问题有了答案：

对变量作任何变换都不能使信息增殖

这是说无论我们对变量作任何逻辑的、数学的处理，它荷带关于另一变量的信息量至多是不减少，而不可能加多。

从一定角度看，上述结论与质量、能量定律几乎是三足鼎立的。图5.9就表现了这种统一性。

图5.9 信息、质量、能量定律的内在统一

这使“信息”完全跻身子与质量、能量完全同等的哲学地位。

这个信息的不增殖原理及在气象预告研究中的巨大价值．在文献[19]中作了充分讨论，这里仅能结合它与熵的关系简略介绍一下。

在信息论中^[20]，除了熵H(x)表示x的不确定性之外，还用H(X|Y)表示“条件熵”，其离散的表达式为

这里p(x|y)表示已知变量y恰为y_j，时x恰为x_i的概率。

变量y已知时所带来的关于x的信息Iy(X)由下式定义

(5.48)

这个式子非常直观：不确知部分的减少就是知识(信息)’的获得。可见信息就是两种熵的差值。

显然根据以上二式，在已知概率分布和条件概率分布之后就可以求出已知y以后会带来多少关于X的信息。

如果中国人的身高x的不确定性为log200，而已知体重后此不确定性变成log20，依上式就说体重知识带来了关于身高的信息，其值为log200-log20=log10=1Hartley(3.2比特)。

这3．2比特就是体重y带来的关于身高x的信息量。现在问，如果对变量y进行某种确定性的变换，使y变成新变量z。那么z带来的关于x的信息Iz(X)会大于原变量y带来的信息Iy(X)吗?

答案是:不。

这就是在变换中(y变成z)信息的不增殖原理。

气象预告员发现今天的气压值p与24小时以后的降水量r有一定关系，即Ip(r)>0，可这个信息量不够大。预告员希望提高信息值，为此他把p作了种种变换(把它开平方、取对数……，这类数学或逻辑花样是无穷尽的)使它成了新变量z，他指望这种数学努力，总可以找到一种变换方式，使新变量带来的信息要多于原变量。

可是简单、冷酷的信息不增殖原理说明变换后的新变量带来的信息至多是与原变量相等，而绝不会多一点!这可写成

Iz(X)≤Iy(X) (5．49)

过去(现在仍有!)曾有不少预告人员为作好天气预告在变量的数学变换上费尽气力，可从信息原理看，这是无效的劳动，我们应当把这种信息原理编入统计预告的教科书中，以免把人力浪费在这种无益活动中。

由于信息是两种熵的差，因而这个信息原理也作为熵原理的组成部分在此一并讨论。

这一信息原理与热力学第二定律还隐含着一层统一性，这就是它们都有非常近似的“可逆与不可逆”问题。

在热力学中的可逆过程是可以使工作介质与环境都从终止状态完全返回初始状态的可逆的热力过程，而不可逆过程则不能全部回到初始状态(详见热力学书)。

在数学中变量y与z的变换也有类似问题。如果每个变量y都能变成唯一的变量z，而且从z又能变回原来的y值，则说它们是互为单值的一一变换(如线性变换等)。如知道了变量z以后有时无法确知原来的y值(如四舍五入处理)则就不是互为单值的一一变换。

在集合论中，把y与z看成两个集合，上述变换对应于一一对应映射。

利用热力学第二定律可以导出，在可逆过程中系统的熵不加大，不可逆过程的则熵加大。

利用信息的不增殖原理可以导出在一一映射中信息不损失，在非一一映射中则信息减少。

一个是物理学定律，一个是数学的，但上述相似性十分引人注目，这内中含有更深层的原因。

限于篇幅，我们不能对此简单又深刻的原理作进一步讨论了。在笔者看来这个原理在科技领域的宣传与应用都有大量的工作尚待人去做。

§5 熵原理讨论

以上4节我们从不同侧面对可以归入熵原理的若干原理作了介绍，这4节的内容有一定的独立性，又有一定的同一性。现在应当讨论一下究竟什么是熵原理了。

热力学第二定律是个物理定律，它可以通俗地表述为热量仅能自发地从高温传向低温，其他的结论都可从这一经验性的事实出发一步步地推证出来。熵自动地达到它的极大值就是其一个重要推论。

信息论中的熵极大原理，从形式上看它几乎仅是一种假设，而不是个原理。这个熵极大假设是十分有用的，它常可帮我们找出多种分布函数。在某些统计问题中，它好象是个很好的“估计理论”。而在某些自然现象中(这涉及面十分宽)如果约束不足以完全决定事物的状态，那么事物自身容有的自由会使“某种”混乱性达到最大的程度，这也就是人们从物理上判知该系统内达到最混乱、最任意、最复杂的程度。这种物理判断讲的抽象一点就是熵极大。而在一定约束条件下，利用熵极大这一假设推求出来的分布函数总是与事实相符，这就导致我们把假设改称原理了。

这样似乎是说熵极大原理用到哪里都正确(如果约束取的妥当、并归纳成分布问题、而且数学求解无误)，可是对上述物理思路认识不深的人又责怪说，一切事物都遵守熵极大原理似乎太抽象、太神秘了。

统计力学思路一节是想消除上述神秘性。在那里撇开熵极大原理而经过较细的物理和数学推证也达到了几乎完全一致的结论。在那里熵极大不是开始时给定的假设，而是推论的结果，它说明了熵极大的含义，可这却使熵极大失去了作为公理的地位。

在推证熵极大与统计力学思路的一致性时，我们是否又用了另外的假设和公理呢?

有。确实在那里用了另外两个前提。一个是假设各个基本状态的出现机会是相等的(如每个雨量筒收集到多少降水的机率是相同的)，这正是统计力学中的一个基本假设。这是个假设，看来也仅能作为一个假设，它符合实际(从推论结果与事实相符来证明)又合情理，它简单，又是最易于想到的假设。

另一个前提就是概率高的事件几乎是必然要出现的事件，这一点在图5.8中显示的很清楚。图中的a，b两种降水分布不是不能出现，但是与d种降水分布比起来，a、b的实现办法仅有1、100种，而d的分布有10⁹⁷种方案都与之对应。一个袋子中有10⁹⁷个白球而仅有1或者100个黑球，任取一球，则抽中白球的概率就几乎是1，即这几乎是必然事件了。

统计力学思路一节的上述两点会推导出与熵极大一致的结论来，这又促使人们思考熵极大的本质。我觉得这里最核-心的思想是第二点，即高概率事件几乎是必然要出现的事件(分布)。

我们知道在统计数学中还有一个最大似然原理，依靠这个原理可以较准确的去估计已知概率分布中的未知参数。最大似然原理的核心思想就是高概率的事件最易于出现，这导致人们定义各概率的连乘积为似然函数，令此函数最大就成了最大似然原理。从似然函数达极大,人们就反求出各统计参数的最大似然估计值。而在一些文献中^{[21, 22]}认为，熵极大与最大似然是一致的，这确实也是十分值得重视的思路。

最大似然原理易于理解，在统计数学中应用了多年，认同它与熵原理的一致性对统一人类的知识当然是有益的。

在笔者已往的有限记忆中，似乎抽象科学(如数学、逻辑)中的原理并不依赖于实验；而物质科学(物理、化学、生物)中的基本原理都依赖于实验；因而尽管我们经常混同地使用抽象科学和物质科学中的规律，但它们的出身不同仍然是不言自明的。可是在讨论什么是熵原理时，它的出身似乎是模糊起来了，它是物理学的定律又似是抽象的数学定律，这种两栖身份在科学定律中显得很特殊。

由于篇幅限制，我们对熵原理的议论就要结束了。我得承认尽管对熵原理作了说明，但尚未给读者一个统一的表述。我认为这是科学界至今没有解决的大事，今天应当首先承认这个问题存在，才会促进人们去解决这个悬了两个世纪的问题。

笔者相信，在各领域广泛应用熵原理的各个组成部分，会促进整个熵原理的科学归纳，而横跨数学与物理，科学与技术、自然与社会的熵原理[2]的再发现会把人类文明提高到一个新的水平。

下一章我们就把问题收回到它如何用于气象学的问题这种较具体事物上去了。

附录 logx≤(x-1)loge的证明

令f(x)=x-1-lnx其中x>0。现分析函数f(x)有无

极值存在，为此将f对x求导数，得

f'(x)=1-(1/x)

令f’(x)=0，得x=1。这是f(x)的仅有的一个极值出现的

位置。再对f’(x)求导数

f’’(x)=1/x²

由此知在x定义域内f’’(x)>0，即f(x)应当在x=1处有一个极小值，此值为／(1)二0。换言之，应当有

f(x)≥0

或 x-1-1nx≥0

即 1nx≤x-1

这就证明了最初的不等式，此式仅在x=1时等号成立。

再者，利用对数换底的关系可将上式变成

logx≤(x-1)loge

此处e为自然数(e=2.71828．．)。

参 考 文 献

[1] 马本昆等，热力学与统计物理，人民教育出版社，30，(1980)。

[2] 唐有祺，统计力学及其在物理化学中的应用，科学出版社，155--159,(1979)。

[3] 赤池弘次，统计与熵，数学译林，(2)141--153，(1984)。

[4] E．T．Jaynes，Information Theory and Statistical Mechanics，Physical Review，106，4，620--630(1957)

[5] E．T．Jaynes，Where do we stand on Maximum Entroy? The Maximum Entropy FormaIism Massa．Inst．Tech.，second print, 15-118，(1981)。

[6] E．T．Jaynes，Clearing up mysteries—The original goal，Maximum Entropy and Bayesian Methods, kluwer Academic Publishers, 1-27，(1989) 。

[7] M．Tribus，Thirty years of lnformation Theory，The Maximum Entropy Formalism, Maesa．Inst．Tech.， second printing，1--14，(1981)。

[8] 祖巴列夫，非平衡态统计热力学，李沅柏等译，高等教育出版社，94--98，(1982)．

[9] H．哈肯，信息与自组织，四川教育出版社，93--100，(1988)。

[10] 廖胜清，最大熵原理及其对统计力学的应用，熵与交叉科学，30--35，(1988)。

[11] G．Silviu，Infoemation Theory with Application, McGRAW-HILL．CO., 293--30l，(1977)．

[12] F.M. Reza, An Introduction To Information Theory，278—282，McGRAW—HiLL．CO.,(1961)．

[13] 近藤次郎，数学模型，张春茹等译，机械工业出版社， 413--414，(1985)。

[14] R．S．Ingarden Information Theory and Variational

Principles in Statistical Theories，Bull Acad. Polon. Sci．Seral．Sci. Math. Astronom．Phys.,11，541—547，(1963)．

[15] S．Kullback and R．A．Leibler，On Information and sufficiency, Ann．Math．statist.，22，79--86，(1951)．

[16] 同(2)，30—34。

[17] F．Mandl，统计物理学，范印哲译，人民教育出版社，66—67，(1981)。

[18] 张学文，气象预告问题的信息分析，科学出版社，24--27，(1981)。

[19] 同[18]，75--90。

[20] 同[18]，15-21。

[21] 同[2]。

[22] 王文军，统计信息熵与充分性度量，熵与交叉科学，气象出版社，65--69，(1988)。

（本稿开始页码与原书相同，结束页码比原书少1页

2007-7-21完成。编辑：张学文，图：胡占永，张学文）

第四章 谱分析

§4.1谱分析的定义

谱分析是一种将模态分析的结果与一个已知的谱联系起来计算模型的位移和应力的分析技术。谱分析替代时间－历程分析，主要用于确定结构对随机载荷或随时间变化载荷(如地震、风载、海洋波浪、喷气发动机推力、火箭发动机振动等等)的动力响应情况。

§4.2什么是谱

谱是谱值与频率的关系曲线，它反映了时间－历程载荷的强度和频率信息。ANSYS的谱分析有三种类型：

·响应谱分析

Ø单点响应谱（Single-point　Response　Spectrum，SPRS）

Ø多点响应谱（Multi-point　Response　Spectrum，MPRS）

·动力设计分析方法（Dynamic　Design　Analysis　Method，DDAM）

·功率谱密度（Power　Spectral　Density，PSD）

在ANSYS/Professional产品中只提供单点响应谱方法。

§4.2.1响应谱分析

一个响应谱代表单自由度系统对一个时间－历程载荷函数的响应，它是一个响应与频率的关系曲线，其中响应可以是位移、速度、加速度、力等。响应谱又分为如下两种形式：

§4.2.1.1单点响应谱

在模型的一个点集上定义一条(或一族)响应谱曲线，例如在所有支撑处，图4-1（a）所示。

ANSYS/LinearPlus　program中只能进行单点响应谱分析。

§4.2.1.2多点响应谱

在模型的不同点集上定义不同的响应谱曲线，图4-1（b）所示。

图4-1单点响应谱和多点响应谱

§4.2.2动力设计分析方法

该法是一种用于分析船用装备抗振性的技术，它所使用的谱是从美国海军研究实验室报告(NRL－1396)中一系列经验公式和振动设计表得来的。

§4.2.3功率谱密度

功率谱密度谱是一种概率统计方法，是对随机变量均方值的量度。一般用于随机振动分析，连续瞬态响应只能通过概率分布函数进行描述，即出现某水平响应所对应的概率。

功率谱密度是结构在随机动态载荷激励下响应的统计结果，是一条功率谱密度值—频率值的关系曲线，其中功率谱密度可以是位移功率谱密度、速度功率谱密度、加速度功率谱密度、力功率谱密度等形式。数学上，功率谱密度值—频率值的关系曲线下的面积就是方差，即响应标准偏差的平方值。

与响应谱分析相似，随机振动分析也可以是单点的或多点的。在单点随机振动分析时，要求在结构的一个点集上指定一个功率谱密度谱；在多点随机振动分析时，则要求在模型的不同点集上指定不同的功率谱密度谱。

§4.2.4确定性分析与概率分析

响应谱和动力设计分析方法都是定量分析技术，因为分析的输入输出数据都是实际的最大值。但是，随机振动分析（功率谱分析PSD）是一种定性分析技术，分析的输入输出数据都只代表它们在确定概率下的可能性发生水平。

§4.3谱分析使用的命令

建立有限元模型和执行谱分析所使用的命令与其它有限元分析完全一样。同样，无论进行那种分析都可选用相似的GUI操作进行建模和求解。

本章节后的“谱分析例题（GUI交互方法和命令与批处理方法）”，讲述在GUI和命令环境下进行谱分析过程。如果要更详细地了解ANSYS的命令，参看《ANSYS命令参考手册》。

下面将详细地探讨两种常用的谱分析方法---单点响应谱（SPRS）和随机振动（PSD）。动力设计方法（DDAM）和多点响应谱（　MPRS）分析将简单地讨论与前两种分析步骤的不同点。

§4.4单点响应谱(SPRS)分析步骤

单点响应谱分析有如下五个步骤：

1.建造模型；

2.获得模态解；

3.获得谱解；

4.扩展模态；

5.合并模态；

6.观察结果。

结构的振型和固有频率是谱分析所必须的数据，因此要先进行模态分析。另外，在扩展模态时，只需扩展到对最后进行谱分析有影响的模态就可以。

§4.4.1建造模型

该步与其它分析类型建立模型的过程相似，即定义工作名、分析的标题、单元类型、单元实常数、材料性质、模型几何形状等。注意以下两点：

·只有线性行为在谱分析中才是有效的。任何非线性单元均作为线性处理。如果含有接触单元，那么它们的刚度始终是初始刚度，不再改变；

·必须定义材料弹性模量（EX）（或其他形式的刚度）和密度（DENS）。材料的任何非线性将被忽略，但允许材料特性是线性的、各向同性或各向异性以及随温度变化或不随温度变化。

§4.4.2获得模态解

结构的模态解（固有频率和振型）是计算谱解所必须的。模态分析的具体过程在《模态分析》中已经阐述过，这里还需注意以下几点：

·使用Block　Lanczos法（缺省）、子空间法或缩减法提取模态。非对称法、阻尼法、QR阻尼法以及PowerDynamics法对下一步谱分析是无效的；

·所提取的模态数目应足以表征在感兴趣的频率范围内结构所具有的响应；

·如果使用GUI交互式方法进行分析，模态分析设置[MODOPT]对话框的扩展模态选项置为NO状态，那么模态计算时将不进行模态扩展，但是可以选择地扩展模态（参看MXPAND命令的SIGNIF输入项的用法）。否则，将扩展模态选项置为YES状态。

·材料相关阻尼必须在模态分析中进行指定；

·必须在施加激励谱的位置添加自由度约束；

·求解结束后退出SOLUTION处理器。

§4.4.3获得谱解

1.进入求解器　(/SOLU命令)

命令：/SOLU

GUI　:Main　Menu　>　Solution

2.定义分析类型及分析选项

ANSYS提供下表所示谱分析选项：

表　4.1分析类型及分析选项

选项	命令	GUI 路径
新分析	^ANTYPE	^Main^Menu　>　Solution　>-Analysis　Type-New　Analysis
分析类型:谱分析	^ANTYPE	^Main^Menu　>　Solution　>-Analysis　Type-New　Analysis
谱分析类型:　SPRS	^SPOPT	^Main^Menu^>^Solution^>^Analysis^Option
求解所需的扩展模态数	^SPOPT	^Main^Menu^>^Solution^>^Analysis^Option

·选择新分析(ANTYPE命令)

只能选择新分析(New　Analysis)。

·选择分析类型(ANTYPE命令)

选择谱分析作为分析类型。

·选择谱分析类型(SPOPT命令)

选择单点响应谱(SPRS)分析。

·定义求解所需的扩展模态数(SPOPT命令)

选择足够的扩展模态数，使足以覆盖谱所决定的频率范围，并足够表征结构的响应特性。求解的精度取决于所使用的模态数：数值越大，精度越高。如果要计算单元应力，　SPOPT命令选项必须置为YES状态。

3.定义载荷步选项。对单点响应谱分析而言，下面的选项是有效的：

表4.2载荷步选项

选项	命令	GUI 路径
谱选项
响应谱类型	SVTYP	^Main^Menu　>　Solution　>　-Load　Step　Opts-Spectrum　>　-Single-Point-Settings
激励方向	SED	^Main^Menu　>　Solution　>　-Load　Step　Opts-Spectrum　>　-Single-Point-Settings
谱值-频率曲线	FREQ,SV	^Main^Menu　>　Solution　>　-Load　Step　Opts-Spectrum　>　-Single-Point-Freq　Table/Spectr　Values
阻尼（动力选项）
b(刚度)阻尼	BETAD	^Main^Menu　>　Solution　>　-Load　Step　Opts-Time/Frequenc　>　Damping
恒定阻尼比	BETAD	^Main^Menu　>　Solution　>　-Load　Step　Opts-Time/Frequenc　>　Damping
模态阻尼	MDAMP	^Main^Menu　>　Solution　>　-Load　Step　Opts-Time/Frequenc　>　Damping

·响应谱的类型[SVTYPE]

谱的类型可以是位移、速度、加速度、力或PSD。除了力谱以外，所有谱都是地震谱，亦即它们都是假定作用于结构的基础上。力谱用F或FK命令定义于非基础节点上，方向通过FX、FY和FZ方向指定。PSD谱[SVTYPE,4]在内部被转换成位移响应谱并被限定为平面窄带谱。更完整的随机振动分析参见“随机振动分析”章节。

·激励方向[SED]

·谱值与频率的关系曲线[FREQ,SV]

SV和FREQ命令用于定义谱曲线，可定义一系列不同阻尼比的谱曲线。用STAT命令列表显示当前谱曲线的对应值。还可以用ROCK命令定义摆动谱。

·阻尼选项(动力分析选项)

如果指定多种阻尼，ANSYS程序将计算出对应每个频率的有效阻尼比。接着，对谱曲线取对数计算出与该有效阻尼比对应的谱值。如果不指定任何阻尼，程序将自动选用阻尼最低的谱曲线。有关各种阻尼的详细介绍参见《瞬态动力学分析》“阻尼”一节。

可以选用的阻尼有如下几种：

1)b(刚度)阻尼[BETAD]

定义频率相关阻尼比。

2)恒定阻尼比[DMPRAT]

指定在所有频率上具有。恒定的阻尼比

3)模态阻尼[MDAMP]

注意：材料相关阻尼比(MP,　DAMP 命令) 可以选用，但只局限于模态分析。MP,　DAMP 命令还可以指定材料相关的恒定阻尼比（不是其它分析中选用的材料相关刚度阻尼）。

4.开始求解计算

Command:　SOLVE

GUI:Main　Menu　>　Solution　>　-Solve-Current　LS

求解的输出结果包括参与系数表。作为打印输出的一部分，参与系数表列出了参与系数、(基于最低阻尼比的)模态系数以及每阶模态的质量分布。模态系数乘以振型就是每阶模态的最大响应。用*GET命令提取模态系数后，将其作为SET命令中的一个比例系数来完成这个过程。

5.若要获得更多的响应谱解，只要重复3、4步即可。注意：这时的求解结果不要写进原来的结果文件。

6.退出求解器

Command:　FINISH

GUI:退出求解器。

§4.4.4扩展模态

1.将“Expansion　Pass　dialog　box”对话框中的扩展模态选项（expansion　pass　option）设置为YES。

Command:　EXPAND

GUI:Main　Menu　>　Solution　>　New　Analysis-Modal-Expansion　pass-On　>　L.S.Opt-Expansion　Pass　>　Single　Expand-Modal-Expand　Modes

2．无论使用Block　Lanczos法、子空间法还是缩减法，都必须进行模态扩展。关于模态扩展，在《动力学分析指南—模态分析》部分中作为独立的求解过程有详细讲述。另外，还需注意以下几点：

·只选择重要的模态进行扩展（参见MXPAND命令的SIGNIF项的用法）。如果使用GUI交互方法进行模态扩展操作，那么在模态分析阶段应将modal　analysis　options对话框[MODOPT]中的mode　expansion选项设置为NO，以便把模态扩展作为一个独立的求解过程，并放在谱分析完成之后进行；

·只有扩展模态才能在以后的模态合并过程中进行模态合并操作；

·如果对谱所产生的应力感兴趣，这时必须进行应力计算。在缺省情况下，模态扩展过程是不包含应力计算的，同时意味着谱分析将不包含应力结果数据。

·如果需要扩展所有模态，只要在模态求解过程中执行MXPAND命令，就同时进行模态扩展过程。如果使用GUI交互方法，并想扩展所有模态，那么在模态分析阶段应将modal　analysis　options对话框（MODOPT命令）中的mode　expansion选项设置为YES。如果只想扩展有明显意义的模态，就必须将模态扩展作为一个独立求解过程放在谱分析之后进行。

注意：即使进行模态扩展分析，模态分析解都将写进结果文件（Jobname.RST ）

§4.5随机振动(PSD)分析步骤

PSD分析包括如下六个步骤：

1．建造模型；

2．求得模态解；

3．扩展模态；

4．获得谱解；

5．合并模态；

6．观察结果。

以上六步中，前两步跟单点响应谱分析一样，后四步将在下面作详细讲解。ANSYS/Professional产品中不能进行随机振动分析。

如果选用GUI交互方法进行分析，模态分析选择对话框（MODOPT命令）中包含有是否进行模态扩展选项（MXPAND命令），将其设置为YES就可以进行下面的：扩展模态。这样，第二步（求得模态解）和第三步（扩展模态）就合并到一个步骤中进行计算。

§4.5.1扩展模态

无论选用子空间法、Block　Lanczos法还是缩减法，都必须进行模态扩展。关于模态扩展，《动力学分析指南—模态分析》部分“扩展模态”一节有详细讲述。另外还需注意以下几点：

·只有扩展后的模态才能在以后的模态合并过程中进行模态合并操作；

·如果对谱所产生的应力感兴趣，这时必须进行应力计算。在缺省情况下，模态扩展过程是不包含应力计算的，这同时意味着谱分析将不包含应力结果数据。

·模态扩展可以作为一个独立的求解过程，也可以放在模态分析阶段；

·在模态扩展结束之后，应执行FINISH命令退出求解器（SOLUTION）。

正如《动力学分析指南—模态分析》部分中讲述的那样，在进行模态分析时执行MXPAND命令就可以将模态求解和模态扩展合并成一步（GUI交互方法和批处理方法）。

§4.6随机振动分析结果应用

§4.6.1随机振动结果与失效计算

随机振动的计算结果（rms,average　frequence等）是各种失效计算的基础。这里仅仅就如何利用随机振动分析结果作出一些讨论，并不代表所有。假定所有随机过程为正态（高斯）平稳随机过程，确定失效的重要统计参数是均值(mean)、均方根（rms）和平均频率（average　frequency）。在ANSYS中，PSD分析均假定均值为零。

失效一般分成两类：一类是可逆的（随着激励的撤消而失效消失），如真空管、激光瞄准器和水晶激发器；第二类是不可逆的（随着激励的撤消而失效仍然存在），如结构件的断裂和机器零部件的疲劳破坏。一般情况下，振动越厉害，失效发生的可能性就越大。

在随机工作状态下，至少存在三种失效方式：（1）一次失效，即研究对象的某物理量在工作时的数值第一次达到确定水平时就会出现失效。例如，当应力达到一定水平时，结构可能发生延展失效。（2）时间失效，即研究对象的某物理量在工作时的数值超过预定水平的一定（寿命）时间百分比之后就会出现失效。这种方法通常用在电器元件，是一种不可逆的行为。（3）累积损伤失效，即每次数值达到某一水平产生微小但可以定义的损伤，所有次损伤累积起来直至发生失效。这就叫疲劳破坏，是随机振动分析最常用的方法。下面仅就随机疲劳失效进行讨论

§4.6.2随机疲劳失效

每次随机振动循环都对累积损伤具有贡献，当总的累积损伤达到100%时就表示发生了失效，并且随机疲劳和确定性疲劳分析在理论上是相同的。对于恒定应力幅的周期载荷作用下的应力历程，波动应力的每次循环都在材料中产生小的变化，产生微小裂缝，并不断发展。应力幅越大，每次循环裂缝增长的越快越多。随着应力循环次数的不断增加，累积损伤不断增多，直至发生失效。对于恒定应力幅，疲劳失效允许的循环次数是按材料的疲劳（S-N）曲线进行确定。一般工程材料的S-N曲线是对数S-N曲线，在坐标图上近似为一条直线，近似的数学公式为：

NS ^b =c

其中：

N=应力幅对应的最大循环次数；

S=应力幅；

b=对于大多数材料为（5~20）的常数；

c=取决于材料的常数。

当应力历程是随机过程而不是固定幅值时，疲劳计算相对比较复杂困难，但已经有好多种处理方法，这里将介绍其中的两种方法：（1）Steinberg提出的比较简化方法，即三区间法，它广泛用于航天电子工业中；（2）Crandall和Mark提出的更精确的方法，即基于随机振动结果的Miner方法。象Miner定律的线性累积损伤假设一样，这两种方法都是基于固定振幅，假定应力幅循环了n次循环时，消耗了材料疲劳寿命的n/N部分，而其它应力水平的循环也以相同方式对材料产生部分损伤。如果应力幅S _i 的循环次数n _i 没有达到对应的许可次数N _i ，那么产生的累积损伤为：

当D=1时，表示疲劳寿命已经耗尽，预测发生了疲劳破坏。考虑到该假设的不精确性和缺点，累积损伤D经常假定为一个小于1的一个数就表示疲劳已经发生。

Crandall和Mark提出的基于随机振动结果的Miner方法的假设前提是随机过程是一个典型的窄带随机过程。随机过程中应力幅每次循环大小不一致，因而对总体累积损伤的贡献也就大小不一致。假定随机过程的时刻T对应的总体累积损伤为D(T)，当T=T _F 时总体损伤达到1并预测疲劳发生。

假定 表示单位时间内正零交点的平均数目（即频率），那么在时刻T的期间内循环的数目为 ，对应的总体损伤期望值为：

其中：

P(a)da=幅值在a和a+da间的循环（部分）次数；

p(a)=定义应力幅值的峰值概率密度；

N(a)　=幅值a所对应的允许循环次数。

期望的损伤是基于峰值的分布p(a)。对于正态分布数据，峰值满足瑞利（Rayleigh）分布，这样将上式中p(a)换成瑞利密度函数，利用S-N关系式替换N(a)，那么总体损伤期望值变成：

其中：

=χ（gamma）函数；

b,c=前面讲述的疲劳曲线公式中的常数。

另外，Steinberg提出了更加简单的基于高斯分布的三区间法，它表示：

68.3%的时间应力值在-1σ~+1σ之间

95.4%的时间应力值在-2σ~+2σ之间

99.73%的时间应力值在-3σ~+3σ之间

因而，在利用Miner定律进行疲劳计算，可以将应力可以处理成3个水平：

应力区间 发生的时间

-1σ~+1σ68.3%的时间

-2σ~+2σ27.1%的时间

-3σ~+3σ 4.33% 的时间

99.73%

该方法的前提是，大于3σ的应力仅仅发生在(100-99.73)0.27%的时间内，假定它们不造成任何损伤。这样，利用Miner定律进行疲劳计算，总体损伤的计算公式是：

其中：

n _{1 σ} =等于或低于1σ水平的实际循环数目（0.6831 ）；

n _{2 σ} =等于或低于2σ水平的实际循环数目（0.271 ）；

n _{3 σ} =等于或低于3σ水平的实际循环数目（0.0433 ）；

N _{1 σ} ,　N _{2 σ} ,　N _{3 σ} =根据疲劳曲线查得的1σ、2σ和3σ应力水平分别对应许可循环的次数。

利用1σ、2σ和3σ应力和统计平均频率计算随机疲劳是一个有效的过程。注意，统计平均频率等于载荷步4(1　sigma速度)除以载荷步3(1　sigma位移)结果的商。这样随机疲劳计算的一般过程是：

(1)计算感兴趣应力分量的统计平均频率(应力速度/应力)；

(2)假定68%的时间处于1σ水平，(95.73-68)27.45%的时间处于2σ水平，(99.73-95.45)4.33%的时间处于3σ水平；

(3)基于期望(工作)寿命和统计平均频率，计算1σ，2σ和3σ水平下的循环次数；

(4)基于S-N曲线计算疲劳寿命使用系数。

注意，上面所述只是为了说明如何利用ANSYS随机振动结果进行随机疲劳计算，仅仅介绍部分应用途径，用户还可以根据其它随机振动疲劳计算理论进行计算，这里不一一罗列。在ANSYS中不能直接实现上述计算过程。

§4.7　DDAM（动力设计分析方法）谱分析

除了以下五点外，DDAM谱分析与单点响应谱分析是一样的。

·所有的输入数据（几何尺寸、材料性质、单元实常数等）必须采用英制单位(英寸(不是英尺)、磅等)。

·选择DDAM而不是SPRS作为谱分析类型（SPOPT命令）。

·使用ADDAM和VDDAM命令，而不是使用SVTYPE、SV和FREQ等命令来定义谱值及其类型。使用SED命令指定激励的总方向。使用ADDAM和VDDAM命令定义计算系数，程序根据《理论部分》17.7.4节的经验方程来计算模态系数。

·NRL求和法（NRLSUM命令)是最适用的模态合并法。模态合并处理方法与单点响应谱分析时是一样的。模态合并要求指定阻尼。

·执行ADDAM和VDDAM命令时已经指定阻尼，在求解中也就毋需定义阻尼。如果定义阻尼，只将其用于模态合并中，在求解过程中将被忽略。

同单点响应谱分析一样，DDAM谱分析也是按照六大步骤进行系统地分析。如果采用批处理方法分析，注意以下两点：

·模态求解和求得DDAM谱解过程可以合并到模态分析求解（ANTYPE，MODAL）过程中。DDAM谱载荷（ADDAM,VDDAM和SED命令）后，

·模态合并命令后，上面的第四步和第五步（即模态扩展和模态合并求解过程）可以合并到模态分析求解（ANTYPE，MODAL且EXPASS，ON）过程中去。

ANSYS/Professional产品中不能进行DDAM谱分析。

§4.8多点响应谱（MPRS)分析

除了以下六点外，多点响应谱分析与随机振动（PSD）分析（包括文件要求）的过程是相同的。

·选用MPRS而不是PSD作为谱分析类型（SPOPT命令）。

·此时，“PSD—频率”关系表对应谱值—频率关系。

·各谱之间不能定义任何程度的相关性（即假定各谱之间是不相关的）。

·只计算（相对于基础激励的）相对结果，不计算绝对结果。

·除了PSDCOM的模态合并方法以外，其他所有模态合并方法都可以选用。

·多点响应谱分析的结果是以POST1的命令格式写入到模态合并文件（　Jobname.MCOM）中。这些命令依据（求解器中模态合并命令指定的）某种方式合并最大的模态响应，最终计算出结构的总响应。总响应包括总位移，如果在模态扩展过程中要求有这些结果，那么还包括总应力、总应变和总的反作用力。在求解器中，执行模态合并命令（即SRSS、CQC、GRP、DSUM和NRLSUM命令）时Label项如果设置成VELO或ACEL，对应的速度和加速度响应也将写入到模态合并文件中。

ANSYS/Professional产品中不能进行DDAM谱分析。

§4.9谱分析的实例（GUI命令流和批处理）

§4.9.3随机振动和随机疲劳分析算例

分析§4.7.1中的梁板结构在基础位移激励谱（如右表所示）作用下结构的响应情况。

谱表

激励谱
频率(Hz）	位移(´10^3 m)
0.5	0.01
1.0	0.02
2.4	0.03
3.8	0.02
17	0.005
18	0.01
20	0.015
32	0.01

§4.9.3.3随机疲劳计算

根据上述计算，可以得到载荷步(Load　Step　)3的第1个子步的1σ的最大Von　Mises应力为66.7Mpa，位置为节点27。按照Steinberg提出了基于高斯分布的三区间法，利用Miner定律进行疲劳计算。

假设结构振动时间T＝3e5，振动平均频率 ＝15次，那么：

n _{1 σ} ＝　0.6831 ＝3.074e6

n _{2 σ} ＝　0.271 ＝1.220e6

n _{3 σ} ＝　0.0433 ＝0.195e6

根据<<机械工程材料性能数据手册>>（机械工业出版社,1994.12）Q235A级钢的P-S-N曲线(图2.1-1)，查得：

应力等于1σ=66.7　Mpa时，N _{1 σ} ＝+∝

应力等于2σ=133.4Mpa时，N _{2 σ} ＝5.4e6

应力等于3σ=200.1Mpa时，N _{3 σ} ＝3.7e5

将上述数值代入总体损伤的计算公式：

＝0.753　<　1.0

通过上述计算，证明设计满足疲劳要求。

第六章 谱分析 Spectral Analysis

到目前为止，时刻变量的数值一般都表示成为一系列随机扰动的函数形式，一般的模型形式为：

我们研究的重点在于，这个结构对不同时点和上的变量和的协方差具有什么样的启示。这种方法被称为在时间域(time domain)上分析时间序列的性质。

在本章中，我们讨论如何利用型如和的周期函数的加权组合来描述时间序列数值的方法，这里表示特定的频率，表示形式为：

上述分析的目的在于判断不同频率的周期在解释时间序列性质时所发挥的重要程度如何。如此方法被称为频域分析(frequency domain analysis)或者谱分析(spectral analysis)。我们将要看到，时域分析和频域分析之间不是相互排斥的，任何协方差平稳过程既有时域表示，也有频域表示，由一种表示可以描述的任何数据性质，都可以利用另一种表示来加以体现。对某些性质来说，时域表示可能简单一些；而对另外一些性质，可能频域表示更为简单。

§6.1 母体谱

我们首先介绍母体谱，然后讨论它的性质。

6.1.1 母体谱及性质

假设是一个具有均值的协方差平稳过程，第个自协方差为：

假设这些自协方差函数是绝对可加的，则自协方差生成函数为：

这里表示复变量。将上述函数除以，并将复数表示成为指数虚数形式，，则得到的结果(表达式)称为变量的母体谱：

注意到谱是的函数：给定任何特定的值和自协方差的序列，原则上都可以计算的数值。

利用De Moivre定理，我们可以将表示成为：

因此，谱函数可以等价地表示成为：

注意到对于协方差平稳过程而言，有：，因此上述谱函数化简为：

利用三角函数的奇偶性，可以得到：

假设自协方差序列是绝对可加的，则可以证明上述谱函数存在，并且是的实值、对称、连续函数。由于对任意，有：，因此是周期函数，如果我们知道了内的所有的值，我们可以获得任意时的值。

§6.2 不同过程下母体谱的计算

假设随机过程服从过程：

这里：

，，

根据前面关于过程自协方差生成函数的推导：

因此得到过程的母体谱为：

例如，对白噪声过程而言，，这时它的母体谱函数是常数：

下面我们考虑过程，

此时：，则母体谱为：

可以化简成为：

显然，当时，谱函数在内是的单调递减函数；当时，谱函数在内是的单调递增函数。

对过程而言，有：

这时只要，则有：，因此谱函数为：

该谱函数的性质为：当时，谱函数在内是的单调递增函数；当时，谱函数在内是的单调递减函数。

一般地，对过程而言：

则母体谱函数为：

如果移动平均和自回归算子多项式可以进行下述因式分解：

则母体谱函数可以表示为：

从母体谱函数中计算自协方差

如果我们知道了自协方差序列，原则上我们就可以计算出任意的谱函数的数值。反过来也是对的：如果对所有在内的，已知谱函数的数值，则对任意给定的整数k，我们也能够计算k阶自协方差。这意味着母体谱函数和自协方差序列包含着相同的信息。其中任何一个都无法为我们提供另外一个无法给出的推断。

下面的命题为从谱函数计算自协方差提供了一个有用的公式：

命题6.1 假设是绝对可加的自协方差序列，则母体谱函数与自协方差之间的关系为：

上述公式也可以等价地表示为：

利用上述谱公式，可以实现谱函数与自协方差函数之间的转换。

解释母体谱函数

假设，则利用命题6.1可以得到时间序列的方差，即，计算公式为：

根据定积分的几何意义，上式说明母体谱函数在区间内的面积就是，也就是过程的方差。

更一般的，由于谱函数是非负的，对任意，如果我们能够计算：

这个积分结果也是一个正的数值，可以解释为的方差中与频率的绝对值小于的成分相关的部分。注意到谱函数也是对称的，因此也可以表示为：

这个积分表示频率小于的随机成分对方差的贡献。

但是，频率小于的随机成分对方差的贡献意味着什么？为了探索这个问题，我们考虑更为特殊一些的时间序列模型：

这里和是零均值的随机变量，这意味着对所有时间t，有。进一步假设序列和是序列不相关和相互不相关的：

，

，对所有的j和k

这时的方差是：

因此，对这个过程来说，具有频率的周期成分对的方差的贡献部分是。如果频率是有顺序的：，则的方差中由频率小于或者等于的周期形成的部分是：。

这种情形下的k阶自协方差为：

因为过程的均值和自协方差函数都不是时间的函数，因此这个过程是协方差平稳过程。但是，可以验证此时的自协方差序列不是绝对可加的。

虽然在上述过程中，我们已经过程的方差分解为频率低于某种程度的周期成分的贡献，我们能够这样做的原因在于这个过程是比较特殊的。对于一般的情形，著名的谱表示定理(the spectral representation theorem)说明：任何协方差平稳过程都可以表示成为不同频率周期成分的和形式。

对任意给定的固定频率，我们定义随机变量和，并假设可以将一个具有绝对可加自协方差的协方差平稳过程表示为：

这里需要对随机变量和的相关性给出更为具体的假设，但是上述公式便是谱表示定理的一般形式。

§6.2 样本周期图 Sample Periodogram

对一个具有绝对可加自协方差的协方差平稳过程，我们已经定义在频率处的谱函数值为：

，

注意到母体谱是利用表示的，而表示的是母体的二阶矩性质。

给定由表示的T个样本，我们可以利用下述公式计算直到阶的样本自协方差：

，

对于给定的，我们可以获得母体谱密度对应的样本情形，我们称其为样本周期图：

样本周期图也可以表示成为如下形式：

类似地，我们可以证明样本周期图下的面积等于样本方差：

样本周期图也是关于原点对称的，因此也有：

更为重要的是，谱表示定理在样本情形也有类似的表示。我们将要说明，对于平稳过程的任意一个容量为的观测值序列，存在频率和系数，，使得期的值可以表示成为：

其中：

当时，与不相关；

当时，与不相关；

对于所有的和，与不相关。

的样本方差是，该方差中可以归因于频率为的周期成分的部分由样本周期图给出。

我们对样本容量是奇数的情形展开讨论上述谱表示模式。这时可以表示成为由个不同频率构成的周期函数，频率如下：

，，……，

因此最高频率为：

我们考虑基于常数项、正弦函数和余弦函数的线性回归：

将这个回归方程表示成为下述方式：

其中：，这是一个具有个解释变量的回归方程，因此解释变量与观测值是一样多的。我们将证明解释变量之间是线性无关的，这意味着基于回归的OLS估计具有惟一解。该回归方程的 系数具有显著的统计意义：表示中可以归因于频率的周期成分的那部分。这就是说，任意观测到的序列，它都可以利用上述周期函数形式表示，并且不同频率的周期成分对方差的贡献都可以在样本周期图中找到。

命题6.2 假设样本容量是奇数，定义，并设定，，假设解释变量为：

则有：

进一步，假设是任意个实数，则下述推断成立：

(a) 过程可以表示为：

这里：

，，

(b) 的样本方差可以表示为：

样本方差可以归因于频率为的周期成分的部分为。

(c) 的样本方差中可以归因于频率为的周期成分的部分还可以表示为：

其中是样本周期图在频率处的值。

上述结果说明，是对角矩阵，这意味着包含在向量中的向量之间是相互正交的。这个命题断言：任何奇数个观测到的时间序列可以表示成为一个常数加上具有个不同频率的个周期成分的加权和。当是偶数整数的时候，类似的结果也是成立的。因此，这个命题给出了类似谱表示定理的有限样本的类似情况。这个命题进一步表明了样本周期图的特征是将的方差按部分分解为不同频率的周期成分的贡献。

注意到解释的方差的频率都落在区间中。为什么不使用负的频率？假设数据确实是由上述过程的一种特殊情形生成的：

这里代表某个特殊的负频率，和是零均值的随机变量，利用三角函数的奇偶性，可以将表示为：

因此，利用上述式子无法从数据中识别数据是从正发频率还是负的频率生成的。这时一种简单的方式是假设数据是从具有正的频率中生成的。

为什么只考虑作为最大的频率呢？假设数据真的是从频率的周期函数中生成的，例如：

这时正弦和余弦函数的周期性质表明，上式可以表示成为：

因此，根据以前的讨论，具有频率的周期在观测值上等价于具有频率的周期。

注意到频率和周期之间的关系，频率对应的周期为。由于我们考虑的最高频率为，因此我们所观测到的能够自己重复的最短阶段是。如果，则周期是每阶段重复自己。但是，如果数据是整数阶段观测的，因此数据可以观测的时间间隔仍然是每4个阶段观测到，这对应着周期频率是。例如，函数和函数在整数的时间间隔上，它们的观测值是一致的。

命题6.2也为计算在频率()上的样本周期图的数值提供了方法。定义：

这里：

，

因此可以得到：

§6.3 估计总本谱 Estimating the Population Spectrum

上面我们介绍了母体谱的意义和性质，下面我们面对的问题是：获得了观测样本以后，如何估计母体谱函数？

样本周期图的大样本性质

一个显然的方法是利用样本周期图去估计母体谱函数。但是，这种方法具有显著的限制。假设对于无限移动平均过程而言：

这里系数是绝对可加的，是具有均值和方差的独立同分布序列，假设是如上定义的母体谱函数，且对所有的，都有。假设是如上定义的样本谱函数，Fuller (1976) 证明了，对和充分大的样本容量，样本周期图与母体谱函数之比的二倍具有下述渐近分布：

进一步，如果，也有：

并且上述两个渐近分布的随机变量是相互独立的。

注意到的均值等于自由度，因此有：

因为是母体数量，不是一个随机变量，因此上式也可以表示成为：

因此，对充分大的样本容量，样本周期函数为母体谱提供了一个渐近无偏估计。

母体谱的参数化估计

假设我们认为数据可以由模型表示：

这里是具有方差的白噪声。这时一个估计母体谱的出色方法是先利用前面介绍的极大似然估计估计参数，具有绝对可加自协方差的协方差平稳过程，我们已经定义在频率处的谱函数

§6.4 谱分析的应用 Uses of Spectral Analysis

我们利用美国制造业生产的数据来说明谱分析的应用。书中给出了联邦储备委员会的季节非调整的月度指数，从1947年1月至1989年11月。其中出现经济衰退的时候出现了生产的下降，大约持续一年左右。数据中出现了显著的季节成分，大约在7月出现下降，而在8月出现复苏。

图6.4给出了原始数据的样本周期图。这里显示的是的函数，这里。

---

[1] 在表5.3中的v1=1，它代表第一档次的降水的平均值，即0到2的平均值，而此处的v1=0是降水的下限。

[2] 1989年在全国第二届熵与交叉科学研讨会上（兰州），笔者在论文中以“最复杂原理”称呼它。（2007注：该文扩充为2003年笔者写成的《组成论》）