什么是稳定性0
Lyapunov stability:
首先要定义一个平衡态,之后,系统会围绕这个状态受到扰动。如果系统返回到平衡状态,则认为稳定的;如果系统偏离平衡状态系统(或者更准确地说,这个特定的平衡状态)被认为是不稳定的。在Lyapunov稳定性的定义中,一个无限的时间范围允许返回到平衡。
短时间的Lyapunov stability:(放大因子或增益)
E(q) denotes the kinetic energy of the perturbation q.
singular value decomposition
任意函数可以展开成vj基函数的形式
The maximum amplitude of output for a given input amplitude?
M的作用如图所示,对圆盘进行拉伸和旋转
SVD截断
A least squares solution x0 that minimises $|M x-b|$.
全局稳定性分析方法
我们来考虑N-S方程的统一动力学形式:
N 是 Navier-Stokes operator.
定义base flow
The temporal evolution of perturbations q’
the linearized Navier-Stokes operator about the base flow, also called the Jacobian of the system
A global linear stability analysis consists of searching modes of the form
reduces N-S equation into the eigenproblem:
The corresponding eigenvalues characterize the stability of $q_b$: if an eigenvalue $\lambda_i$ has a positive real part $\delta_i = real(\lambda_i)$ (positive growth rate), then the base flow is unstable. Imaganry part means the oscillating mode at a frequency $w_i$.
Resolvent analysis for the of pseudo-resonance
另一种全局方法,被称为预解分析,包括研究外部作用力对baseflow的影响流动而不是专注于自我维持的振荡。类似在机械振动领域的谐波激励下的稳态响应。这个提法最初Trefethen引入,1993.
Resolvent analysis最早应用是在稳定性和转捩中的nonmetal analysis中运用。例如在寻找线性方程组的矩阵时间指数的范数,对transient growth的研究,e-pseudospectra, 以及定义线性NS在对外部forcing最大反馈的比例。
基本形式可以写成:
代入上式可得:
R 被称为resolvent operator, 或者叫放大因子。
The optimal forcing which maximizes the energy gain G
Mean flow stability vs base flow stability analysis
In some cases (e.g. the cylinder flow), the frequency selection process is strongly driven by nonlinear effects. 因此,线性基流稳定性分析不能正确地预测流体的动力学行为,而这产生不满足RPIF(((real part positive) with a frequency (the imaginary part))属性的情况.
最终的结论是mean flow 比base flow 更有效。
整体程序框架
考虑因素:stability analysis、resolvent analysis、Arbitrary Lagrangian Eulerian methods 可考虑流固耦合的统一框架、implicit newton approach for non-linearity、选用成熟的FVM或者FEM开源C++求解器、小扰动假设,主要关心的问题在于流体或者固体起振之初,比如失速起始、无量纲化
固体域$\Omega_s$:
其中,$M_s$是质量矩阵,P为Kirchhhoff 应力张量.
另外,对于可压固体,我们还需要对其密度进行修正:
流体域:
暂时考虑不可压的N-S方程:
$\sigma$ 为Cauchy 应力张量,跟流体的形变和雷诺数又关系。
- 固体域使用Lagrangian框架、流体域使用Eulerian框架,这样的好处是表达更加自然美观。
Arbitrary Lagrangian Eulerian methods
在于考虑相互的交界面:
因此,得到extension displacement in the reference fluid domain:
链式法则,
最终,ALE formulation of N-S equation:
对比之前的N-S,我们可以发现
弱解形式,应用FEM求解
freefem++ 参考代码:https://github.com/ERC-AEROFLEX/ALE
2.1.2公式转化为稳定性问题:
简化矩阵:
问题? 如何衔接到n-periodic的模型上面去?
Linear Stability Analysis 更深一层次的理解
这是一个Duffing oscillator动力学系统,用相空间的形式表示。
如图,可以看到三个fixed point, 或者称为equilibrium points $\mathbf{X}^{*}$。
找到equilibrium points之后,我们需要来判断这个状态是stable还是unstable的。 Lyapunov stability的解释是,过了很长时间能回到equilibrium points。
首先,假设系统的perturbation满足如下动力学方程:
根据Taylor 公式在平衡点展开,
$\mathcal{A}$ 为Jacobian matrix of F.
也就是说,假设给定初始条件x_0, x可以解析得到:
我们对A进行spectral decomposition:
得到x的另外一种表达形式:
可以发现,其特征值如果按照实部从大到小排列,
也就是说,$\lambda_1$就决定了系统的稳定性,如果实部大于0,那么系统就会呈指数增长,即线性不稳定;如果实数小于0,那么系统就会指数衰减,即线性稳定。实数等于0的情况较为特殊,无法直接判断,需要再继续判断A的二阶Taylor展开。
这里主要针对长时间的稳定性问题,短时间内需要后面再进行理论的探讨。
Non-modal Stability Analysis 更深一层次的理解
给一个例子,
其特征值为
根据特征值,可以判断长时间内的衰减趋势。但是短时间内,两者都大幅度的指数放大。这是在之前Stability Analysis无法解释的。
这是由于矩阵A的non-normality所带来的:
$\mathcal{A}^{\dagger} $为伴随矩阵,它具有非常好的性质。关于伴随矩阵和逆矩阵的关系,和其性质,参考链接。
As a result of this non-normality, the eigenvectors of A do not form an orthonormal set of vectors。
长时间的过程,都会趋向于线性稳定。但是,non-modal在短时间内可能会造成瞬时复制增大。本质是应为non-modal的矩阵性质,分解得到非正交的特征,物理上表现在更初始值的选取有关。但一旦逼近平衡点,就会按照平衡点的性质稳定下来。
从图中解释,目前这两个特征值是非正交的,初始的选取位置不一样,造成短期的增长效应。但如果是正交的,那么类似于垂直的网格布置,两个特征值的互不影响,这样也就是说,直接取消于稳定点的趋势。
Optimal Perturbation Analysis 紧接着引出
也就是,即使是stable的点逼近,但是短时间内也可能产生大的放大响应,找到最大的放大点是一个非常有意思的问题。
两个策略:
求优化问题
并转化为:
SVD
引入SVD分解:
那么,最大的特征值即为最大能量响应:
其对应的初始x_0 由对应的右奇异特征值给出。
Resolvent Analysis
the response of the system to the forcing f(t) 表现为卷积的形式:
Consider linear stable systems,influence of the forcing on the current state decays exponentially according to the least stable eigenvalue。
根据卷积和傅立叶变换的性质,可以很容易的得到响应矩阵,这个矩阵也叫做Resolvent operator:
类似,寻找最大的响应放大因子:
同样,可以利用SVD寻找到最大的奇异值:
线化不可压N-S方程
将均匀项和脉动量分开
fluaction equation:
RANS equation: (w=0)
消除压力项:
pressure Poisson equation
分解得到:
mean component
frequency component
代入到前面的fluaction equation,消去压力项:
前面的推导的具体过程参考5
最后总结框图如下:
首先,我们将速度和压力分解,
预解算子:
对管道坐标系来说:
见代码:https://github.com/mluhar/resolvent
McKeon & Sharma (2010), 做的是在管道里面。
预解算子的边界条件
硬壁面无滑移边界条件
对于软壁面边界条件的思路🤔:
Complex wall admittance linking pressure and velocity。2
Impose boundary conditions for Navier-Stokes Resolvent This accounts for a simple compliant surface model, which essentially leads to a harmonic-motion relationship between pressure and velocity。
探究谐波响应激起的最大压力响应,而非速度 1
possion equation:
boundary condition at the wall:
构建格林公式满足:
可以解析得到压力对格林公式的表达式:
加强筋边界-channel
Solid obstructions within the fluid domain—riblets in our case—are modeled as a spatially varying permeability function Kx; y; z,
which appears as an additional body force in the momentum equations. 3
【方程1】
在流体域,K趋向于∞;在固体域,K趋向于0;
【方程2】
同样的,meanflow也会呈现类似的周期性特征:
【方程3】
在此之上,脉动小量被提取出来:可以看到,均匀流与x的方向无关
【方程4】
代入动量方程,均匀量满足关系:
【方程5】
脉动量满足:
【方程6】
将总结的均匀流的周期性特性【方程2、3】代入到【方程5】,得:
【方程7】
脉动量方程6的分解形式:
【方程8】
方程8可以整理成矩阵的形式:
【方程9】
可以发现,k和k_s空间波数耦合。k_s在压气机中类似叶片扫过的均匀流。nk_s 表示其多阶谐波。
将方程9整理成输入和输出的形式:
最终,写成简洁的形式:
SVD分解
对预解算子重新整理:
f_k仅仅于速度有关系,因此也可以写成:
svd分解:
对应的force mode:
对应的response mode:
基函数满足正交性质:
响应的放大程度表现在奇异值的大小。
整体目标:把握信号的特点和机理的关系,构建模型
验证模型1:
2d,周期性变换,均匀流,环形
优点:定向分析,目前论文分析流肋条大小对减阻对效果
缺点:未考虑湍流的影响
验证模型2:
循环矩阵,应用到上述的模型1中,DNS(湍流模型)
优点:模型1的局部可以得到涡流的细节,利用peter schmid 循环矩阵的思路,再进行周期化,并做预解分析
验证模型3:
https://su2code.github.io/tutorials/Unsteady_NACA0012/
https://su2code.github.io/tutorials/Inc_Turbulent_NACA0012/
https://su2code.github.io/tutorials/Unsteady_Shape_Opt_NACA0012/
可压的resolvent anlaysis N-S model
参考教材6
可压LES model (Mach > 0.2,Mach < 6)
基本方程推导
推导
能量方程可以写成焓的表达形式:
viscous dissipation
能量方程可以写成internal energy的表达形式:
能量方程可以写成pressure的表达形式:
能量方程可以写成entropy的表达形式:
LES的滤波概念:尺度分离
假设湍流各向同性、均匀性。滤波器特性与坐标系无关。非各向同性可以利用非中心滤波器。滤波本质是卷积的思想。无非类似空间的格林函数。
其中,cut-off scale in space and time,$\Delta$ ,$\tau_{c}$.
该方程可写成时域的简化形式:
转化为频域,
或者:
non-resolved part:
常用的三种滤波器类型:
利用taylor公式:
代入公式后,得
我们的目标当然是整理成矩阵的离散形式:
Commutation Error
这里是建模的核心。
LES 方程
利用滤波器Favre Filtering
分别分解为低频项和高频项
主控方程组:
【方程1】
【方程1滤波形式】:方程整体滤波,卷积核函数操作,不可避免引入Commutation Error,通常密度最稳定,分离误差较小。
【方程2-焓的形式】
viscous dissipation
【方程2-焓的滤波形式】:
最终写成(此处漏掉了一些推导过程):
【方程2-温度的形式】
【方程2-温度的滤波形式】
简化成:
【方程2-压力的形式】
【方程2-压力的滤波形式】
简化:
【方程2-熵的形式】
【方程2-熵的滤波形式】
【方程3-状态方程】
【方程3-状态方程的滤波形式】
第二项没有办法直接计算,因此,需要引入模型。
【模型】:还有一些变参数模型,不具体介绍
等价:
和
【方程4-Total Energy】
【方程4-Total Energy滤波形式】
简化:
还有一些类似的模型,具体见📖。
【方程5-动量方程】
【方程5-动量方程滤波形式】
这里,模型还略去了压力的脉动量,也就是说,声学小扰动项不在les作为重点考虑。
Resolvent analysis of compressible N-S equation
可压N-S方程:
strain and stress tensor:
steady version of compressible Navier- Stokes equations:
脉动量:
上述方程是基于baseflow的。一些论文结果表明,在线性稳定性分析过程中,其结果不理想。
接下来,用真正的meanflow再分析一波:周期性满足下,对t的导数为0;
代入
可以写成:
Weak form of incompressible N-S equation
其弱解形式:
整个过程需要考虑边界条件:
最终,整理出易于编写代码的形式:(引入分步积分)
通过Newton Iteration求解base-flow
首先猜测初始值,然后不断的求梯度逼近,
代入上面的式子,类似taylor分解操作:
同样整理成弱解形式:
其中,
最后整理出如下形式:
Linear Stability
将参数分解为baseflow和脉动量的表达形式,导入N-S方程
其弱解形式也可以写出:
最终整理出离散的矩阵形式:
Adjoint Eigenvalue Problem
Connection between adjoint and resolvent
Adjoint Eigenvalue Problem (resolvent anlaysis的再次不同理解)
前面已经提到线性N-S方程其弱解形式也可以写出:
Adjoint 形式:
代入脉动量,最后得到频域表达式:
伴随特征值和特征值互为复共轭关系。
若
则
Iterative Methods for Eigenvalue Computations
Compressible N-S equation
baseflow solution:
eigenvalue:
Meanflow
在非线性域,baseflow引入的线性小扰动假设不再成立。因此,有必要重新定义新的”baseflow“,称为:meanflow
Ajoint 9
q plus叫做ajoint或者Langrange multiplier.
优化问题 转化为求极值问题:
这三项分别可以写成:
Ajoint优化问题的求解思路(梯度下降):
第1个方程: u_guess, 然后求解;将q代入到第2方程分别求解q_plus;然后再代到第3个方程里面,即可得到gradient,或者叫做sentisitivity。作为优化问题,我们明确了迭代的方向,通过优化u,再重新循环迭代。直到,梯度为0。表示达到了极值条件。
在深度学习思路中,同样有类似的cost fucntion, 最终殊途同归,利用梯度下降法,求极值问题。
根据此思路,我们可以了解系统的stability, receptivity, sensitivity.
The role of ajoint
- 有自己独立的系统方程
how does a perturbution of the governing equations influence the output measure? (引入小扰动)
外部扰动(加扰动)
内部扰动(比如,将雷诺数上调一些,或其他参数或形状稍微变化一些)
将其代入到上面的变分方程中,得到的类似的线性小扰动方程(线性化过程):
外部扰动
内部扰动
总结:表明包含 ajoint系统的sensitivity信息,可以用来判断:流场哪些区域,或者哪些参数对系统比较敏感,哪些对系统不敏感。
- 流动中最敏感区域是什么?
- 当我们改变平衡点、几何结构或参数时,稳定性/接受性等是如何变化的?
- 不稳定的“根源”在哪里?
- 我们应该在哪里放置控制元素来有效地操纵流?
比如临近失速的压气机系统,哪些位置的扰动会是系统立刻恶化,该位置用于布置传感器进行失速的检测。也可以对系统进行控制或者扩稳。
adjoint执行
在执行中,我们必须要获得Jabobian matrix
对于incompressible,已经有很多这方面的工作。但对于compressible,大家都还在努力。
ajoint可以采用automatic differentiation, 但目前效率可能不高。
adjoint的好处:10
正向:100个参数需要求解100次10000维 整个流场(n-s equ);逆向:2个参数需要求解2次10000维 整个流场(ajoint方程);两者得到同样的结果。
一个论文里面经常出现的案例:圆柱扰流
SU2中的adjoint11
可压RANS equation:
写出Langrange multiplier
变分,小扰动
第一项为原始状态下几何结构和标量函数值的变化,而第二项则依赖于原始几何和变形引起的标量函数的变化。
阻力或升力优化问题:
最终,变分可以写出:
第二项线性化:
linearized Navier-Stokes equations,
代回到原来的式子
利用分步积分,整理
如果不考虑边界条件的影响
最终,整理得到
至此,adjoint N-S equation可以得到:
0. Analysis of fluid systems: stability, receptivity, sensitivity ↩
1. On the structure and origin of pressureuctuations in wall turbulence: predictionsbased on the resolvent analysis ↩
2. A framework for studying the eect of compliant surfaces on wall turbulence ↩
3. Resolvent Analysis for Turbulent Channel Flow with Riblets ↩
4. On the existence of steady-state resonant waves in experiments ↩
5. Elements of resolvent methods in uid mechanics: notes for an introductory short course v0.3 ↩
6. Large Eddy Simulation for Compressible Flows ↩
7. phd thesis: Fluid dynamic instabilities in complexflow systems ↩
8. Computation of the blu -body sound generation by a self-consistent mean flow formulation ↩
9. https://www.youtube.com/watch?v=JFPTQ-ALDcM ↩
10. http://www2.eng.cam.ac.uk/~mpj1001/MJ_ABCD.html ↩
11. Viscous Continuous Adjoint Approach for the Design of Rotating Engineering Applications ↩