激波管的解法——有限差分方法

方程及状态参数

进排气管内一维变截面非定常控制方程为:
{ ∂ ( ρ A ) ∂ t + ∂ ( ρ u A ) ∂ x = 0 ∂ ( ρ u A ) ∂ t + ∂ ( ρ u 2 A + p A ) ∂ x = − ρ G A + p ∂ A ∂ x ∂ ( ρ e 0 A ) ∂ t + ∂ ( ρ u h 0 A ) ∂ x = ρ q ˙ A (1) \left\{ \begin{array}{l} \frac{ {\partial \left( {\rho A} \right)}}{ {\partial t}} + \frac{ {\partial \left( {\rho uA} \right)}}{ {\partial x}} = 0\\ \frac{ {\partial \left( {\rho uA} \right)}}{ {\partial t}} + \frac{ {\partial \left( {\rho {u^2}A + pA} \right)}}{ {\partial x}} = - \rho GA + p\frac{ {\partial A}}{ {\partial x}}\\ \frac{ {\partial \left( {\rho {e_0}A} \right)}}{ {\partial t}} + \frac{ {\partial \left( {\rho u{h_0}A} \right)}}{ {\partial x}} = \rho \dot qA \end{array} \right. \tag{1} t(ρA)+x(ρuA)=0t(ρuA)+x(ρu2A+pA)=ρGA+pxAt(ρe0A)+x(ρuh0A)=ρq˙A(1)
其中总热力学能为静热力学能和动能之和: e 0 = e + 1 2 u 2 = R g T γ − 1 + 1 2 u 2 = p ρ ( γ − 1 ) + 1 2 u 2 {e_0} = e + \frac{1}{2}{u^2} = \frac{ { {R_g}T}}{ {\gamma - 1}} + \frac{1}{2}{u^2} = \frac{p}{ {\rho \left( {\gamma - 1} \right)}} + \frac{1}{2}{u^2} e0=e+21u2=γ1RgT+21u2=ρ(γ1)p+21u2.
U = [ ρ A ρ u A ρ e 0 A ] {\bf{U}} = \left[ \begin{array}{l} \rho A\\ \rho uA\\ \rho {e_0}A \end{array} \right] U=ρAρuAρe0A F = [ ρ u A ( ρ u 2 + p ) A ρ h 0 u A ] {\bf{F}} = \left[ \begin{array}{l} \rho uA\\ \left( {\rho {u^2} + p} \right)A\\ \rho {h_0}uA \end{array} \right] F=ρuA(ρu2+p)Aρh0uA S = [ 0 p d A d x − ρ G A ρ q ˙ A ] {\bf{S}} = \left[ \begin{array}{l} 0\\ p\frac{ {dA}}{ {dx}} - \rho GA\\ \rho \dot qA \end{array} \right] S=0pdxdAρGAρq˙A,则式(1)化为:
∂ U ∂ t + ∂ F ( U ) ∂ x = S ( U ) (2) \frac{ {\partial {\bf{U}}}}{ {\partial t}} + \frac{ {\partial {\bf{F}}\left( {\bf{U}} \right)}}{ {\partial x}} = {\bf{S}}\left( {\bf{U}} \right)\tag{2} tU+xF(U)=S(U)(2)
在初始化网格热力状态时,一般给的都是 u , p , T u,p,T u,p,T和广义组分表征量 α \alpha α,然后求 U \bf{U} U,这里定义一个类Node,Node中的一个方法为 U {\bf{U}} U的初始化函数:

    def Uinit(self, u, p, T, AFAK=1.e8, A=1):
        from GasProperty import Rg
        from GasProperty import k_Justi
        self.p = p
        self.T = T
        self.rho = p / (Rg(AFAK) * T)
        self.u = u

        self.k = k_Justi(self.T, AFAK)
        self.e = self.p / self.rho / (self.k - 1.)
        self.e0 = self.e + self.u ** 2 / 2.
        self.h0 = self.e0 + self.p / self.rho
        self.a = pow(self.k * self.p / self.rho, 1. / 2.)

        from numpy import array
        self.U = array([self.rho * A, self.rho * self.u * A,
                        (self.p / (self.k - 1) + 1. / 2. * self.rho * self.u ** 2) * A]).transpose()

F ( U ) {\bf{F}}(\bf{U}) F(U)的表达式:
F ( U ) = [ U 2 U 2 2 U 1 + ( γ − 1 ) ( U 3 − U 2 2 2 U 1 ) U 2 ( U 3 + ( γ − 1 ) ( U 3 − U 2 2 2 U 1 ) ) U 1 ] (3) {\bf{F}({\bf{U}})} = \left[ \begin{array}{l} {U_2}\\ \frac{ { {U_2}^2}}{ { {U_1}}} + \left( {\gamma - 1} \right)\left( { {U_3} - \frac{ { {U_2}^2}}{ {2{U_1}}}} \right)\\ \frac{ { {U_2}\left( { {U_3} + \left( {\gamma - 1} \right)\left( { {U_3} - \frac{ { {U_2}^2}}{ {2{U_1}}}} \right)} \right)}}{ { {U_1}}} \end{array} \right]\tag{3} F(U)=U2U1U22+(γ1)(U32U1U22</

  • 1
    点赞
  • 7
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值