索杆结构由受拉索和受压杆组成，包括张拉整体、索穹顶和索桁等一系列结构形式^[1]。索杆结构具有重量轻、跨度大、可扩展性强等特点，在土木工程、航空航天、机器人等领域应用广泛^{[2 − 3]}。在结构服役过程中，轻量化和低阻尼的索杆结构在地震作用、风荷载、冲击荷载等外部激励作用下，会产生高水平的振动响应^[4]。因此，有必要采取有效的控制方法来降低结构响应，提高结构的抗振能力和防灾性能。

结构振动控制方法分为主动控制、半主动控制、被动控制和混合控制^{[5, 6]}。主动控制在结构中布置传感器监测结构状态，由控制器引导嵌入结构系统的作动器输入最优控制力，抵消结构系统在外部激励下的振动能量，可有效降低结构的振动水平^{[7 − 8]}。索杆结构具有良好的可扩展性，可通过嵌入作动器实现结构振动响应的主动控制^[9]。XU等^[10]采用线性二次型控制算法(Linear Quadratic Regulator, LQR)与瞬时最优控制算法对索穹顶结构进行了主动控制模拟，验证了张力结构主动控制的可行性。YANG等^[11]结合非线性自适应控制算法与线性H∞算法实现张拉整体结构主动控制。PENG等^[12]提出了分散模型预测控制算法，将张拉整体结构主动控制问题分解为一系列的线性互补问题，实现张拉整体结构的分级主动控制。这些控制算法根据特定性能指标优化或迭代求解最优控制力，存在计算时间较长的问题。同时，由于结构模型的复杂性、控制器设计参数的多样性以及算法鲁棒性设计需求等客观因素的存在，主动控制算法的计算成本日益增长^[13]。因此，降低主动控制算法的计算成本十分有必要。

降阶模型是缓解动力学求解资源限制的一种有效的解决方案，有望降低主动控制算法计算时间^{[14 − 15]}。在结构动力学求解常用的模型降阶方法包括模态综合法、准静态里兹向量法、Krylov子空间法、本征正交分解法(Proper Orthogonal Decomposition，POD)等。ZHANG等^[16]利用二阶Krylov子空间法提升了瞬态激励下的结构极端变形求解效率。AYOUB等^[17]将POD方法扩展到多层框架结构的非线性动力求解中，大幅降低计算成本。Wang等^[18]将准静态里兹向量法应用于瞬态动力学方程求解，优化方程迭代计算效率。其中，POD方法由奇异值分解、主成分分析和Karhunen-Loève分解组成，是一种强大而有效的特征提取、模型降阶和数据压缩技术^[19]。POD方法从数据中提取本征正交模态，关注系统的空间特征，生成投影子空间用于模型降阶^[20]。

基于POD方法的降阶模型方法已在结构动力响应计算中得到广泛应用。AZEEZ等^[21]采用POD方法从时间序列数据中提取正交模态，建立悬臂梁的降阶模型，重构悬臂梁的振动冲击响应。AZAM等^[22]研究了基于POD算法的动力结构系统降阶建模性能，结果表明POD算法优于标准模态分析方法。CORIGLIANO等^[23]基于POD方法提出了一种适用于高度非线性结构的模型降阶方法，使得求解弹塑性问题所需的计算时间减少约50%，同时保持了与完整模型相当的精度。由上述研究可以发现，基于POD方法的降阶模型在结构动力响应分析方面展现出了强大的能力。与此相比，在降低振动主动控制算法的计算成本方面，POD方法的巨大潜力还未被挖掘。

本文提出了基于POD方法的索杆结构模型降阶和POD-LQR振动控制器设计方法。推导了基于奇异值分解从结构位移、速度或加速度数据提取投影子空间的计算公式。基于能量精度指标建立了投影子空间阶数的选择标准，实现了降阶模型阶数的自动化选择。建立了控制器设计参数在全阶空间和投影子空间之间的转换关系，确保POD-LQR控制器的控制性能指标不受模型阶数变化的影响。基于降阶模型设计POD-LQR控制器的时长大幅降低。

1   基于本征正交分解的索杆结构模型降阶

1.1   索杆结构有限元模型建立

