留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

高zeta势下Phan-Thien-Tanner(PTT)流体的电渗微推进器

郑佳璇 梁韵笛 菅永军

郑佳璇, 梁韵笛, 菅永军. 高zeta势下Phan-Thien-Tanner(PTT)流体的电渗微推进器[J]. 应用数学和力学, 2023, 44(10): 1213-1225. doi: 10.21656/1000-0887.430346
引用本文: 郑佳璇, 梁韵笛, 菅永军. 高zeta势下Phan-Thien-Tanner(PTT)流体的电渗微推进器[J]. 应用数学和力学, 2023, 44(10): 1213-1225. doi: 10.21656/1000-0887.430346
ZHENG Jiaxuan, LIANG Yundi, JIAN Yongjun. Electroosmotic Micro Thrusters of Phan-Thien-Tanner (PTT) Fluid at High Zeta Potential[J]. Applied Mathematics and Mechanics, 2023, 44(10): 1213-1225. doi: 10.21656/1000-0887.430346
Citation: ZHENG Jiaxuan, LIANG Yundi, JIAN Yongjun. Electroosmotic Micro Thrusters of Phan-Thien-Tanner (PTT) Fluid at High Zeta Potential[J]. Applied Mathematics and Mechanics, 2023, 44(10): 1213-1225. doi: 10.21656/1000-0887.430346

高zeta势下Phan-Thien-Tanner(PTT)流体的电渗微推进器

doi: 10.21656/1000-0887.430346
基金项目: 

国家自然科学基金项目 12262026

内蒙古自治区自然科学基金项目 2021MS01007

内蒙古自治区高校创新科研团队计划项目 NMGIRT2323

内蒙古自治区“草原英才”工程 12000-12102013

详细信息
    作者简介:

    郑佳璇(1992—),女,讲师,博士(E-mail: 31536015@mail.imu.edu.cn)

    梁韵笛(1996—),女,硕士生(E-mail: 31936003@mail.imu.edu.cn)

    通讯作者:

    菅永军(1974—),男,教授,博士生导师(通讯作者. E-mail: jianyj@imu.edu.cn)

  • 中图分类号: O357.1; O361.4

