基于偏微分方程模型降阶方法的最优控制

田容雨; 朱慧

doi:10.13718/j.cnki.xsxb.2019.01.017

基于偏微分方程模型降阶方法的最优控制

田容雨¹,
朱慧²

1.
山东大学(威海) 继续教育学院, 山东威海 264209

2.
山东大学(威海) 教务处, 山东威海 264209

详细信息

作者简介:
田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究 .

中图分类号: TP273

Optimal Control Based on Reduced Order Models of Partial Differential Equations

Rong-yu TIAN¹,
Hui ZHU²

1.
Continuing Education College, Shandong University, Weihai Shandong 264209, China

2.
Office of Teaching Affairs, Shandong University(Weihai Campus), Weihai Shandong 264209, China

摘要: 为了快速准确地求解含有偏微分方程约束（PDE）的优化问题，提出了一种基于偏微分方程模型降阶的最优控制问题求解方法.含有偏微分方程约束会使得优化问题的求解耗费大量的时间，难以满足现有控制与优化的需求.在研究了偏微分方程性质的基础上，得出了一种新的模型降阶方法.通过使用奇异值分解法来提取原模型的主要特性，得到低维空间的基函数，再使用伽辽金投影法，将原模型投影到现有基函数构成的低维空间中，从而达到降低模型阶次来快速计算PDE优化问题的目的.实验结果表明在降阶模型阶次较低的情况下，依然能对原模型有较好的逼近效果.该方法用于快速准确地求解含有偏微分方程约束的优化问题是可行的、有效的.
- 偏微分方程约束 /
- 奇异值分解 /
- 伽辽金方法 /
- 模型降阶
Abstract: In order to solve constrained optimization problems containing partial differential equations (PDE) quickly and accurately, a model has been proposed on the basis of the optimal control problem solved in partial differential order reduction method. Because the issue of partial differential equations containing constraint solving optimization, problems will be time-consuming, which is difficult to meet the needs of existing control and optimization. For such complex issues in the study of the nature of partial differential equations on the basis, we have proposed a new model reduction method. By extracting the main characteristics of the original model in singular value decomposition method to obtain a low-dimensional space of basic functions, and in the Galerkin projection method, low-dimensional space of the original model has been projected onto the existing base function constituted, which can reduce model the purpose of the order to quickly calculate PDE optimization problems. The simulation experiment results show that under the reduced-order model order times lower case, the original model is still able to have a better approximation results. It follows that the method for quickly and accurately solving constrained optimization problems containing partial differential equations is feasible and effective.
- partial differential equation constraints /
- singular value decomposition /
- Galerkin method /
- reduced order model .

偏微分方程被用来描述很多物理现象和自然现象^[1]，在工业过程控制中有很多模型都是基于偏微分方程约束的控制问题和优化问题.但基于偏微分方程的最优控制问题，往往需要大量的时间来求解^[2-3]，无法满足现有优化与控制的需求，因此如何快速求解偏微分方程模型或者对偏微分方程模型进行降阶显得十分必要.

目前常见的偏微分方程的求解方法为有限差分法和有限体积法，都是对偏微分方程进行时间与空间的离散^[4-5]，将偏微分方程化简为一系列代数方程来求解.但是随着系统复杂度的日益提升，离散所得的代数方程的维数往往很高，计算时间很长.文献[6-7]基于有限元的方法来求解偏微分方程最优控制问题，虽然提高了计算精度，依然存在着计算时间长的不足.文献[8-10]提出了基于伽辽金方法来求解偏微分方程的思想，大大降低了该问题的计算量，但是并没有给出一种较好的基函数选取方法，因此求解的准确性相对较差.文献[11-12]提出了基于奇异值分解的方法来降低代数方程的阶次，从而达到简化计算的目的.

本文在结合了伽辽金方法的基础上引入正交基的概念，将原系统空间的能量投影到低维空间中，先采用有限差分法离散偏微分方程，然后根据这些离散的递归高阶代数方程生成快照.在此基础之上，本文提出了偏微分方程的模型降阶方法，能够提取高维状态空间模型的特性来获得一个低维状态空间模型.该方法能够快速准确地求解偏微分方程.本文以一个偏微分方程的最优控制问题为例，进行了仿真实验.实验表明该方法可以快速准确地求解偏微分方程，满足优化问题求解的需求.

1. 一般PDE模型的离散方法说明

1.1 PDE模型的描述

偏微分方程一般用来描述常见的物理现象，模型采用常见的抛物型偏微分方程，同时含有对流项与扩散项，如式(1)

$ \begin{array}{l} \frac{{\partial y}}{{\partial t}} = - v\frac{{\partial y}}{{\partial x}} + D\frac{{{\partial ^2}y}}{{\partial {x^2}}}\\ y\left( {x,0} \right) = {y_0}\;\;\;\;\;\;y \in \left[ {0,L} \right]\\ vy\left( {{0^ - },t} \right) = vy\left( {{0^ + },t} \right) - \beta \frac{{\partial y\left( {x,t} \right)}}{{\partial x}}\left| {_{z = {0^ + }}} \right.\\ \frac{{\partial y\left( {x,t} \right)}}{{\partial x}}\left| {_{z = L}} \right. = 0 \end{array} $