索杆结构在平衡状态下的动力响应受几何非线性的影响较小，围绕平衡构型的线性化动力学模型可有效表征索杆结构动力行为^{[23 − 24]}。在外部激励作用下，采用有限元方法分析索杆结构运动，并考虑结构阻尼影响，其运动方程可表示为：

式中：${{\boldsymbol{M}}} \in {\mathbb{R}^{n \times n}}$为质量矩阵；${{\boldsymbol{E}}} \in {\mathbb{R}^{n \times n}}$为阻尼矩阵；${{\boldsymbol{K}}} \in {\mathbb{R}^{n \times n}}$为刚度矩阵；${\ddot {\boldsymbol{a}}} \in {\mathbb{R}^n}$、$\dot {\boldsymbol{a}} \in {\mathbb{R}^n}$和${\boldsymbol{a}} \in {\mathbb{R}^n}$分别为全局节点加速度、速度和位移向量；${\boldsymbol{f}} \in {\mathbb{R}^n}$为外力向量；$n$为结构的自由度数。

基于主动控制算法计算作动器最优控制力输入，以实现索杆结构的振动主动控制。本文通过将控制力输入到结构相应节点上，对作动器进行仿真，得到索杆结构在外部动态激励和作动器输入控制力共同作用下的运动方程：

式中：${{{\boldsymbol{B}}}_{\rm{s}}} \in {\mathbb{R}^{n \times l}}$为作动器的位置矩阵；${\boldsymbol{u}} \in {\mathbb{R}^l}$为控制力输入向量；$l$为作动器的数量。

将式(2)所示的运动方程转换为状态空间形式：

式中：${{\boldsymbol{A}}_{\rm{c}}} = \left[ \begin{matrix} {0}&{{\boldsymbol{I}}} \\ { - {{{\boldsymbol{M}}}^{ - {1}}}{{\boldsymbol{K}}}}&{ - {{{\boldsymbol{M}}}^{ - {1}}}{{\boldsymbol{E}}}} \end{matrix} \right]$；${\boldsymbol{x}} = \left[ \begin{matrix} {\boldsymbol{a}} \\ {\dot {\boldsymbol{a}}} \end{matrix} \right]$；$\dot {\boldsymbol{x}} = \left[ \begin{matrix} {\dot {\boldsymbol{a}}} \\ {\ddot {\boldsymbol{a}}} \end{matrix} \right]$；${{{\boldsymbol{B}}}_{\rm{c}}} = \left[ \begin{matrix} {0} \\ {{{{\boldsymbol{M}}}^{ - 1}}{{{\boldsymbol{B}}}_{\rm{s}}}} \end{matrix} \right]$；${{{\boldsymbol{H}}}_{\rm{c}}} = \left[ \begin{matrix} {0} \\ {{{{\boldsymbol{M}}}^{ - 1}}} \end{matrix} \right]$；${{\boldsymbol{y}}}$为观测向量；${{\boldsymbol{C}}}$为观测矩阵，当采用全状态测量时，有${{\boldsymbol{C}}} = {{\boldsymbol{I}}}$，$N = 2n$是索杆结构状态空间模型的维数。

1.2   模型降阶

结构位移向量${\boldsymbol{a}}$、速度向量$\dot {\boldsymbol{a}}$和加速度向量${\ddot {\boldsymbol{a}}}$投影过程可表示为：

式中：${{{\boldsymbol{P}}}_{\rm{a}} } \in {\mathbb{R}^{n \times {r_{\rm{a}} }}}$为投影矩阵；${\tilde {\boldsymbol{a}}} \in {\mathbb{R}^{{r_{\rm{a}} }}}$、${\tilde {\dot {\boldsymbol{a}}}} \in {\mathbb{R}^{{r_{\rm{a}} }}}$和${\tilde {\ddot {\boldsymbol{a}}}} \in {\mathbb{R}^{{r_{\rm{a}} }}}$分别为投影子空间中的位移向量、速度向量和加速度向量；${r_{\rm{a}} }$为投影子空间的维数。由于投影矩阵${{{\boldsymbol{P}}}_{\rm{a}} }$仅与位置有关，因此对由位移向量${\boldsymbol{a}}$、速度向量$\dot {\boldsymbol{a}}$或加速度向量${\ddot {\boldsymbol{a}}}$组成的数据矩阵进行奇异值分解均可得到投影矩阵${{{\boldsymbol{P}}}_{\rm{a}} }$。本文以结构位移数据为例进行理论推导，根据索杆结构特点适配与扩展本征正交分解，获得索杆结构的降阶状态空间模型。结构位移响应数据的矩阵形式${{\boldsymbol{X}}_{\rm{a}} }$表示如下：