Electroosmotic Micro Thrusters of Phan-Thien-Tanner (PTT) Fluid at High Zeta Potential

  • 摘要: 在高的壁面zeta电势下, 考查了Phan-Thien-Tanner(PTT)黏弹性流体在平行板微通道中的电渗推进器问题. 在没有考虑Debye-Hückel线性近似的条件下, 求解了非线性Poisson-Boltzmann方程, 得到了高zeta电势下电势的解析解. 通过求解PTT流体满足的Cauchy动量方程, 获得了Navier滑移条件下微推进器速度的数值解. 进而通过数值积分得到了电渗微推进器的性能分布, 包括比冲、推力、效率和推力-功率比. 最后, 详细分析了黏弹性参数、壁面zeta电势、滑移系数和双电层厚度对速度分布及推进器性能的影响. 结果表明, 与Newton流体相比, PTT流体作为推进剂有利于推进器性能的提高, 比如, 流体速度随着黏弹性参数的增大而增大, 导致推进器性能也呈增大的趋势. 此外, 当前推进器比冲为800~1 000 ms时,推力可达0~250 μN, 效率为6%~12%, 推力-功率比为0~20 mN/W.
  • 近年来, 随着科学技术的进步与航天事业的发展, 复杂多样的空间活动日益增加, 促使微纳推进器广泛被关注. 与传统卫星相比, 微型卫星或航行器具有轻便、成本低、精度高等优势, 它们的小型化已成为现代航空事业发展的主要趋势. 其中, 微纳推进器作为小型卫星动力系统中空间推进的重要组成部分, 也逐渐向低质量和小体积的方向发展. 比如, 遥感卫星在绕地旋转时会形成一个旋转的轨道, 卫星实现自身功能的同时还有姿态调整的任务, 那么航天器的精准控制和轨道的校正就可以用微纳推进器去完成[1-2]. 此外, 电推进技术也已被广泛研究并应用于许多相关的空间任务中[3-7]. 而电渗就是一个基本的电动现象, 它能有效驱动和控制小型化系统中流体的流动. 在外加电场的作用下,流体被产生的电渗体力所驱动, 即电渗流(EOF)[8-9], 这与微纳管道内的双电层(EDL)[10]密切相关. 微纳管道内电渗流动的问题已被广泛探讨, 比如, 电渗与压力混合驱动的幂律流体特性[11]、电磁电渗驱动的流体流动[12]以及电渗驱动溶液浓度的扩散[13]等.

    基于电渗流具有实现微纳推进器高效率和低功耗的优势, 许多专家和学者将其作为电渗推进器中微/纳流体的运输流而研究[14-15]. Diez等[16]首先将纳米管道内流体的电动力学特征引入空间推进的研究中, 使微纳流体力学进入一个新的领域, 为之后研究微推进器提供了理论基础和研究参考. Huang等[17]将纳米电动推进器与两层流体的电渗流动相结合, 当存在壁面滑移条件时, 对电渗推进器的推进性能进行了研究. Zheng和Jian[18-19]考虑了在柔性纳米管道内空间电渗推进器的问题, 得到了推进器性能的解析表达式, 这些性能包括推力、比冲、效率以及推力-功率比, 它们是衡量推进器是否满足不同空间推进任务的关键因素. 随后, 他们数值研究了steric效应对柔性纳米管道内电渗推进器的影响[20], 并且发现由于Joule热效应的存在, 在考虑离子大小的实际情况下, 虽然推进器的推力可以提高, 但效率为30%~70%, 推力-功率比的峰值也会逐渐地消失. 然而, 在以上研究中, 推进器内的流体认为都是Newton流体, 因此, 非Newton黏弹性流体还没有在电渗推进器的研究中被考虑.

    大部分的生物流体是以非Newton黏弹性行为为特征的,具有长链分子结构的溶液, 常见的有血液、唾液、蛋白质溶液以及DNA溶液等. 具有模拟流体流动黏弹性的模型之一就是简化的Phan-Thien-Tanner(s-PTT)流体模型[21], 它与原始的PTT流体模型相比, 区别在于前者仅考虑分子的仿射运动, 这意味着简化的PTT流体在较小程度上显示剪切变薄的特性[22]. 近年来, 利用简化的PTT模型, 非Newton黏弹性流体的流变学特性被学者们广泛研究[15, 23-26]. Hashemabadi等[27]调查了通过两平行板间的简化PTT流体流动的问题, 并且考虑了本构方程中应力系数分为线性和指数形式的两种情况. 他们发现, 对于线性应力的情况, 弹性系数对流体的速度以及速度的极值有强烈影响. 进一步地, 在线性和非线性的Navier滑移边界条件下, Ferrás等[28]增加了二次型PTT黏弹性流体模型的数值研究, 并且他们获得了线性PTT流体速度的完全解析解. 此外, 一个黏弹性流体的电渗流模型被Sarma等[29]考虑, 他们利用简化的PTT本构关系描述了弹性流体的流变学行为, 并且在壁面处带有高zeta电势和水动力滑移的条件. 他们发现,较高的zeta电势和滑移系数导致体积流率明显增强.

    本文首次选取PTT流体, 即典型的黏弹性非Newton流体, 作为推进器中的推进剂, 并在高zeta电势下, 利用数值方法调查PTT流体的电渗微推进器的性能分布. 通过求解非线性的Poisson-Boltzmann方程, 得到完整的电势解析解. 采用广义的Cauchy动量方程, 根据Navier线性滑移边界条件数值地求解电渗流的速度分布. 从而, 由流体电势和速度所产生的微推进器的推力、比冲、效率和推力-功率比就可以数值地获得. 最后, 对系统的各个影响参数和推进器的性能进行分析, 主要目的是考察PTT流体在高zeta电势下对电渗微推进器性能的影响.

    在推进器的微管道中, 我们考虑一个稳定且完全发展的不可压缩黏弹性流体的电渗流动. 如图 1所示, 微管道的长和宽分别为LW, 高为2H, 并假设它的长和宽远远大于高, 即$L \gg 2 H, W \gg 2 H$, 因此, 管道内的流动可看作是在两无限平行的平板间的流动. 在通道中心建立直角坐标系, 并且假设上下壁面处的双电层不重叠[30]. 当在x轴方向上施加强度为Ex的外部电场时, 双电层中的离子被驱动而发生移动, 流体微团由于黏性作用而一起移动, 从而在整个微管道中形成电渗流. 此外, 在壁面处具有均匀的表面电荷, 即壁面处存在均匀的zeta电势. 本文采用简化的PTT模型来描述流体的流变行为, 并用线性Navier滑移定律作为壁面边界条件, 其中uw是滑移速度. 研究中仅考虑温度均匀且不施加外界压力的纯电渗流的流动.

    图  1  平行板微管道中黏弹性流体的纯电渗流
    Figure  1.  The pure electroosmotic flow of the viscoelastic fluid in a parallel plate microchannel

    控制流体流动的动力学方程是连续性方程和动量方程. 连续性方程表示为

    $$ \nabla \cdot \boldsymbol{u}=0, $$ (1)

    Cauchy动量方程为

    $$ \rho \frac{\mathrm{D} \boldsymbol{u}}{\mathrm{D} t}=\nabla \cdot \boldsymbol{\tau}+\rho_{e} \boldsymbol{E}, $$ (2)

    其中,u是速度矢量, ρρe分别是流体的密度和净电荷密度, τ是聚合物额外应力张量, E是电场.

    在研究电渗微推进器的性能之前, 首先应该考虑微管道内的电势分布. 电势Φ是两种电势贡献之和, 一部分为位于入口和出口由电极产生的电势ϕ, 由Laplace方程$\nabla^{2} \phi=0$描述, 另一部分是管道壁面附近自发诱导的电势ψ, 由如下Poisson方程给出[30-32]:

    $$ \nabla^{2} \psi=-\frac{\rho_{e}}{\varepsilon_{\mathrm{r}} \varepsilon_{0}}, $$ (3)

    其中, εrε0分别是聚合物溶液的相对介电常数和真空介电常数, 净电荷密度ρe服从Boltzmann分布[33]:

    $$\rho_{e}=-2 n_{0} e z \sinh \left(\frac{e z}{k_{\mathrm{b}} T} \psi\right), $$ (4)

    其中, n0是离子密度(n0 = CNA为聚合物溶液中离子的体积数浓度, C为离子的摩尔浓度, NA为Avogadro常数), kb是Boltzmann常数, T是绝对温度, e是基本点电荷, z是离子化合价. 由于流动是完全发展的、不重叠的、稳定的电渗流, 故方程(3)可简化为

    $$\frac{\mathrm{d}^{2} \psi}{\mathrm{d} y^{2}}=-\frac{\rho_{e}}{\varepsilon_{\mathrm{r}} \varepsilon_{0}}. $$ (5)

    联立方程(4)和(5), 有

    $$\frac{\mathrm{d}^{2} \psi}{\mathrm{d} y^{2}}=\frac{2_{0} e z}{\varepsilon_{\mathrm{r}} \varepsilon_{0}} \sinh \left(\frac{e z}{k_{\mathrm{b}} T} \psi\right). $$ (6)

    对应的边界条件为

    $$ \left.\frac{\mathrm{d} \psi}{\mathrm{d} y}\right|_{y=0}=0, $$ (7a)
    $$ \psi=\left.\zeta\right|_{y=H}. $$ (7b)

    对方程(6)积分, 并利用边界条件(7), 可得电势的解析解. 需要指出的是, 由于存在另一隐含边界条件$\psi=\left.\zeta\right|_{y=-H}$, 给出整个通道上完整的电势解析解为[16]

    $$ \psi=\frac{4 k_{\mathrm{b}} T}{e z}\left[\operatorname{arctanh}\left(\tanh \left(\frac{e z \zeta}{4 k_{\mathrm{b}} T}\right) \exp (-\kappa y-\kappa H)\right)+\operatorname{arctanh}\left(\tanh \left(\frac{e z \zeta}{4 k_{\mathrm{b}} T}\right) \exp (\kappa y-\kappa H)\right)\right], $$ (8)

    其中, $\kappa=\left(2 n_{0} e^{2} z^{2} /\left(\varepsilon_{\mathrm{r}} \varepsilon_{0} k_{\mathrm{b}} T\right)\right)^{1 / 2}$定义为双电层厚度的倒数. 为了便于计算, 我们引入以下无量纲量:

    $$y^{*}=\frac{y}{H}, \kappa^{*}=\kappa H, \left(\psi^{*}, \zeta^{*}\right)=\frac{(\psi, \zeta)}{k_{\mathrm{b}} T /(e z)} . $$ (9)

    故方程(8)的电势解析解可用无量纲的形式给出:

    $$ \psi^{*}=4\left[\operatorname{arctanh}\left(\tanh \left(\frac{\zeta^{*}}{4}\right) \exp \left(-\kappa^{*} y^{*}-\kappa^{*}\right)\right)+\operatorname{arctanh}\left(\tanh \left(\frac{\zeta^{*}}{4}\right) \exp \left(\kappa^{*} y^{*}-\kappa^{*}\right)\right)\right]. $$ (10)

    简化PTT流体的本构方程可以写为[34]

    $$ f(\operatorname{tr}(\boldsymbol{\tau})) \boldsymbol{\tau}+\lambda \stackrel{\nabla}{\boldsymbol{\tau}}=\eta\left(\nabla \boldsymbol{u}+\nabla \boldsymbol{u}^{\mathrm{T}}\right), $$ (11)

    其中, f是应力张量τ的迹函数, λ是弛豫时间(即聚合物分子在流动时重新分配所需的时间), η是黏度系数, $\stackrel{\nabla}{\tau}$是上随体导数,表示如下:

    $$\stackrel{\nabla}{\boldsymbol{\tau}}=\frac{\partial \boldsymbol{\tau}}{\partial t}+\boldsymbol{u} \cdot \nabla \boldsymbol{\tau}-\left[(\nabla \boldsymbol{u})^{\mathrm{T}} \cdot \boldsymbol{\tau}+\boldsymbol{\tau} \cdot \nabla \boldsymbol{u}\right]. $$ (12)

    函数f可以用线性和指数的形式表示, 本文选择线性形式[34], 具体如下:

    $$f(\operatorname{tr}(\boldsymbol{\tau}))=1+\frac{\varepsilon \lambda}{\eta} \operatorname{tr}(\boldsymbol{\tau}). $$ (13)

    其中, ε是拉伸参数.

    对于二维区域内完全发展的简化PTT流体的单向流动(即u只有x轴方向的分量), 本构方程简化为

    $$ f\left(\tau_{x x}+\tau_{y y}\right) \tau_{x x}=2 \lambda \tau_{x y}\left(\frac{\mathrm{d} u}{\mathrm{~d} y}\right), $$ (14)
    $$f\left(\tau_{x x}+\tau_{y y}\right) \tau_{y y}=0, $$ (15)
    $$ f\left(\tau_{x x}+\tau_{y y}\right) \tau_{x y}=\eta\left(\frac{\mathrm{d} u}{\mathrm{~d} y}\right)+\lambda \tau_{y y}\left(\frac{\mathrm{d} u}{\mathrm{~d} y}\right) . $$ (16)

    从方程(13)中可以看出, f(τxx+ τyy)=0或者τyy=0, 但如果前者等于零, 将得到不现实的结果, 因此唯一可能的结果是τyy=0[35]. 结合方程(13)—(16), 得到以下速度梯度与剪应力的关系:

    $$ \tau_{x y}+\frac{2 \varepsilon \lambda^{2}}{\eta^{2}}\left(\tau_{x y}\right)^{3}=\eta \frac{\mathrm{d} u}{\mathrm{~d} y}. $$ (17)

    此外, 本文是忽略压力梯度的纯电渗流问题, 因此Cauchy动量方程(2)化为

    $$ \frac{\mathrm{d} \tau_{x y}}{\mathrm{~d} y}-\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x} \frac{\mathrm{d}^{2} \psi}{\mathrm{d} y^{2}}=0. $$ (18)

    对方程(16)积分, 得

    $$ \tau_{x y}=\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x} \frac{\mathrm{d} \psi}{\mathrm{d} y}+A. $$ (19)

    假设双电层不重叠, 则微管道上下两壁面处的电势和滑移速度相同, 即边界条件对称, 故

    $$ \left.\frac{\mathrm{d} \psi}{\mathrm{d} y}\right|_{y=0}=0, \left.\tau_{x y}\right|_{y=0}=0. $$ (20)

    根据方程(19)和相应的边界条件(20), 得A=0. 将A代入方程(19), 有

    $$\tau_{x y}=\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x} \frac{\mathrm{d} \psi}{\mathrm{d} y}. $$ (21)

    结合方程(13)、(16)、(17)和(19), 得到流场的速度梯度为

    $$ \frac{\mathrm{d} u}{\mathrm{~d} y}=\frac{\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x}}{\eta} \frac{\mathrm{d} \psi}{\mathrm{d} y}+2 \varepsilon \lambda^{2}\left(\frac{\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x}}{\eta}\right)^{3}\left(\frac{\mathrm{d} \psi}{\mathrm{d} y}\right)^{3} . $$ (22)

    引入以下无量纲变量:

    $$ \tau_{x y}^{*}=\frac{\tau_{x y}}{\eta \kappa u_{\mathrm{sh}}}, u^{*}=\frac{u}{u_{\mathrm{sh}}}, $$ (23)

    其中, $u_{\mathrm{sh}}=-\varepsilon_{\mathrm{r}} \varepsilon_{0} E_{x} \zeta / \eta$是Helmholtz-Smoluchowski速度[36-38], 定义$D e=\lambda \kappa u_{\mathrm{sh}}$是基于双电层厚度的Deborah数. 将无量纲变量(9)和(23)代入方程(22), 确定微管道中流体的无量纲速度梯度:

    $$\frac{\mathrm{d} u^{*}}{\mathrm{~d} y^{*}}=-\frac{1}{\zeta^{*}} \frac{\mathrm{d} \psi^{*}}{\mathrm{~d} y^{*}}-\frac{2 \varepsilon D e^{2}}{\kappa^{* 2} \zeta^{* 3}}\left(\frac{\mathrm{d} \psi^{*}}{\mathrm{~d} y^{*}}\right)^{3} . $$ (24)

    在微纳尺度的管道中, 当双电层存在时, 流体在壁面处的滑移条件是不容忽视的, 壁面存在滑移有利于电渗流速度的增强, 这种增强源于Navier滑移边界条件[39]. 因此, 在非Newton流体的研究中, 选择线性Navier滑移定律来描述流体在固液接触面的滑移条件, 该定律定义的壁面滑移速度uw与壁面剪切应力$\tau_{y x, \mathrm{w}}$的关系如下[40-42]

    $$ u_{\mathrm{w}}=k_{\mathrm{nl}}\left(\mp \tau_{y x, \mathrm{w}}\right) \text {, } $$ (25)

    其中$k_{\mathrm{nl}}$是滑移系数, 符号“-”表示图 1中上壁面, “+”表示下壁面. 引入无量纲化的$u_{\mathrm{w}}^{*}=u_{\mathrm{w}} / u_{\mathrm{sh}}$和$k_{\mathrm{nl}}^{*}=$ $k_{\mathrm{nl}} \eta / H$, 则无量纲滑移速度可写为

    $$u_{\mathrm{w}}^{*}=k_{\mathrm{nl}}^{*}\left(\mp \kappa^{*} \tau_{y x, \mathrm{w}}^{*}\right) \text {. } $$ (26)

    因此, 在边界条件(26)下, 用数值计算方法可求得速度分布.

    基于几何形状的对称性, 本节只对微管道上半部分进行讨论, 即0≤yH.

    推进器的比冲定义为推进器每单位重量消耗推进剂所产生的反应脉冲量[43]. 假设流动和传热是稳定的, 且不考虑不稳态效应的情况下, 推进器的比冲直接由平均喷射速度决定, 而本文将推进剂排气速度视为平均流速, 故比冲可以用排气速度除以重力加速度常数来表示[18, 43-46]

    $$I_{\mathrm{sp}}=\frac{1}{g_{0} H} \int_{0}^{H} u(y) \mathrm{d} y. $$ (27)

    以上方程中的比冲通过$I_{\mathrm{sp}}^{*}=g_{0} I_{\mathrm{sp}} / u_{\mathrm{sh}}$无量纲化, 有

    $$ I_{\mathrm{sp}}^{*}=\int_{0}^{1} u^{*}\left(y^{*}\right) \mathrm{d} y^{*}. $$ (28)

    推进器是将任何形式的能量转化为机械能从而产生推力的装置. 在外加电场的驱动下, 流体产生速度使得推进器产生推力. 推力由流体的平均流速与质量流之间的关系给出, 可表示为速度的平方对y从0到H积分[18, 44-45]

    $$ T_{\mathrm{h}}=2 W \int_{0}^{H} \rho u^{2}(y) \mathrm{d} y. $$ (29)

    推力的无量纲表达式可写为

    $$ T_{\mathrm{h}}^{*}=\int_{0}^{1} u^{* 2}\left(y^{*}\right) \mathrm{d} y^{*}, $$ (30)

    其中$\quad T_{\mathrm{h}}^{*}=T_{\mathrm{h}} /\left(2 W H \rho u_{\mathrm{sh}}^{*}\right)$.

    推进器效率是指发动机提供的推进功率与产生的总功率之比. 在我们研究的纯电渗力驱动的情况下, 推进剂产生的动能(K)提供推进动力, 除动能外, 总能量还包括Joule热效应(Pj)、黏性耗散(Pv)和壁面处摩擦生热(Pf), 这里我们忽略流体内部的摩擦, 仅考虑固液接触面处的摩擦. 综上, 总功率写为

    $$ P_{\mathrm{in}}=K+P_{\mathrm{j}}+P_{\mathrm{v}}+P_{\mathrm{f}}, $$ (31)

    动能、Joule热效应、黏性耗散和摩擦热分别写为[18, 29, 32, 47]

    $$ K=\int_{A_{\text {out }}} \frac{1}{2} \rho u^{3} \mathrm{~d} A_{\text {out }}, $$ (32)
    $$ P_{\mathrm{j}}=\int_{V} \sigma E_{x}^{2} \mathrm{~d} V, $$ (33)
    $$P_{\mathrm{v}}=\int_{V} \tau_{x y} \frac{\mathrm{d} u}{\mathrm{~d} y} \mathrm{~d} V, $$ (34)
    $$ P_{\mathrm{f}}=\int_{A_{\text {wall }}} \mu(y=H) u(y=H) \frac{\mathrm{d} u}{\mathrm{~d} y}(y=H) \mathrm{d} A_{\text {wall }}, $$ (35)

    其中, Aout, AwallV分别代表微管道的横截面积、表面积和总体积, μ=τxy/(du/dy)是黏度剖面[29]. 为了使分析简单化, 以下是各种能量的无量纲表达式:

    $$ K^{*}=\int_{0}^{1} \frac{1}{2} u^{* 3} \mathrm{~d} y^{*}, $$ (36)
    $$ P_{\mathrm{j}}^{*}=\delta, $$ (37)
    $$P_{\mathrm{v}}^{*}=\gamma \int_{0}^{1} \tau_{x y}^{*} \frac{\mathrm{d} u^{*}}{\mathrm{~d} y^{*}} \mathrm{~d} y^{*}=\gamma P_{\mathrm{v}}^{\prime}, $$ (38)
    $$ P_{\mathrm{f}}^{*}=\gamma\left[\tau_{y x, \mathrm{w}}^{*} u^{*}\left(y^{*}=1\right)\right]=\gamma P_{\mathrm{f}}^{\prime}, $$ (39)
    $$P_{\mathrm{in}}^{*}=K^{*}+\delta+\gamma\left(P_{\mathrm{v}}^{\prime}+P_{\mathrm{f}}^{\prime}\right), $$ (40)
    $$ \delta=\frac{L \sigma E_{x}^{2}}{\rho u_{\mathrm{sh}}^{3}}, \gamma=\frac{L \eta \kappa^{*}}{\rho H^{2} u_{\mathrm{sh}}}, $$ (41)

    其中, $\left[P_{\mathrm{in}}^{*}, K^{*}, P_{\mathrm{j}}^{*}, P_{\mathrm{v}}^{*}, P_{\mathrm{f}}^{*}\right]=\left[P_{\mathrm{in}}, K, P_{j}, P_{\mathrm{v}}, P_{\mathrm{f}}\right] /\left(2 \rho H W u_{\mathrm{sh}}^{3}\right)$. 效率是动能与总输入功率的比值, 表示为

    $$ \eta_{\mathrm{t}}=\frac{K}{P_{\mathrm{in}}}=\frac{K^{*}}{P_{\mathrm{in}}^{*}}. $$ (42)

    在研究推进器性能的过程中, 另一个重要物理量是推力-功率比, 其定义为[44]

    $$ \xi=\frac{T_{\mathrm{h}}}{P_{\mathrm{in}}} \sim \frac{\rho H W v_{\mathrm{sh}}^{2}}{K / \eta_{\mathrm{t}}}=\eta_{\mathrm{t}} \frac{\rho H W v_{\mathrm{sh}}^{2}}{\rho H W v_{\mathrm{sh}}^{3} / 2}=\frac{2 \eta_{\mathrm{t}}}{v_{\mathrm{sh}}}, $$ (43)

    上式表示推力-功率比与效率成正比. 利用ξ*=ushξ将推力-功率比的无量纲形式写为如下形式:

    $$ \xi^{*}=\frac{T_{\mathrm{h}}^{*}}{P_{\mathrm{in}}^{*}} . $$ (44)

    本节中, 式(28)、(30)、(36)和(38)可通过数值积分方法求得.

    在第2节中我们得到了纯电渗作用下黏弹性PTT流体速度的数值解, 并利用数值积分方法得到推进器性能分布. 我们将对推进器的速度、比冲、推力、效率和推力-功率比进行讨论, 调查Deborah数、滑移系数和zeta电势对推进器性能的影响. 微管道的半高、宽度和长度分别设置为H=5 μm [15, 48], W=200 μm和L=1 mm[43], 流体密度为ρ=103 kg/m3[15], 电导率为σ=10-3 S/m[48], 无量纲滑移系数为0≤knl*≤1[32], 无量纲zeta电势为1≤ζ*≤10[29]. 其中, 含有Deborah数的εDe2称为黏弹性基团, 它是描述黏性和弹性行为的重要物理参数, 取值范围是0≤εDe2 ≤1[29]. 当对整个系统施加500 V电压时, x轴方向的电场强度为Ex=5×105 V/m. 其他物理量使用以下典型值: 10≤κ* ≤100, ε0=8.854×10-12 C/(V · m), εr=80, T=300 K, kb=1.38×10-23, z=1, e=1.6×10-19 C[29].

    为了验证本文速度数值解的正确性, 图 2首先给出当前研究的速度数值解与Sarma等[29]的速度解析解的比较, 显然本文速度的数值解与Sarma等[29]的解析解非常吻合. 图 3表示电渗微推进器在黏弹性参数εDe2取不同值时所对应的有量纲速度分布,其中κ=4 μm-1, knl=5×10-3 m/(Pa · s), ζ=150 mV. 当εDe2=0时, 所示结果为Newton流体的速度分布, 当εDe2小于0.1时, 速度仅有微小变化, 也可近似视为Newton流体[49], 比如水. 当εDe2为0.5或1时, 所示结果为非Newton PTT流体的速度分布, 在图 3中流体速度明显增加, 这是因为随着黏弹性参数εDe2的增加, 流体的剪切变薄行为增强, 即黏弹性流体的黏度降低, 使得在相同条件下, PTT流体的流速比Newton流体速度更大. 因此,与Newton流体相比,PTT流体作为推进剂,其速度的增加有助于推进器性能的提高(图 4).

    图  2  速度的数值解与Sarma等[29]速度的解析解对比
    Figure  2.  The current numerical velocity compared to the analytical velocity of Sarma et al. [29]
    图  3  对于不同的εDe2, 电渗微推进器的速度分布
    Figure  3.  The velocity distributions of the electroosmotic microthruster for different values of εDe2
    图  4  对于不同的εDe2值, 推进器性能随κ的变化情况
    Figure  4.  The thruster performance changes with κ for different values of εDe2

    图 4中, 当εDe2=0,0.1,0.5,1时, 我们绘制了微推进器的性能分布随双电层厚度倒数κ的变化图像, 这些性能为比冲Isp、推力Th、效率ηt以及推力-功率比ξ(此处及之后的讨论中, 要求κ >2 μm-1κ* >10, 否则不满足薄的双电层的假设). 其中knl=5×10-3 m/(Pa · s), ζ=150 mV. 当固定双电层厚度的倒数κ时, 推进器的比冲、推力、效率和推力-功率比都随着黏弹性参数值εDe2的增大而增大. 在图 4中, 当εDe2小于0.1时, 视为Newton流体的性能曲线, 显然PTT流体的性能高于Newton流体的性能, 其中效率和推力-功率比的增长更为明显. 这是由于εDe2的增大, 意味着黏性和弹性的比值变小, 所受阻力减小, 使得流速增大(图 3), 而比冲、推力、动能分别与速度、速度的二次方、速度的三次方成正比, 效率与动能成正比, 推力-功率比又与效率成正比,故推进器的性能与速度变化趋势保持一致. 另一方面, 当固定εDe2时,随着κ的增大, 双电层厚度变小, 电荷密度增大, 致使离子受到较大的电场力而使电渗透速度增大, 进而使得比冲、推力、效率和推力-功率比增加.

    图 5显示的是不同壁面zeta电势ζ对电渗微推进器速度分布的影响,其中, κ=2 μm-1, knl=1×10-2 m/(Pa · s), εDe2=0.5. 正如我们所期望的那样, 固定其他参数的值, 当壁面zeta电势不断增大时电渗流速度也增大. 这是因为当电势变大时, 电双层内的电荷密度变大, 离子浓度增高, 流体粒子具有较高的电渗驱动力, 这实际上增加了流体的速度[29]. 在不同zeta电势下, 图 6为推进器的性能随双电层厚度的倒数κ的变化分布,其中knl=1×10-2 m/(Pa · s), εDe2=0.5. 如图 6所示, 当固定κ时, 随着壁面zeta电势ζ的增大, 推进器的比冲、推力、效率和推力-功率比都明显提高. 这归因于图 5所示的随着zeta电势的增大流体速度变大, 从而比冲、推力、动能也随zeta电势的增大而增大. 同样地, 上面已经讨论过推进器效率与推力-功率比成正相关, 故它也随zeta电势的增加而增加.

    图  5  在不同的zeta电势下, 电渗微推进器的速度分布
    Figure  5.  The velocity distributions of the electroosmotic microthruster for different zeta values
    图  6  在不同zeta电势下, 推进器性能随κ的变化情况
    Figure  6.  The thruster performance changes with κ for different zeta values

    图 7描述了滑移系数knl和双电层厚度的倒数κ对微管道内速度分布的影响,其中κ=3 μm-1, ζ=125 mV, εDe2=0.5. 所有参数值都在所绘图中给出. 壁面摩擦力的大小可量化为流体在壁面处的Navier滑移条件, 当滑移系数为0即壁面是无滑移时, 壁面处的摩擦力较大导致流体流动的阻力强, 故流速最小, 如图 7所示. 随着滑移系数的增大, 流体速度明显增大, 这意味着界面滑移的存在减小了流动阻力使得流速增加. 图 8绘制了在滑移系数knl取不同值时推进器性能随κ的变化. 从图 8可以看出, 当固定参数κ时, 随着滑移系数knl的增大, 推进器的比冲、推力、效率和推力-功率比也明显增大, 这主要归因于速度的增加. 其中ζ=125 mV, εDe2=0.5. 特别地, 当knl=0即壁面是无滑移时, 推进器的比冲、推力、效率和推力-功率比都近乎为0, 这主要因为无滑移时流体速度很小所导致的,因此壁面存在滑移有利于推进器性能的提高.

    图  7  当滑移系数knl取不同值时, 电渗微推进器的速度分布
    Figure  7.  The velocity distributions of the electroosmotic microthruster for different slip coefficients
    图  8  当滑移系数knl取不同值时, 推进器的性能随κ的变化图象
    Figure  8.  The thruster performance changes with κ for different slip coefficients

    本文基于PTT流体模型描述流体的流变行为, 从理论上研究电渗微推进器在空间推进中的应用. 求解了完整的Poisson-Boltzmann方程和Cauchy动量方程, 同时界面水动力滑移使用线性Navier滑移定律建模, 得到电渗微推进器电势的完整解析解和速度的数值解. 利用电势解析解和数值速度揭示双电层厚度、黏弹性参数、zeta电势和滑移系数对推进器性能的影响, 即考查推进器比冲、推力、效率、推力-功率比随各个参数的化变情况. 在所考虑的微推进器模型中, 与Newton流体相比, PTT流体作为推进剂使得推进器的性能有所提高, 其中效率和推力-功率比的增长较为明显,这是由于黏弹性参数的增加产生了更大的流速. 此外, 增加壁面zeta电势和滑移系数, 流速和推进器的各项性能显著提高. 在所选取的参数范围内, 推进器比冲为800~ 1 000 ms时, 推力可达0~250 μN, 效率为6%~12%, 推力-功率比为0~20 mN/W. 我们认为, 本分析为目前黏弹性流体作为推进剂的微推进器提供了深入的理论见解, 而真正应用到太空任务的操作与运行中, 仍需要根据实际情况对每个参数进行优化.

    致谢: 本文作者衷心感谢内蒙古师范大学基本科研业务费(2022JBYJ022)对本文的资助.
  • 图  1  平行板微管道中黏弹性流体的纯电渗流

    Figure  1.  The pure electroosmotic flow of the viscoelastic fluid in a parallel plate microchannel

    图  2  速度的数值解与Sarma等[29]速度的解析解对比

    Figure  2.  The current numerical velocity compared to the analytical velocity of Sarma et al. [29]

    图  3  对于不同的εDe2, 电渗微推进器的速度分布

    Figure  3.  The velocity distributions of the electroosmotic microthruster for different values of εDe2

    图  4  对于不同的εDe2值, 推进器性能随κ的变化情况

    Figure  4.  The thruster performance changes with κ for different values of εDe2

    图  5  在不同的zeta电势下, 电渗微推进器的速度分布

    Figure  5.  The velocity distributions of the electroosmotic microthruster for different zeta values

    图  6  在不同zeta电势下, 推进器性能随κ的变化情况

    Figure  6.  The thruster performance changes with κ for different zeta values

    图  7  当滑移系数knl取不同值时, 电渗微推进器的速度分布

    Figure  7.  The velocity distributions of the electroosmotic microthruster for different slip coefficients

    图  8  当滑移系数knl取不同值时, 推进器的性能随κ的变化图象

    Figure  8.  The thruster performance changes with κ for different slip coefficients

  • [1] SMITH R S, HADAEGH F Y. Control of deep-space formation-flying spacecraft; relative sensing and switched information[J]. Journal of Guidance, Control, and Dynamics, 2005, 28(1): 106-114. doi: 10.2514/1.6165
    [2] 徐晓辉. 深空探测中的微型推进器技术[J]. 机电产品开发与创新, 2007, 20(6): 29-31. https://www.cnki.com.cn/Article/CJFDTOTAL-JDCP200706013.htm

    XU Xiaohui. Micro thruster technology in deep-space exploring[J]. Development & Innovation of Machinery & Electrical Products, 2007, 20(6): 29-31. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-JDCP200706013.htm
    [3] BLANCO A, ROY S. Rarefied gas electro jet (RGEJ) micro-thruster for space propulsion[J]. Journal of Physics D: Applied Physics, 2017, 50(45): 455201. doi: 10.1088/1361-6463/aa8c47
    [4] MARTINEZ-SANCHEZ M, POLLARD J E. Spacecraft electric propulsion-an overview[J]. Journal of Propulsion and Power, 1998, 14(5): 688-699. doi: 10.2514/2.5331
    [5] TAKAHASHI T, TAKAO Y, ERIGUCHI K, et al. Microwave-excited microplasma thruster: a numerical and experimental study of the plasma generation and micronozzle flow[J]. Journal of Physics D: Applied Physics, 2008, 41(19): 194005. doi: 10.1088/0022-3727/41/19/194005
    [6] KEIDAR M, ZHUANG T, SHASHURIN A, et al. Electric propulsion for small satellites[J]. Plasma Physics and Controlled Fusion, 2014, 57(1): 014005.
    [7] 朴思扬, 朱春艳, 褚福运, 等. 基于等效梁模型的运载火箭动力学特性仿真预示研究[J]. 应用数学和力学, 2020, 41(3): 280-291. doi: 10.21656/1000-0887.400256

    PIAO Siyang, ZHU Chunyan, CHU Fuyun, et al. Dynamic characteristics prediction of launch vehicles based on the equivalent beam model[J]. Applied Mathematics and Mechanics, 2020, 41(3): 280-291. (in Chinese) doi: 10.21656/1000-0887.400256
    [8] LYKLEMA J. Fundamentals of Interface and Colloid Science: Soft Colloids[M]. Elsevier, 2005.
    [9] QIAO R. Control of electroosmotic flow by polymer coating: effects of the electrical double layer[J]. Langmuir, 2006, 22(16): 7096-7100. doi: 10.1021/la060883t
    [10] 解智勇. 两层微流体系统中的电动流动及传热分析[D]. 博士学位论文. 呼和浩特: 内蒙古大学, 2019.

    XIE Zhiyong. The analysis of electrokinetic flow and heat transfer characteristic through a two-layer microfluidic system[D]. PhD Thesis. Hohhot: Inner Mongolia University, 2019. (in Chinese)
    [11] 罗艳, 李鸣, 杨大勇. 微通道内电渗压力混合驱动幂律流体流动模拟[J]. 应用数学和力学, 2016, 37(4): 373-381. doi: 10.3879/j.issn.1000-0887.2016.04.005

    LUO Yan, LI Ming, YANG Dayong. Simulation of mixed electroosmotic and pressure-driven flows of power-law fluids in microchannels[J]. Applied Mathematics and Mechanics, 2016, 37(4): 373-381. (in Chinese) doi: 10.3879/j.issn.1000-0887.2016.04.005
    [12] 王爽, 菅永军. 周期壁面电势调制下平行板微管道中的电磁电渗流动[J]. 应用数学和力学, 2020, 41(4): 396-405. doi: 10.21656/1000-0887.400151

    WANG Shuang, JIAN Yongjun. Magnetohydrodynamic electroosmotic flow in zeta potential patterned micro-parallel channels[J]. Applied Mathematics and Mechanics, 2020, 41(4): 396-405. (in Chinese) doi: 10.21656/1000-0887.400151
    [13] 张凯, 林建忠, 李志华. 电渗驱动微通道流中的扩散[J]. 应用数学和力学, 2006, 27(5): 512-518. http://www.applmathmech.cn/article/id/714

    ZHANG Kai, LIN Jianzhong, LI Zhihua. Diffusion in the micro-channel flow driven by electroosmosis[J]. Applied Mathematics and Mechanics, 2006, 27(5): 512-518. (in Chinese) http://www.applmathmech.cn/article/id/714
    [14] JI J, QIAN S, LIU Z. Electroosmotic flow of viscoelastic fluid through a constriction microchannel[J]. Micromachines, 2021, 12(4): 417. doi: 10.3390/mi12040417
    [15] MARTÍNEZ L, BAUTISTA O, ESCANDÓN J, et al. Electroosmotic flow of a Phan-Thien-Tanner fluid in a wavy-wall microchannel[J]. Colloids and Surfaces A: Physicochemical and Engineering Aspects, 2016, 498: 7-19.
    [16] DIEZ F J, HERNAIZ G, MIRANDA J J, et al. On the capabilities of nano electrokinetic thrusters for space propulsion[J]. Acta Astronautica, 2013, 83: 97-107. doi: 10.1016/j.actaastro.2012.09.020
    [17] HUANG K H, HUANG H F. Two-liquid electroosmotic thrusters for micro propulsion applications[J]. Physics of Fluids, 2019, 31(12): 122003. doi: 10.1063/1.5128274
    [18] ZHENG J, JIAN Y. Electroosmotic thrusters in soft nanochannels for space propulsion[J]. Physics of Fluids, 2020, 32(12): 122005. doi: 10.1063/5.0033436
    [19] ZHENG J, JIAN Y. Space electroosmotic thrusters in ion partitioning soft nanochannels[J]. Micromachines, 2021, 12: 777. doi: 10.3390/mi12070777
    [20] ZHENG J, JIA B, JIAN Y. Steric effects on space electroosmotic thrusters in soft nanochannels[J]. Mathematics, 2021, 9(16): 1916. doi: 10.3390/math9161916
    [21] PHAN-THIEN N, TANNER R I. A new constitutive equation derived from network theory[J]. Journal of Non-Newtonian Fluid Mechanics, 1977, 93(2/3): 353-365.
    [22] OLIVEIRA P J, PINHO F T. Analytical solution for fully developed channel and pipe flow of Phan-Thien-Tanner fluids[J]. Journal of Non-Newtonian Fluid Mechanics, 1999, 387: 271-280. http://web.fe.up.pt/~fpinho/pdfs/jfma.pdf
    [23] PINHO F T, OLIVEIRA P J. Analysis of forced convection in pipes and channels with the simplified Phan-Thien-Tanner fluid[J]. International Journal of Heat and Mass Transfer, 2000, 43(13): 2273-2287. doi: 10.1016/S0017-9310(99)00303-8
    [24] PINHO F T, OLIVEIRA P J. Axial annular flow of a nonlinear viscoelastic fluid: an analytical solution[J]. Journal of Non-Newtonian Fluid Mechanics, 2000, 93(2/3): 325-337.
    [25] CRUZ D O A, PINHO F T. Skewed Poiseuille-Couette flows of sPTT fluids in concentric annuli and channels[J]. Journal of Non-Newtonian Fluid Mechanics, 2004, 121(1): 1-14. doi: 10.1016/j.jnnfm.2004.03.007
    [26] 许晓阳, 彭严, 邓方安. 三维PTT黏弹性液滴撞击固壁面问题的改进SPH模拟[J]. 应用数学和力学, 2015, 36(6): 616-627. doi: 10.3879/j.issn.1000-0887.2015.06.006

    XU Xiaoyang, PENG Yan, DENG Fang'an. Numerical simulation of 3D PTT droplet impact onto solid surface with an improved smoothed particle hydrodynamics method[J]. Applied Mathematics and Mechanics, 2015, 36(6): 616-627. (in Chinese) doi: 10.3879/j.issn.1000-0887.2015.06.006
    [27] HASHEMABADI S H, ETEMAD S G, THIBAULT J, et al. Analytical solution for dynamic pressurization of viscoelastic fluids[J]. International Journal of Heat and Fluid Flow, 2003, 24(1): 137-144. doi: 10.1016/S0142-727X(02)00204-7
    [28] FERRÁS L L, AFONSO A M, ALVES M A, et al. Annular flow of viscoelastic fluids: analytical and numerical solutions[J]. Journal of Non-Newtonian Fluid Mechanics, 2014, 212: 80-91. doi: 10.1016/j.jnnfm.2014.07.004
    [29] SARMA R, DEKA N, SARMA K, et al. Electroosmotic flow of Phan-Thien-Tanner fluids at high zeta potentials: an exact analytical solution[J]. Physics of Fluids, 2018, 30(6): 062001. doi: 10.1063/1.5033974
    [30] PARDON G, VAN DER WIJNGAART W. Modeling and simulation of electrostatically gated nanochannels[J]. Advances in Colloid and Interface Science, 2013, 199/200: 78-94. doi: 10.1016/j.cis.2013.06.006
    [31] FERRÁS L L, AFONSO A M, ALVES M A, et al. Electro-osmotic and pressure-driven flow of viscoelastic fluids in microchannels: analytical and semi-analytical solutions[J]. Physics of Fluids, 2016, 28(9): 093102. doi: 10.1063/1.4962357
    [32] ANAND V. Effect of slip on heat transfer and entropy generation characteristics of simplified Phan-Thien-Tanner fluids with viscous dissipation under uniform heat flux boundary conditions: exponential formulation[J]. Applied Thermal Engineering, 2016, 98: 455-473. doi: 10.1016/j.applthermaleng.2015.12.025
    [33] BRUUS H. Theoretical Microfluidics[M]. Oxford: Oxford University Press, 2008.
    [34] PHAN-THIEN N. A nonlinear network viscoelastic model[J]. Journal of Rheology, 1978, 22(3): 259-283. doi: 10.1122/1.549481
    [35] FERRÁS L L, NÓBREGA J M, PINHO F T. Analytical solutions for channel flows of Phan-Thien-Tanner and Giesekus fluids under slip[J]. Journal of Non-Newtonian Fluid Mechanics, 2012, 171: 97-105.
    [36] SMOLUCHOWSKI M. Versuch einer mathematischen theorie der koagulationskinetik kolloider lösungen[J]. Zeitschrift für Physikalische Chemie, 1918, 92(1): 129-168.
    [37] PARK H M, LEE W M. Helmholtz-Smoluchowski velocity for viscoelastic electroosmotic flows[J]. Journal of Colloid and Interface Science, 2008, 317(2): 631-636. doi: 10.1016/j.jcis.2007.09.027
    [38] PROBSTEIN R F. Physicochemical Hydrodynamics[M]. Wiley Interscience, 2003.
    [39] EIJKEL J C T. Liquid slip in micro- and nanofluidics: recent research and its possible implications[J]. Lab on a Chip, 2007, 7(3): 299-301. doi: 10.1039/b700364c
    [40] SCHOWALTER W R. The behavior of complex fluids at solid boundaries[J]. Journal of Non-Newtonian Fluid Mechanics, 1988, 29: 25-36. doi: 10.1016/0377-0257(88)85048-1
    [41] NAVIER C. Mémoire sur les lois du mouvement des fluides[J]. Mémoires de l'Académie Royale des Sciences de l'Institut de France, 1823, 6(1823): 389-440.
    [42] FERNANDES C, FERRÁS L L, HABLA F, et al. Implementation of partial slip boundary conditions in an open-source finite-volume-based computational library[J]. Journal of Polymer Engineering, 2019, 39(4): 377-387.
    [43] 张雪儿, 李得天, 张天平. 电推进三种比冲的定义及其工程应用[J]. 真空与低温, 2020, 26(6): 486-493. https://www.cnki.com.cn/Article/CJFDTOTAL-ZKDW202006011.htm

    ZHANG Xueer, LI Detian, ZHANG Tianping. Three types of specific impulses for electric propulsions: their definitions and engineering applications[J]. Vacuum & Cryogenics, 2020, 26(6): 486-493. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-ZKDW202006011.htm
    [44] GOEBEL D M, KATZ I. Fundamentals of Electric Propulsion: Ion and Hall Thrusters[M]. John Wiley & Sons, 2008.
    [45] HEY F G. Micro Newton Thruster Development: Direct Thrust Measurements and Thruster Downscaling[M]. Springer, 2018.
    [46] JONCQUIERES V, VERMOREL O, CUENOT B. A fluid formalism for low-temperature plasma flows dedicated to space propulsion in an unstructured high performance computing solver[J]. Plasma Sources Science and Technology, 2020, 29(9): 095005.
    [47] ESCANDÓN J P, BAUTISTA O, MÉNDEZ F, et al. Theoretical conjugate heat transfer analysis in a parallel flat plate microchannel under electro-osmotic and pressure forces with a Phan-Thien-Tanner fluid[J]. International Journal of Thermal Sciences, 2011, 50(6): 1022-1030.
    [48] ESCANDÓN J, SANTIAGO F, BAUTISTA O, et al. Hydrodynamics and thermal analysis of a mixed electromagnetohydrodynamic-pressure driven flow for Phan-Thien-Tanner fluids in a microchannel[J]. International Journal of Thermal Sciences, 2014, 86: 246-257.
    [49] DHINAKARAN S, AFONSO A M, ALVES M A, et al. Steady viscoelastic fluid flow between parallel plates under electro-osmotic forces: Phan-Thien-Tanner model[J]. Journal of Colloid and Interface Science, 2010, 344(2): 513-520.
  • 期刊类型引用(3)

    1. 长龙,布仁满都拉,娜仁,孙艳军,菅永军. pH调节的纳米平行通道中Powell-Eyring流体的电渗流动. 应用数学和力学. 2025(01): 72-83 . 本站查看
    2. 长龙,布仁满都拉,孙艳军,刘全生,菅永军. 具有正弦波纹的环形微通道中电渗流动. 内蒙古大学学报(自然科学版). 2025(01): 26-36 . 百度学术
    3. 长龙,布仁满都拉,孙艳军,菅永军. 具有正弦波纹的平行板微通道中Jeffrey流体周期电渗流动. 应用数学和力学. 2024(05): 622-636 . 本站查看

    其他类型引用(0)

  • 加载中
图(8)
计量
  • 文章访问数:  520
  • HTML全文浏览量:  198
  • PDF下载量:  54
  • 被引次数: 3
出版历程
  • 收稿日期:  2022-10-31
  • 修回日期:  2023-03-13
  • 刊出日期:  2023-10-31

目录

/

返回文章
返回