$ (1) $

式中y(x，t)表示温度场，L表示长度，D表示扩散系数.

1.2 PDE模型的离散

一般选择有限差分法作为偏微分方程的离散化方法对式(1)进行离散，可以得到式(2)

$ \frac{{y_i^{N + 1} - y_i^N}}{\tau } = - v\frac{{{y_{i + 1}} - {y_{i - 1}}}}{{2\Delta x}} + D\frac{{{y_{i + 1}} - 2{y_i} + {y_{i - 1}}}}{{\Delta {x^2}}} $

$ (2) $

同理对边界条件离散化之后，带入式(2)可以得到式(3)

$ \begin{array}{*{20}{c}} {y_i^{N + 1} = \left( {\alpha + \beta } \right)y_{i - 1}^N + \left( {1 - 2\beta - {\rm{\Phi }}} \right)y_i^N + \left( {\beta - \alpha } \right)y_{i + 1}^N}&{i \in \left\{ {2,3, \cdots M} \right\}} \end{array} $

$ (3) $

对于边界上的点有所不同，可得式(4)

$ \begin{array}{*{20}{c}} {y_1^{N + 1} = \left[ {1 - 2\beta - \left( {\alpha + \beta } \right)\gamma - {\rm{\Phi }}} \right]y_1^N + 2\beta y_2^N + \left( {\alpha + \beta } \right)m}\\ {y_{M + 1}^{N + 1} = 2\beta y_M^N + \left( {1 - 2\beta - {\rm{\Phi }}} \right)y_{M + 1}^N} \end{array} $

$ (4) $

其中$\alpha = \frac{{v\tau }}{{2\Delta x}}$，$\beta = D\frac{{\Delta t}}{{\Delta {x^2}}}$，$\gamma = \frac{{2v\Delta x}}{D}$，由此易知可以将得到离散之后的递推代数方程式整理成状态空间的形式

$ x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right),y\left( k \right) = cx\left( k \right) $

$ (5) $

易知

$ A = \left| {\begin{array}{*{20}{c}} {{k_1}}&{2\beta }&{}&{}&{}\\ {\alpha + \beta }&{1 - 2\beta }&{\beta - \alpha }&{}&{}\\ {}&{}& \ddots&\ddots&\ddots \\ {}&{}&{\alpha + \beta }&{1 - 2\beta }&{\beta - \alpha }\\ {}&{}&{}&{2\beta }&{1 - 2\beta } \end{array}} \right| $

$ (6) $

其中k₁=1－2β－(α+β)Φ，从上式分析可知，矩阵A与维数和空间离散的步长有关，如果是大规模的问题，A的维数会随之大大增加，此时对该偏微分方程的计算会耗费很多时间.

2. PDE模型降阶处理与误差估计

2.1 PDE模型降阶

对上面所示的状态空间方程进行处理，首先可知

$ y = \left| {\begin{array}{*{20}{c}} {y\left( {{x_1},{t_1}} \right)}&{y\left( {{x_1},{t_2}} \right)}& \cdots &{y\left( {{x_1},{t_k}} \right)}\\ {y\left( {{x_2},{t_1}} \right)}&{y\left( {{x_2},{t_2}} \right)}& \cdots &{y\left( {{x_2},{t_k}} \right)}\\ \vdots&\vdots &{}&{}\\ {y\left( {{x_{M + 1}},{t_1}} \right)}&{y\left( {{x_{M + 1}},{t_2}} \right)}& \cdots &{y\left( {{x_{M + 1}},{t_k}} \right)} \end{array}} \right| $

$ (7) $

从集合y的k个列向量中选取L个列向量构成子集{y}^L，通常情况下L≪M，将这L个列向量构成的集合称为瞬像.在解决实际问题的时候，该瞬像集合通常由对以前的实验或者模拟结果组成，例如在解数值天气预报方程时，可以利用以前的天气结果来构成瞬像集合.

在得到该方程的瞬像之后，利用奇异值分解的方法来研究该偏微分方程的差分格式.首先将瞬像的集合表示成矩阵的形式

$ z = \left| {\begin{array}{*{20}{c}} {z_1^0}&{z_1^1}& \cdots &{z_1^L}\\ {z_2^0}&{z_2^1}& \cdots &{z_2^L}\\ \vdots&\vdots&\ddots&\vdots \\ {z_M^0}&{z_M^1}& \cdots &{z_M^L} \end{array}} \right| $

$ (8) $

对于矩阵z，进行奇异值分解^[13]

$ z = U\left( {\begin{array}{*{20}{c}} S&0\\ 0&0 \end{array}} \right){V^{\rm{T}}} $

$ (9) $

上式中U∈R^M×M，V∈R^L×L都是正交矩阵，S是一个对角矩阵，S=diag(σ₁，σ₂，σ_L+1)，其对角线的元素就是矩阵Z的正奇异值.