式中，$m$为所有位移向量的个数。

对位移响应矩阵${{\boldsymbol{X}}_{\rm{a}} }$进行奇异值分解：

式中：${{{\boldsymbol{U}}}_{\rm{a}} } \in {\mathbb{R}^{n \times n}}$；${{{\boldsymbol{\varSigma}} }_{\rm{a}} } \in {\mathbb{R}^{n \times n}}$；${{{\boldsymbol{V}}}_{\rm{a}} } \in {\mathbb{R}^{m \times n}}$；${{\tilde {\boldsymbol{U}}}_{\rm{a}} } \in {\mathbb{R}^{n \times {r_{\rm{a}} }}}$；${{\tilde {\boldsymbol{\varSigma}} }_{\rm{a}} } \in {\mathbb{R}^{{r_{\rm{a}} } \times {r_{\rm{a}} }}}$；${\tilde {\boldsymbol{V}}}_{\rm{a}} ^ \top \in {\mathbb{R}^{{r_{\rm{a}} } \times m}}$。左奇异矩阵${{{\boldsymbol{U}}}_{\rm{a}} }$集合了结构所有的本征正交模态，奇异值的平方${\boldsymbol{\varSigma}} _{{\rm{a}} ,ii}^2$表示相应本征正交模态的最大能量。因此，可以指定能量精度指标${p_1}$，根据以下条件自动选择投影子空间的维数${r_{\rm{a}} }$：

取${{{\boldsymbol{P}}}_{\rm{a}}} = {{\tilde {\boldsymbol{U}}}_{\rm{a}} }$，并将式(4)代入式(2)可得到索杆结构的降阶运动方程：

式中：$\tilde {\boldsymbol{M}} = {{\boldsymbol{P}}}_{\rm{a}}^ \top {\boldsymbol{M}}{{{\boldsymbol{P}}}_{\rm{a}}}$；${\tilde {\boldsymbol{E}}} = {{\boldsymbol{P}}}_{\rm{a}}^ \top {{\boldsymbol{E}}}{{{\boldsymbol{P}}}_{\rm{a}}}$；${\tilde {\boldsymbol{K}}} = {{\boldsymbol{P}}}_{\rm{a}}^ \top {{\boldsymbol{K}}}{{{\boldsymbol{P}}}_{\rm{a}}}$；${\tilde{\boldsymbol{ f}}} = {{\boldsymbol{P}}}_{\rm{a}}^ \top {\boldsymbol{f}}$；${{\tilde {\boldsymbol{B}}}_{\rm{s}}} = {{\boldsymbol{P}}}_{\rm{a}}^ \top {{{\boldsymbol{B}}}_{\rm{s}}}$。将式(8)转换为状态空间形式可得：

式中：${\tilde {\boldsymbol{A}}_{\rm{c}}} = \left[ \begin{matrix} {0}&{{\boldsymbol{I}}} \\ { - {{{\tilde {\boldsymbol{M}}}}^{ - {1}}}{\tilde {\boldsymbol{K}}}}&{ - {{{\tilde {\boldsymbol{M}}}}^{ - {1}}}{\tilde {\boldsymbol{E}}}} \end{matrix} \right]$；$\tilde {\boldsymbol{x}} = \left[ \begin{matrix} {{\tilde {\boldsymbol{a}}}} \\ {\tilde {\dot {\boldsymbol{a}}}} \end{matrix} \right]$；$\tilde {\dot {\boldsymbol{x}}} = \left[ \begin{matrix} {\tilde {\dot {\boldsymbol{a}}}} \\ {\tilde {\ddot {\boldsymbol{a}}}} \end{matrix} \right]$；${{\tilde {\boldsymbol{B}}}_{\rm{c}}} = \left[ \begin{matrix} {0} \\ {{{{\tilde {\boldsymbol{M}}}}^{ - 1}}{{{\tilde {\boldsymbol{B}}}}_{\rm{s}}}} \end{matrix} \right]$；${{\tilde {\boldsymbol{H}}}_{\rm{c}}} = \left[ \begin{matrix} {0} \\ {{{{\tilde {\boldsymbol{M}}}}^{ - 1}}} \end{matrix} \right]$。