若令U=(φ₁，φ₂，…φ_M)，V=(ψ₁，ψ₂，…ψ_L).定义矩阵的范数为‖A‖_α，β=sup$\frac{{{{\left\| {Ax} \right\|}_\alpha }}}{{{{\left\| x \right\|}_\beta }}}$，令Z_M=$\sum\limits_{i = 1}^M {{\sigma _i}{\varphi _i}{\psi _i}} $，如果M < r=rank(Z)，有如下表达式成立

$ \mathop {\min }\limits_{rank\left( B \right) \le M} {\left\| {Z - B} \right\|_{2,2}} = {\left\| {Z - {Z_M}} \right\|_{2,2}} = {\sigma _{M + 1}} $

$ (10) $

这就表明了Z_M是Z的最优逼近.令Z的L+1个列向量为a_M^L=(z₁^l，z₂^l，…z_M^l)，利用矩阵范数和向量范数的相容性有

$ {\left\| {a_M^L - P\left( {a_M^l} \right)} \right\|_2} = \left\| {\left( {Z - {Z_M}} \right)\varepsilon _{L + 1}^l} \right\| \le \left\| {\left( {Z - {Z_M}} \right)} \right\|\left\| {\varepsilon _{L + 1}^l} \right\| = {\sigma _{M + 1}} $

$ (11) $

式中P_M(a_M^l)=$\sum\limits_{i = 1}^M $(φ_i，a_M^l)φ_i(a_M^l，φ_i)是向量φ_i和向量a_M^l的标准内积，在对Z_M构造的过程中发现Φ=(φ₁，φ₂，…φ_M)是一组最优基.

$ y = \sum\limits_{i = 1}^K {{\varphi _i}{\lambda _i}\psi _i^{\rm{T}}} $

$ (12) $

在定义了矩阵范数的情况下，求Φ={φ₁，φ₂，…φ_m}即最小化目标函数

$ J\left( {\rm{\Phi }} \right) = {\left\| {y - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}y} \right\|_2} $

$ (13) $

Φ={φ₁，φ₂，…φ_m}可以理解为标准化正交基的集合

$ \varphi _i^{\rm{T}}{\varphi _j} = \left\{ {\begin{array}{*{20}{c}} 0&{i \ne j}\\ 1&{i = j} \end{array}} \right\} $

$ (14) $

因此可得

$ y\left( k \right) = \varphi a\left( k \right),a\left( k \right) \approx {\varphi ^{\rm{T}}}x\left( k \right) $

$ (15) $

将式(15)带入到式(5)中可得

$ \begin{array}{l} {\varphi ^{\rm{T}}}\varphi a\left( {k + 1} \right) = {\varphi ^{\rm{T}}}A\varphi a\left( k \right) + {\varphi ^{\rm{T}}}Bu\left( k \right)\\ y\left( k \right) = c\varphi a\left( k \right) \end{array} $

$ (16) $

有φ^Tφ=I，由此可得

$ \begin{array}{l} a\left( {k + 1} \right) = {A_r}a\left( k \right) + {B_r}u\left( k \right)\\ y\left( k \right) = {c_r}\varphi a\left( k \right) \end{array} $

$ (17) $

式中A_r=φ^TAφ，B_r=φ^TB，c_r=cφ.降阶模型的解可以由式(18)表示

$ \tilde y\left( {{x_k},t} \right) = \sum\limits_{i = 1}^M {{a_i}\left( t \right)\varphi \left( {{x_k}} \right)} ,k = 1,2, \cdots N + 1 $

$ (18) $

通过使用伽辽金投影方法，可知该模型由M+1阶降到了r阶.采用降阶模型可以大大简化在差分求解中的运算量.在实际的使用过程中如果瞬像的能量达到了真实模型的99.5%以上，就认为该降阶模型可以很好地去逼近原模型.

2.2 解的误差估计

下面对降阶模型的解的误差进行简要分析，设Z^K为有限差分的真实解，Z^*K为降阶模型的解，根据式(16)可得

$ \begin{array}{l} {\left\| {{Z^K} - {Z^{ * K}}} \right\|_2} = {\left\| {\left( {{Z^K} - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}{Z^K}} \right) + \left( {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}{Z^K} - {Z^{ * K}}} \right)} \right\|_2} \le \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\left\| {{Z^K} - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}{Z^K}} \right\|_2} + {\left\| {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}{Z^K} - {Z^{ * K}}} \right\|_2} \end{array} $

$ (19) $

对于上式的前半部分，瞬像集合是由如下空间构成

$ Z_M^k \in span\left\{ {Z_M^0,Z_M^1,Z_M^2, \cdots ,Z_M^L} \right\} $

$ (20) $

于是存在如下的常数r_i，使得式(21)成立

$ Z_M^k = \sum\limits_{i = 0}^L {{r_i}Z_M^i} $

$ (21) $

式(19)的后半部分有

$ Z_M^k - Z_M^{ * k} = A\left( {Z_{M - 1}^k - Z_{M - 1}^{ * k}} \right) $