根据式(3)得到降阶状态空间模型的另一种方法是取投影矩阵${{{\boldsymbol{P}}}_{{\rm{A}}} }$和投影关系如下：

将式(10)~式(11)代入式(3)，可得降阶后的状态空间模型，表示如下：

式中：${\tilde {\boldsymbol{A}}_{\rm{c}}} = {{\boldsymbol{P}}}_{{\rm{A}}} ^ \top {{\boldsymbol{A}}_{\rm{c}}}{{{\boldsymbol{P}}}_{{\rm{A}}} }$；${\tilde {\boldsymbol{B}}_{\rm{c}}} = {{\boldsymbol{P}}}_{{\rm{A}}} ^ \top {{\boldsymbol{B}}_{\rm{c}}}{{{\boldsymbol{P}}}_{{\rm{A}}} }$；${{\tilde {\boldsymbol{H}}}_{\rm{c}}} = {{\boldsymbol{P}}}_{{\rm{A}}} ^ \top {{{\boldsymbol{H}}}_{\rm{c}}}$。此时降阶后的状态空间模型的维数为${r_{{\rm{A}}} } = 2{r_{\rm{a}} }$。

2   LQR控制器设计

在外部激励作用下，系统偏离平衡状态，线性二次型最优控制的目标是在不消耗过多能量的情况下，施加控制力使系统趋近于平衡状态。对于线性结构系统，选取系统状态和控制输入的二次型函数作为控制性能指标函数^[25]。传统的FEM-LQR控制器直接根据式(3)所示的基于有限元方法得到的状态空间模型进行设计，本文提出的POD-LQR控制器根据式(12)所示的基于POD方法得到的降阶状态空间模型进行设计。其控制性能指标${Z_{{{\rm{L}}{\rm{Q}}{\rm{R}}} }}$和${\tilde Z_{{\rm{L}}{\rm{Q}}{\rm{R}}}}$定义为：

式中：${{\boldsymbol{x}}^ \top }{\boldsymbol{Q}}x$和${\tilde {\boldsymbol{x}}^ \top }\tilde {\boldsymbol{Q}}\tilde {\boldsymbol{x}}$为结构振动能量的大小，分别记为${Z_1}$和${\tilde Z_1}$；${ {\boldsymbol{u}} ^ \top }{{\boldsymbol{R}}u}$为控制能量输入的大小，记为${Z_2}$；${\boldsymbol{Q}}$为关于状态$x$的半正定加权矩阵；$\tilde {\boldsymbol{Q}}$为关于状态$\tilde {\boldsymbol{x}}$的加权矩阵；${{\boldsymbol{R}}}$为关于向量${\boldsymbol{u}}$的正定加权矩阵。

取$\tilde {\boldsymbol{Q}} = {{\boldsymbol{P}}}_{{\rm{A}}} ^ \top {\boldsymbol{Q}}{{{\boldsymbol{P}}}_{{\rm{A}}} }$可以保证POD-LQR控制器的控制性能指标与FEM-LQR控制器的控制性能指标相当，不受降阶模型阶数${r_{{\rm{A}}} }$取值的影响，证明如下：

LQR算法提供的最优控制力输入${\boldsymbol{u}}$、反馈增益矩阵${{{\boldsymbol{G}}}_{\rm{c}}}$和${{\tilde {\boldsymbol{G}}}_{\rm{c}}}$为：

式中，${{{\boldsymbol{S}}}_{\rm{c}}}$和${{\tilde {\boldsymbol{S}}}_{\rm{c}}}$为下式所示Riccati代数方程的解：

索杆结构的几何刚度是由预应力提供的，主动控制过程中应避免预应力松弛，预应力${\boldsymbol{f}}_{\rm{p}}^{\rm{e}}$和输入${\boldsymbol{u}}$应满足的约束关系为：