$ (22) $

由式(19)和式(22)可知

$ (23) $

对上式整理可得

$ \begin{array}{l} \left\| {Z_M^k - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_{M - 1}^k} \right\| + \left\| {Z_M^{ * k} - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_{M - 1}^k} \right\| \le \\ {\left\| {\sum\limits_{i = 0}^L {r_i^k\left( {Z_M^i - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^i} \right)} } \right\|_2} + {\left\| {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}A} \right\|_{2,2}}{\left\| {\left( {Z_M^{k - 1} - Z_M^{ * k - 1}} \right)} \right\|_2} \le \\ \sum\limits_{i = 0}^L {\left| {{r_i}} \right|{{\left\| {Z_M^i - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^i} \right\|}_2}} + {\left\| {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}A} \right\|_{2,2}}{\left\| {\left( {Z_M^{k - 1} - Z_M^{ * k - 1}} \right)} \right\|_2} \end{array} $

$ (24) $

如上式所示，记Z中的L+1个列向量为

$ Z_M^l = \left( {z_1^l,z_1^2, \cdots ,z_M^l} \right) $

$ (25) $

ε_L+1ⁱ是除第L个分量为1，其余为0的单位列向量，同时令$\sum\limits_{i = 0}^L $ |r_i^k|=B_k，‖ΦΦ^TA‖_2，2=C，从上一小节的推导中已经证明了Z_M是Z的最优逼近，因此有

$ \mathop {\min }\limits_{rank\left( {{Z_M}} \right) \le M} {\left\| {Z - {Z_M}} \right\|_{2,2}} = {\sigma _{\left( {M + 1} \right)}} $

$ (26) $

然后根据矩阵范数和向量范数的相容性，有如下表达式成立