式中：${\boldsymbol{f}}_{\rm{p}}^{\rm{e}}$和${\boldsymbol{u}}$受拉时取正值，受压时取负值；$\beta$为安全系数，其值小于1。引入$\beta$是为了考虑环境荷载(如温度荷载或雪荷载)对预应力的影响，可依据环境荷载的统计数据来确定其值。

值得注意的是，使用基于本征正交分解方法的降阶状态空间模型显著节省了求解Riccati代数方程的计算时间，因为其计算复杂度从$O({N^3})$降低到$O(r_{{\rm{A}}} ^3)$(${r_{{\rm{A}}} } \ll N$)。因此，POD-LQR控制器的计算效率显著提高。FEM-LQR控制器和POD-LQR控制器设计流程如图1所示。

图  1  FEM-LQR和POD-LQR控制器设计流程图

Figure  1.  FEM-LQR and POD-LQR controller design flowchart

下载: 全尺寸图片 幻灯片

3   方法验证

将本文提出的POD-LQR方法应用于索杆结构振动主动控制，并与FEM-LQR方法进行对比，说明其在控制效果和控制器设计时间上的优势。具体来说，考虑了两种类型的索杆结构。一种是Levy索穹顶，利用位移响应数据获取降阶状态空间模型，设计POD-LQR控制器；另一种是30杆球形张拉整体，利用加速度响应数据获取降阶状态空间模型，设计POD-LQR控制器。

3.1   算例1：Levy索穹顶

1) 模型降阶和LQR控制器设计

Levy索穹顶结构(图2)由96根拉索和17根压杆组成，跨度为20 m，约束为周边8个节点固接，共113个单元根据对称性分为11组。具体结构参数如表1所示。采用二次奇异值分解方法确定Levy索穹顶结构整体可行预应力，预应力设置详见表2。基于有限元方法的索杆结构状态空间模型的阶数为$N = 204$。选择外圈的4个单元布置作动器，具体位置如图3所示。${{\boldsymbol{Q}}}$和R为LQR理论的结构主动控制系统设计的关键控制参数，增大${{\boldsymbol{Q}}}$和减小${{\boldsymbol{R}}}$均能提升LQR算法的控制效果。${{\boldsymbol{Q}}}$的取值设置为$\left[ \begin{matrix} {{{{\boldsymbol{I}}}_{102 \times 102}}}&{{{0}_{102 \times 102}}} \\ {{{0}_{102 \times 102}}}&{{{0}_{102 \times 102}}} \end{matrix} \right]$，在满足式(22)表示的索杆结构主动控制约束条件，${{\boldsymbol{R}}}$取$1 \times {10^{ - 14}}{{{\boldsymbol{I}}}_{4 \times 4}}$以实现目标控制效果。

图  2  Levy索穹顶

Figure  2.  Levy cable dome

下载: 全尺寸图片 幻灯片

表  1  Levy索穹顶结构参数

Table  1.  Structure parameters of Levy cable dome

结构参数	参数取值
总节点数	42
自由节点数	34
单元数	113
节点集中质量/kg	5
拉索弹性模量/GPa	181
压杆弹性模量/GPa	206
拉索等效截面面积/m2	6×10−4
压杆等效截面面积/m2	1×10−3
材料密度/(kg/m3)	7850
阻尼比/(%)	5
安全系数	0.95

下载: 导出CSV 

| 显示表格

表  2  Levy索穹顶结构参数

Table  2.  Prestress settings of Levy cable dome

单元分组编号	预应力/N		单元分组编号	预应力/N
1	7220		7	5839
2	3388		8	12129
3	6968		9	−11327
4	3286		10	−2636
5	16890		11	−5775
6	8807		—	—

下载: 导出CSV 

| 显示表格

图  3  Levy索穹顶作动器位置示意图(红色突出显示)

Figure  3.  Levy cable dome embedded with actuators (highlighted in red)

下载: 全尺寸图片 幻灯片