$ {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^i = {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z\varepsilon _{L + 1}^i $

$ (27) $

由此可得

$ {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z\varepsilon _{L + 1}^i = {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}\left( {\sum\limits_{i = 1}^M {{\sigma _i}{\varphi _i}\psi _i^{\rm{T}}} } \right)\varepsilon _{L + 1}^i $

$ (28) $

将上式展开可得

$ \begin{array}{l} {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}\left( {\sum\limits_{i = 1}^M {{\sigma _i}{\varphi _i}\psi _i^{\rm{T}}} } \right)\varepsilon _{L + 1}^i = \\ \left[ {\begin{array}{*{20}{c}} {\varphi _1^{\rm{T}}}\\ \vdots \\ {\varphi _M^{\rm{T}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{\sigma _1}}& \cdots &0& \cdots &0\\ \vdots&\ddots&\vdots&\ddots&\vdots \\ 0& \cdots &{{\sigma _M}}& \cdots &0 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\psi _1^{\rm{T}}}\\ \vdots \\ {\psi _M^{\rm{T}}} \end{array}} \right] = \\ \left[ {\begin{array}{*{20}{c}} {\varphi _1^{\rm{T}}}\\ \vdots \\ {\varphi _M^{\rm{T}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{\sigma _1}}& \cdots &0\\ \vdots&\ddots&\vdots \\ 0& \cdots &{{\sigma _M}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\psi _1^{\rm{T}}}\\ \vdots \\ {\psi _M^{\rm{T}}} \end{array}} \right]\varepsilon _{L + 1}^i = \\ \left( {\sum\limits_{i = 1}^M {{\sigma _i}{\varphi _i}\psi _i^{\rm{T}}} } \right)\varepsilon _{L + 1}^i = \\ {Z_M}\varepsilon _{L + 1}^i \end{array} $

$ (29) $

由式(24)和式(29)可知

$ \begin{array}{l} \sum\limits_{i = 0}^L {\left| {{r_i}} \right|{{\left\| {Z_M^i - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^i} \right\|}_2}} + {\left\| {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}A} \right\|_{2,2}}{\left\| {\left( {Z_M^{k - 1} - Z_M^{ * k - 1}} \right)} \right\|_2} \le \\ {B_k}{\left\| {Z_M^i - {Z_M}\varepsilon _{L + 1}^i} \right\|_2} + C{\left\| {\left( {Z_M^{k - 1} - Z_M^{ * k - 1}} \right)} \right\|_2} \le \\ {B_k}{\left\| {Z_M^i - {Z_M}\varepsilon _{L + 1}^i} \right\|_2} + C{\left\| {\left( {Z_M^{k - 1} - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^{k - 1}} \right) + \left( {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^{k - 1} - Z_M^{ * k - 1}} \right)} \right\|_2} \le \\ {B_k}{\left\| {Z - {Z_M}} \right\|_{2,2}}{\left\| {\varepsilon _{L + 1}^i} \right\|_2} + C\left\| {Z_M^{k - 1} - {\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^{k - 1}} \right\| + C\left\| {{\rm{\Phi }}{{\rm{\Phi }}^{\rm{T}}}Z_M^{k - 1} - Z_M^{ * k - 1}} \right\| \le \\ \left( {{B_k} + C{B_{k - 1}} + \cdots + {C^{k - 1}}B} \right){\sigma _{M + 1}} + {C^k}{\left\| {{\rm{\Phi }}\beta _M^0 - Z_M^{ * 0}} \right\|_2} = \\ {{\bar D}_{{\sigma _{M + 1}}}} \end{array} $

综上所述，可以得到

$ {\left\| {{Z^k} - {Z^{ * k}}} \right\|_2} \le {{\bar D}_{{\sigma _{M + 1}}}} $

$ (30) $

由此可以证明，采用模型降阶方法得到的近似解和真实解之间的距离是有界的.在适当选取其中参数的情况下，可以保证近似解能够较好地逼近真实解.

3. 仿真实验

考虑式(31)所示的一个对流扩散的优化问题

$ \begin{array}{l} \min J = \int_0^\infty {{{\left[ {y\left( t \right) - {y_d}} \right]}^2}{\rm{d}}t} + {\alpha ^2}\int_0^\infty {{{\left[ {u\left( t \right) - {u_d}} \right]}^2}{\rm{d}}t} \\ \frac{{\partial y}}{{\partial t}} = - v\frac{{\partial y}}{{\partial x}} + D\frac{{{\partial ^2}y}}{{\partial {x^2}}} + R\left( {x,t} \right)\\ y\left( {x,0} \right) = {y_0}\;\;\;\;\;\;\;\;\;y \in \left[ {0,L} \right]\\ vy\left( {{0^ - },t} \right) = vy\left( {{0^ + },t} \right) - \beta \frac{{\partial y\left( {x,t} \right)}}{{\partial x}}\left| {_{z = {0^ + }}} \right.\\ \frac{{\partial y\left( {x,t} \right)}}{{\partial x}}\left| {_{z = L}} \right. = 0 \end{array} $

$ (31) $

首先采用有限差分法对该优化问题进行离散化处理，由此可以得到式(32)所示的一个最优控制问题.

$ \begin{array}{l} \min J = \sum\limits_{n = 0}^\infty {\left[ {{x^{\rm{T}}}Qx + {u^{\rm{T}}}Ru} \right]} \\ x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\ y\left( k \right) = Cx\left( k \right) \end{array} $

$ (32) $

通过线性二次调节器的反馈控制设计可知

$ u = - Kx,K = {\left( {{B^{\rm{T}}}SB + R} \right)^{ - 1}}\left( {{B^{\rm{T}}}SA} \right) $

$ (33) $

S是通过求解黎卡提方程得到的

$ {A^{\rm{T}}}SA - S - \left( {{A^{\rm{T}}}SB} \right){\left( {{B^{\rm{T}}}SB + R} \right)^{ - 1}}\left( {{B^{\rm{T}}}SA} \right) + Q = 0 $

$ (34) $

在求解该问题的过程中，分别采用不同阶次的简化模型和真实模型进行对比.

根据经验值得到表 1中实验参数，并按照表 1中参数进行实验. 图 1给出了采用本文方法中近似模型不同阶次r与总能量的关系.

图 1 近似模型阶次与总能量的关系

下载: 全尺寸图片幻灯片

表 1 实验参数

参数	数值
v	0.5
L	1
D	0.01
N	200
时间步长	0.001
空间步长	0.005

| Show Table

DownLoad: CSV

由此可以得到当近似模型的阶次达到10以上时，总能量已经解决100%了.由一般的工程经验可知，若瞬像的总能量达到了原系统总能量的99.5%以上，那么就认为有较好的近似效果.

图 2、图 3和图 4分别是近似模型阶次依次选取10，15，20时，真解与近似解的误差曲线.从分析可以看出，当r=10时真解与近似解的误差最大，r=20时真解与近似解的误差最小.但总体而言，近似模型阶次若大于10，则可以较好地逼近原模型.

图 2 r=10时真解与近似解的误差

下载: 全尺寸图片幻灯片

图 3 r=15真解与近似解的误差

下载: 全尺寸图片幻灯片

图 4 r=20真解与近似解的误差

下载: 全尺寸图片幻灯片

图 5是近似模型阶次依次选取10，15，20时，优化求解的控制量的大小变化曲线，从图 5中分析可知，控制量大小几乎是重合的. 图 6是近似模型阶次依次选取10，15，20时，优化求解的控制量的大小变化曲线和控制量一样，状态量在不同阶次模型中几乎重合.

图 5 不同近似模型阶次控制量变化曲线

下载: 全尺寸图片幻灯片

图 6 不同近似模型阶次状态量变化曲线

下载: 全尺寸图片幻灯片

如图 7所示，分别使用近似模型和真实模型的解析解对比，从图 7中分析可知，使用近似模型求解的结果和真实模型解析解的结果几乎一致，由此说明了该方法的有效性.

图 7 解析解与近似模型比较

下载: 全尺寸图片幻灯片

4. 结论

本文在分析了几种偏微分方程求解方法的基础上，基于伽辽金投影法，通过合理地选取正交基函数来近似高阶代数方程，从而达到简化计算的目的.在理论上证明了该方法的近似解可以收敛于真实解，并且以一个最优控制问题为例进行了仿真实验，证明了该方法的有效性.本文接下来的研究方向是考虑将优化的一些思想运用到降阶模型中，建立一个更加准确的降阶模型.

图 1 近似模型阶次与总能量的关系

下载: 全尺寸图片幻灯片

图 2 r=10时真解与近似解的误差

下载: 全尺寸图片幻灯片

图 3 r=15真解与近似解的误差

下载: 全尺寸图片幻灯片

图 4 r=20真解与近似解的误差

下载: 全尺寸图片幻灯片

图 5 不同近似模型阶次控制量变化曲线

下载: 全尺寸图片幻灯片

图 6 不同近似模型阶次状态量变化曲线

下载: 全尺寸图片幻灯片

图 7 解析解与近似模型比较

下载: 全尺寸图片幻灯片

表 1 实验参数

参数	数值
v	0.5
L	1
D	0.01
N	200
时间步长	0.001
空间步长	0.005

下载: 导出CSV

[1]	doi: http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=bb7f19b5539e4b38b5a1012b800b968b ZIENIUK E, SAWICKI D, BOLTUC A. Parametric Integral Equations Systems in 2D Transient Heat Conduction Analysis[J].International Journal of Heat and Mass Transfer, 2014, 78(12):571-587.
[2]	周宏宇, 王爱民.非线性偏微分方程数值求解的自适应方法研究[J].计算机应用与软件, 2012, 47(20):38-40, 84. doi: http://d.old.wanfangdata.com.cn/Periodical/jsjgcyyy201120011
[3]	杨健, 赖晓霞.辅助函数法求解非线性偏微分方程精确解[J].计算机技术与发展, 2017, 27(10):196-200. doi: http://d.old.wanfangdata.com.cn/Periodical/wjfz201711042
[4]	doi: http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=gcsxxb201002022 FU H, RUI H.Characteristic-Mixed Methods on Dynamic Finite Element Spaces for Two-Phase Miscible flow in Porous Media[J].Chinese Journal of Engineering Mathematics, 2013, 27(2):369-374.
[5]	李志华, 喻军, 杨红光.偏微分方程与微分代数方程的一致求解方法[J].中国机械工程, 2015, 26(4):441-445. doi: 10.3969/j.issn.1004-132X.2015.04.003
[6]	doi: http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=fd05976f49add4be0afa944f2c0fe52b FENG X, NEILAN M.Finite Element Approximations of General Fully Nonlinear Second Order Elliptic Partial Differential Equations Based on the Vanishing Moment Method[J].Computers & Mathematics with Applications, 2014, 68(12):2182-2204.
[7]	李明, 崔向照, 赵金娥.求解高次有限元方程的外推瀑布型多重网格法[J].西南师范大学学报(自然科学版), 2016, 41(1):20-23. doi: http://d.old.wanfangdata.com.cn/Periodical/xnsfdxxb201601003
[8]	JIN B, LAZAROV R, LIU Y K, et al.The Galerkin Finite Element Method for a Multi-Term Time-Fractional Diffusion Equation[J].Journal of Computational Physics, 2015, 281:825-843. doi: 10.1016/j.jcp.2014.10.051
[9]	TAN B T, GHATTAS O.A PDE-Constrained Optimization Approach to the Discontinuous Petrov-Galerkin Method with a Trust Region Inexact Newton-CG Solver[J].Computer Methods in Applied Mechanics & Engineering, 2014, 278(6):20-40.
[10]	李洋洋.基于自由度缩减的局部无网格伽辽金方法研究[D].长沙: 湖南大学, 2013.http://cdmd.cnki.com.cn/article/cdmd-10532-1014119579.htm
[11]	陈军胜.基于四元数矩阵分解的线性方程组解的存在性判断研究[J].西南师范大学学报(自然科学版), 2015, 40(5):34-38. doi: http://d.old.wanfangdata.com.cn/Periodical/xnsfdxxb201505008
[12]	谭伟杰, 冯西安, 张杨梅.基于Hankel矩阵分解的互素阵列高分辨目标定向[J].西南大学学报(自然科学版), 2016, 38(7):191-198. doi: http://d.old.wanfangdata.com.cn/Periodical/xnnydxxb201607030
[13]	李葵, 范玉刚, 吴建德.基于二次奇异值分解和最小二乘支持向量机的轴承故障诊断方法[J].计算机应用, 2014, 34(8):2438-2441. doi: http://d.old.wanfangdata.com.cn/Periodical/jsjyy201408061

期刊类型引用(0)

其他类型引用(3)

图( 7) 表( 1)

计量

文章访问数: 1404
HTML全文浏览数: 1113
PDF下载数: 112
施引文献: 3

全文HTML

偏微分方程被用来描述很多物理现象和自然现象^[1]，在工业过程控制中有很多模型都是基于偏微分方程约束的控制问题和优化问题.但基于偏微分方程的最优控制问题，往往需要大量的时间来求解^[2-3]，无法满足现有优化与控制的需求，因此如何快速求解偏微分方程模型或者对偏微分方程模型进行降阶显得十分必要.

目前常见的偏微分方程的求解方法为有限差分法和有限体积法，都是对偏微分方程进行时间与空间的离散^[4-5]，将偏微分方程化简为一系列代数方程来求解.但是随着系统复杂度的日益提升，离散所得的代数方程的维数往往很高，计算时间很长.文献[6-7]基于有限元的方法来求解偏微分方程最优控制问题，虽然提高了计算精度，依然存在着计算时间长的不足.文献[8-10]提出了基于伽辽金方法来求解偏微分方程的思想，大大降低了该问题的计算量，但是并没有给出一种较好的基函数选取方法，因此求解的准确性相对较差.文献[11-12]提出了基于奇异值分解的方法来降低代数方程的阶次，从而达到简化计算的目的.

本文在结合了伽辽金方法的基础上引入正交基的概念，将原系统空间的能量投影到低维空间中，先采用有限差分法离散偏微分方程，然后根据这些离散的递归高阶代数方程生成快照.在此基础之上，本文提出了偏微分方程的模型降阶方法，能够提取高维状态空间模型的特性来获得一个低维状态空间模型.该方法能够快速准确地求解偏微分方程.本文以一个偏微分方程的最优控制问题为例，进行了仿真实验.实验表明该方法可以快速准确地求解偏微分方程，满足优化问题求解的需求.

1. 一般PDE模型的离散方法说明

1.1. PDE模型的描述

偏微分方程一般用来描述常见的物理现象，模型采用常见的抛物型偏微分方程，同时含有对流项与扩散项，如式(1)

式中y(x，t)表示温度场，L表示长度，D表示扩散系数.

1.2. PDE模型的离散

一般选择有限差分法作为偏微分方程的离散化方法对式(1)进行离散，可以得到式(2)

同理对边界条件离散化之后，带入式(2)可以得到式(3)

对于边界上的点有所不同，可得式(4)

其中$\alpha = \frac{{v\tau }}{{2\Delta x}}$，$\beta = D\frac{{\Delta t}}{{\Delta {x^2}}}$，$\gamma = \frac{{2v\Delta x}}{D}$，由此易知可以将得到离散之后的递推代数方程式整理成状态空间的形式

易知

其中k₁=1－2β－(α+β)Φ，从上式分析可知，矩阵A与维数和空间离散的步长有关，如果是大规模的问题，A的维数会随之大大增加，此时对该偏微分方程的计算会耗费很多时间.

2. PDE模型降阶处理与误差估计

2.1. PDE模型降阶

对上面所示的状态空间方程进行处理，首先可知

从集合y的k个列向量中选取L个列向量构成子集{y}^L，通常情况下L≪M，将这L个列向量构成的集合称为瞬像.在解决实际问题的时候，该瞬像集合通常由对以前的实验或者模拟结果组成，例如在解数值天气预报方程时，可以利用以前的天气结果来构成瞬像集合.

在得到该方程的瞬像之后，利用奇异值分解的方法来研究该偏微分方程的差分格式.首先将瞬像的集合表示成矩阵的形式

对于矩阵z，进行奇异值分解^[13]

上式中U∈R^M×M，V∈R^L×L都是正交矩阵，S是一个对角矩阵，S=diag(σ₁，σ₂，σ_L+1)，其对角线的元素就是矩阵Z的正奇异值.

若令U=(φ₁，φ₂，…φ_M)，V=(ψ₁，ψ₂，…ψ_L).定义矩阵的范数为‖A‖_α，β=sup$\frac{{{{\left\| {Ax} \right\|}_\alpha }}}{{{{\left\| x \right\|}_\beta }}}$，令Z_M=$\sum\limits_{i = 1}^M {{\sigma _i}{\varphi _i}{\psi _i}} $，如果M < r=rank(Z)，有如下表达式成立

这就表明了Z_M是Z的最优逼近.令Z的L+1个列向量为a_M^L=(z₁^l，z₂^l，…z_M^l)，利用矩阵范数和向量范数的相容性有

式中P_M(a_M^l)=$\sum\limits_{i = 1}^M $(φ_i，a_M^l)φ_i(a_M^l，φ_i)是向量φ_i和向量a_M^l的标准内积，在对Z_M构造的过程中发现Φ=(φ₁，φ₂，…φ_M)是一组最优基.

在定义了矩阵范数的情况下，求Φ={φ₁，φ₂，…φ_m}即最小化目标函数

Φ={φ₁，φ₂，…φ_m}可以理解为标准化正交基的集合

因此可得

将式(15)带入到式(5)中可得

有φ^Tφ=I，由此可得

式中A_r=φ^TAφ，B_r=φ^TB，c_r=cφ.降阶模型的解可以由式(18)表示

通过使用伽辽金投影方法，可知该模型由M+1阶降到了r阶.采用降阶模型可以大大简化在差分求解中的运算量.在实际的使用过程中如果瞬像的能量达到了真实模型的99.5%以上，就认为该降阶模型可以很好地去逼近原模型.

2.2. 解的误差估计

下面对降阶模型的解的误差进行简要分析，设Z^K为有限差分的真实解，Z^*K为降阶模型的解，根据式(16)可得

对于上式的前半部分，瞬像集合是由如下空间构成

于是存在如下的常数r_i，使得式(21)成立

式(19)的后半部分有

由式(19)和式(22)可知

对上式整理可得

如上式所示，记Z中的L+1个列向量为

ε_L+1ⁱ是除第L个分量为1，其余为0的单位列向量，同时令$\sum\limits_{i = 0}^L $ |r_i^k|=B_k，‖ΦΦ^TA‖_2，2=C，从上一小节的推导中已经证明了Z_M是Z的最优逼近，因此有

然后根据矩阵范数和向量范数的相容性，有如下表达式成立

由此可得

将上式展开可得

由式(24)和式(29)可知

综上所述，可以得到

由此可以证明，采用模型降阶方法得到的近似解和真实解之间的距离是有界的.在适当选取其中参数的情况下，可以保证近似解能够较好地逼近真实解.

4. 结论

本文在分析了几种偏微分方程求解方法的基础上，基于伽辽金投影法，通过合理地选取正交基函数来近似高阶代数方程，从而达到简化计算的目的.在理论上证明了该方法的近似解可以收敛于真实解，并且以一个最优控制问题为例进行了仿真实验，证明了该方法的有效性.本文接下来的研究方向是考虑将优化的一些思想运用到降阶模型中，建立一个更加准确的降阶模型.

参考文献 (13)

姓名

邮箱

手机号码

标题

留言内容

验证码
	验证码错误！

留言板

基于偏微分方程模型降阶方法的最优控制

1.
山东大学(威海) 继续教育学院, 山东威海 264209

2.
山东大学(威海) 教务处, 山东威海 264209

作者简介:
田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究 .

Optimal Control Based on Reduced Order Models of Partial Differential Equations

1.
Continuing Education College, Shandong University, Weihai Shandong 264209, China

2.
Office of Teaching Affairs, Shandong University(Weihai Campus), Weihai Shandong 264209, China

1. 一般PDE模型的离散方法说明

1.1 PDE模型的描述

1.2 PDE模型的离散

2. PDE模型降阶处理与误差估计

2.1 PDE模型降阶

2.2 解的误差估计

3. 仿真实验

4. 结论

期刊类型引用(0)

其他类型引用(3)

计量

基于偏微分方程模型降阶方法的最优控制

作者简介: 田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究
1. 山东大学(威海) 继续教育学院, 山东威海 264209

2. 山东大学(威海) 教务处, 山东威海 264209

English Abstract

Optimal Control Based on Reduced Order Models of Partial Differential Equations

全文HTML

1.1. PDE模型的描述

1.2. PDE模型的离散

2.1. PDE模型降阶

2.2. 解的误差估计

目录

全文HTML

1.1. PDE模型的描述

1.2. PDE模型的离散

2.1. PDE模型降阶

2.2. 解的误差估计

留言板

基于偏微分方程模型降阶方法的最优控制

1. 山东大学(威海) 继续教育学院, 山东 威海 264209 2. 山东大学(威海) 教务处, 山东 威海 264209

作者简介: 田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究 .

Optimal Control Based on Reduced Order Models of Partial Differential Equations

1. Continuing Education College, Shandong University, Weihai Shandong 264209, China 2. Office of Teaching Affairs, Shandong University(Weihai Campus), Weihai Shandong 264209, China

1. 一般PDE模型的离散方法说明

1.1 PDE模型的描述

1.2 PDE模型的离散

2. PDE模型降阶处理与误差估计

2.1 PDE模型降阶

2.2 解的误差估计

3. 仿真实验

4. 结论

期刊类型引用(0)

其他类型引用(3)

计量

出版历程

基于偏微分方程模型降阶方法的最优控制

作者简介: 田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究 1. 山东大学(威海) 继续教育学院, 山东 威海 264209 2. 山东大学(威海) 教务处, 山东 威海 264209

English Abstract

Optimal Control Based on Reduced Order Models of Partial Differential Equations

全文HTML

1.1. PDE模型的描述

1.2. PDE模型的离散

2.1. PDE模型降阶

2.2. 解的误差估计

目录

全文HTML

1.1. PDE模型的描述

1.2. PDE模型的离散

2.1. PDE模型降阶

2.2. 解的误差估计

1.
山东大学(威海) 继续教育学院, 山东威海 264209

2.
山东大学(威海) 教务处, 山东威海 264209

作者简介:
田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究 .

1.
Continuing Education College, Shandong University, Weihai Shandong 264209, China

2.
Office of Teaching Affairs, Shandong University(Weihai Campus), Weihai Shandong 264209, China

作者简介: 田容雨(1980-), 男, 硕士, 实验师, 主要从事数据挖掘研究
1. 山东大学(威海) 继续教育学院, 山东威海 264209

2. 山东大学(威海) 教务处, 山东威海 264209