利用基于有限元法的Levy索穹顶模型生成结构位移响应数据，构建基于POD的降阶状态空间模型。将均值为0，方差为100的随机荷载输入到结构的所有节点上，使结构产生振动，结构位移数据采集时间为5 s，采样频率为100 Hz。将能量精度指标${p_1}$设为99.5%，得到基于POD的降阶状态空间模型阶数为${r_{{\rm{A}}} } = 36$。矩阵$\tilde {\boldsymbol{Q}}$取为$\tilde {\boldsymbol{Q}} = {{\boldsymbol{P}}}_{{\rm{A}}} ^ \top {\boldsymbol{Q}}{{{\boldsymbol{P}}}_{{\rm{A}}} }$。

2) 荷载设置和振动控制效果

考虑风荷载和地震作用两种类型的外部动态激励，验证所提出的POD-LQR方法的优越性。风荷载使用自回归(Autoregressive，AR)模型进行模拟，AR模型阶数设置为5，水平阵风采用Davenport谱，基本风压取0.45 ${{\mathrm{kN}}} /{{\mathrm{m}}^2}$，风荷载持续时间为60 s。将风速转换为作用于上层节点的外力施加在结构上：

式中：${C_{\rm{p}}}$为风压分布系数，为简便取${C_{\rm{p}}} = - 1$；${\rho _{\rm{p}}}$=1.225 ${{\mathrm{kg}}} /{{\mathrm{m}}^3}$为空气密度；${{{\boldsymbol{A}}}_{\rm{p}}}$为作用面积向量；${{{\boldsymbol{V}}}_{\rm{p}}}$为风速向量。

地震作用采用Kobe地震波，时间为20 s，峰值加速度为3.5 g，并将地震波转换为外力施加在结构上：

式中：${{\boldsymbol{M}}}$为质量矩阵；${{{\boldsymbol{I}}}_{g}}$为影响向量；${{\ddot x}_{g}}$为地面运动加速度向量。

风荷载作用下FEM-LQR控制器和POD-LQR控制器的振动控制效果如图4所示。其中，图4(a)~图4(b)为节点1的竖向位移，图4(c)~图4(d)为结构振动能量。图5对比了FEM-LQR控制器和POD-LQR控制器在地震作用下的振动控制效果。从图4和图5可以看出，POD-LQR控制器可以降低结构在风荷载和地震作用下的振动水平，与FEM-LQR控制器的振动控制效果相近。

图  4  风荷载作用下的振动控制效果(Levy索穹顶)

Figure  4.  Vibration control effect under wind load (Levy cable dome)

下载: 全尺寸图片 幻灯片

图  5  地震作用下的振动控制效果(Levy索穹顶)

Figure  5.  Vibration control effect under seismic load (Levy cable dome)

下载: 全尺寸图片 幻灯片

为了定量比较两种控制器的控制效果，定义结构振动能量降低率${\eta _1}$为：

上式表明结构振动能量降低率${\eta _1}$越高，振动控制效果越好。两种控制器的振动能量降低率和控制器设计时间对比如图6所示，可以得到，POD-LQR控制器的控制效果略优于FEM-LQR控制器。具体而言，FEM-LQR控制器和POD-LQR控制器在风荷载作用下的振动能量降低率分别为49%和51%，在地震作用下的结构振动能量降低率分别为51%和53%。将LQR算法求得的矩阵${{{\boldsymbol{S}}}_{\rm{c}}}$代入式(20)左侧，得到Riccati方程的误差指标${p_2}$，定义为：

图  6  振动控制效果与控制器设计时间的比较(Levy索穹顶)

Figure  6.  Comparison of control effect and controller design time (Levy cable dome)

下载: 全尺寸图片 幻灯片

式中，$\Vert ·\Vert$为${l_2}$范数。误差指标${p_2}$越接近零，Riccati方程的解就越精确。在本算例中，FEM-LQR控制器和POD-LQR控制器的误差指标${p_2}$分别为${p_{2,{\rm{F}}{\rm{E}}{\rm{M}} - {\rm{L}}{\rm{Q}}{\rm{R}}}} = 6.83 \times {10^{ - 11}}$和${p_{2,{{\rm{P}}{\mathrm{O}}{\rm{D}}} - {\rm{L}}{\rm{Q}}{{\rm{R}}} }} = 1.47 \times {10^{ - 13}}$。这一结果可归因于POD方法生成的投影子空间来自于结构响应数据，同时降阶状态空间模型的阶数较低，使得求解Riccati方程的精度更高，POD-LQR控制器的控制效果略好于FEM-LQR控制器。此外，FEM-LQR控制器和POD-LQR控制器的设计时间分别为0.60 s和0.02 s，这也表明模型降阶显著提高了POD-LQR控制器的设计效率。

3) 数据源影响分析

值得注意的是，数据源对基于POD获得的降阶模型的质量至关重要。因此，本节研究数据源对模型阶数选择和振动控制效果的影响。基于40组随机荷载和40组脉冲荷载来设计POD-LQR控制器，具体而言，脉冲荷载幅值为1 kN，随机选取3个~10个节点施加脉冲；随机荷载均值为0，方差为100施加到结构的所有节点。数据源对POD-LQR控制器的控制效果的影响如图7所示，可以得到，即使数据源不同，大多数基于位移数据设计的POD-LQR控制器的控制效果都优于FEM-LQR控制器，少量POD-LQR控制器的控制效果略低于FEM-LQR控制器，可能因为结构没有得到充分的激励，得到的降阶模型与原模型误差较大。此外，可以观察到，在相同的能量精度指标${p_1}$要求下，脉冲荷载产生的降阶状态空间模型的阶数小于随机荷载产生的降阶状态空间模型的阶数。

图  7  数据源对模型阶数选择和振动控制效果的影响(Levy索穹顶)

Figure  7.  The influence of data sources on model order selection and vibration control effect (Levy cable dome)

下载: 全尺寸图片 幻灯片

3.2   算例2：30杆球形张拉整体

(1) 模型降阶和LQR控制器设计

30杆球形张拉整体结构(图8)由60个节点和120个单元组成，其中包含30根压杆和90根拉索。压杆和拉索的长度分别为4 m和1.43 m。30杆球形张拉整体的底部5个节点固接。结构参数及预应力设置见表3。基于有限元方法的状态空间模型的阶数为$N = 330$。选择图9中突出显示的10个单元布置作动器。参照上一节的方法选取LQR算法的关键参数取值，具体为${{\boldsymbol{Q}}} = \left[ \begin{matrix} {{{{\boldsymbol{I}}}_{165 \times 165}}}&{{{0}_{165 \times 165}}} \\ {{{0}_{165 \times 165}}}&{{{0}_{165 \times 165}}} \end{matrix} \right]$和${{\boldsymbol{R}}} = 2 \times {10^{ - 12}}{{{\boldsymbol{I}}}_{10 \times 10}}$。

图  8  30杆球形张拉整体

Figure  8.  3D view of 30-strut spherical tensegrity

下载: 全尺寸图片 幻灯片

表  3  30杆球形张拉整体结构参数

Table  3.  Structural parameters settings of 30-strut spherical tensegrity

结构参数	参数取值
总节点数	60
自由节点数	55
单元数	120
节点集中质量/kg	1
拉索弹性模量/GPa	181
压杆弹性模量/GPa	206
拉索等效截面面积/m2	2×10−5
压杆等效截面面积/m2	1×10−4
材料密度/(kg/m3)	7850
阻尼比/(%)	2
安全系数	0.95
拉索预应力/kN	7.45
压杆预应力/kN	−8.00

下载: 导出CSV 

| 显示表格

图  9  30杆球形张拉整体作动器位置示意图(红色突出显示)

Figure  9.  30-strut spherical tensegrity embedded with actuators (highlighted in red).

下载: 全尺寸图片 幻灯片

同样，利用基于有限元法的30杆球形张拉整体模型生成结构加速度响应数据，构建基于POD的降阶状态空间模型。将均值为0，方差为100的随机荷载输入到结构的所有节点上，使结构产生振动，加速度响应数据采集时间为5 s，采样频率为100 Hz。将能量精度指标${p_1}$设为98%，得到基于POD的降阶状态空间模型阶数为${r_{{\rm{A}}} } = 74$。

(2) 荷载设置和振动控制效果

风荷载和地震作用的设置与3.1节一致。FEM-LQR控制器和POD-LQR控制器在风荷载和地震作用下的振动控制效果如图10和图11所示。图12比较了两种控制器的振动能量降低率和控制器设计时间。对比结果表明，POD-LQR控制器的振动控制效果略优于FEM-LQR控制器，产生这一优势的原因是基于POD方法的降阶模型求解Riccati方程的精度较高。具体而言，FEM-LQR控制器和POD-LQR控制器在风荷载作用下的振动能量降低率分别为40%和43%，在地震作用下的结构振动能量降低率分别为53%和56%。FEM-LQR控制器和POD-LQR控制器的误差指标分别为${p_{2,{\rm{F}}{\rm{E}}{\rm{M}} - {\rm{L}}{\rm{Q}}{\rm{R}}}} = 1.72 \times {10^{ - 8}}$和${p_{2,{{\rm{P}}{\mathrm{O}}{\rm{D}}} - {\rm{L}}{\rm{Q}}{{\rm{R}}} }} = 1.59 \times {10^{ - 13}}$。此外，POD-LQR控制器的设计时间为0.05 s，显著低于设计FEM-LQR控制器的2.40 s。

图  10  风荷载作用下的振动控制效果(30杆球形张拉整体)

Figure  10.  Vibration control effect under wind load (30-strut spherical tensegrity)

下载: 全尺寸图片 幻灯片

图  11  地震作用下的振动控制效果(30杆球形张拉整体)

Figure  11.  Vibration control effect under seismic load (30-strut spherical tensegrity)

下载: 全尺寸图片 幻灯片

图  12  振动控制效果与控制器设计时间的比较(30杆球形张拉整体)

Figure  12.  Comparison of control effect and controller design time (30-strut spherical tensegrity).

下载: 全尺寸图片 幻灯片

同样，通过40组随机荷载和40组脉冲荷载对30杆球形张拉整体进行激励，获得设计POD-LQR控制器所需的加速度数据。数据源对POD-LQR控制器的控制效果的影响结果如图13所示。可以看出，大多数基于加速度数据设计的POD-LQR控制器的控制效果超过FEM-LQR控制器，少量POD-LQR控制器的控制效果低于FEM-LQR控制器。此外，对于给定的能量精度指标${p_1}$，可以发现基于脉冲荷载构建的降阶状态空间模型的阶数小于基于随机荷载构建的降阶状态空间模型的阶数。

图  13  数据源对模型阶数选择和振动控制效果的影响(30杆球形张拉整体)

Figure  13.  The influence of data sources on model order selection and vibration control effect (30-strut spherical tensegrity)

下载: 全尺寸图片 幻灯片

4   结论

本文提出了基于本征正交分解的索杆结构模型降阶和POD-LQR振动控制器设计新方法。通过Levy索穹顶和30杆球形张拉整体两个数值算例，将所提出的POD-LQR控制器用于结构在风荷载和地震作用下的振动主动控制。最后，深入探讨了不同数据源对降阶模型阶数选择和振动控制效果的影响。主要结论如下：

(1) 数据驱动的本征正交分解方法从结构位移、速度或加速度响应数据中提取投影子空间。降阶模型的应用提高了Riccati方程的求解精度，其模型阶数可以基于能量精度指标自动选择得到。结果表明，所提出的POD-LQR控制器有效地降低了结构在外部荷载作用下的振动响应，取得了稍好于FEM-LQR控制器的振动控制效果。

(2) 基于本征正交分解方法设计POD-LQR控制器用时更少，Riccati方程的计算复杂度大大降低。对于Levy索穹顶而言，POD-LQR控制器的设计时间为0.02 s，显著低于设计FEM-LQR控制器的0.60 s；对于30杆球形张拉整体而言，POD-LQR控制器的设计时间为0.05 s，显著低于设计FEM-LQR控制器的2.40 s。

(3) POD-LQR控制器的振动控制效果和降阶模型的阶数选择会受数据源影响。大多数基于脉冲荷载和随机荷载的结构响应数据所设计的POD-LQR控制器的振动控制效果都超过FEM-LQR控制器；在相同的能量精度指标要求下，基于脉冲荷载构建的降阶状态空间模型的阶数小于基于随机荷载构建的降阶状态空间模型的阶数。