引言

在临床治疗中,植入自体细胞（考虑到组织的相容性）进行组织修复和重建常被认为是最佳的治疗方法。考虑到自身组织移植会造成二次创伤,获取大量自体细胞最直接的途径就是进行细胞的体外培养。一些研究^[1,2]发现适宜的物理微环境是细胞生长不可或缺的因素,在体外培养时生物灌流反应器（如图1）可以方便地模拟体内多种物理微环境。许多实验室引入了不同的生物反应器进行细胞培养,并证明了液流有利于组织的构建^[3]。相比静态培养,微流控培养腔中的细胞在移植后会对体内的微流动环境具有更强的适应性。人体的一些组织结构（例如骨）在受到外力时会产生变形,进而在组织内微液流环境中形成孔隙压力和流动电位,这些力、电信号微环境进一步调控当地组织细胞的生长、分化,最终使组织（重/塑建）适应外部载荷环境。在体外细胞培养生物反应器上也可施加压力和电场进行驱动来模拟体内细胞所处的力-电微环境。很多学者的工作^[1,4]也已经证明压力和电场的干预可以调控细胞的生长、分化,虽然明确的作用机理（通路）还不是很清楚,但已明确的事实是压力梯度和电场驱动液体流动产生的流体切应力对细胞生长、分化起到了至关重要的作用。

细胞可以感受不同大小的切应力并能做出相应的响应^[5,6,7],能较好的将这些细胞能感受的力电信号和宏观的外载荷联系起来,微流动场被认为是信号传导的最佳媒介。Chang等^[8]得到了牛顿型液体在电场驱动下的流速分布,Ganguly等^[9]考虑了电和磁对压力驱动流传热特性的影响,Movahed等^[10]对比了压力和电场驱动小型设备中液流的优劣。这些研究背景均基于微流控技术,以讨论外力和流速的关系为主,微流控技术以其用料少、操作性强和便于观察分析等优点,近年来受到众多领域学者的青睐。在生物医学领域,Pappas^[11]在对癌细胞的生长和扩散进行分析时,利用微流控技术建立可控规模的反应区评估癌细胞和器官的相互作用。Uzel等^[12]在微流控芯片上设定不同的生化梯度,用于模拟体内细胞所处的化学微环境。Sun^[13]则将微流控芯片应用到细胞的筛选工作,并通过电场实现微阀控制。这些研究均以实验和有限元建模为主要手段,并未建立系统的理论模型进行具体量化分析,并探明如下调控机制：外部物理驱动场→微流体流动行为（流体切应力等）→细胞的生长/分化。

在细胞动态流动培养腔中,细胞多用粘附分子（整合素）将自身吸附在培养腔壁上^[14]（避免被液流带走）,这样使得液流对细胞的力学刺激作用集中在腔壁附近。Becquart等^[15]研究发现,处于液体环境中的细胞能感知静水压力及流体剪切力的刺激信号,并且发现流体剪切力比静水压力对引起相关基因表达的诱导作用更大。Stavenschi等^[6]在实验中利用压力驱动对比了不同大小和频率的切应力下细胞的基因表达效果,实验表明2Pa、2Hz的切应力刺激最有利于成骨基因的表达。Zhang等^[16]也通过实验发现了人类间充质干细胞（MSCs）基因的表达和流体切应力有极大的关系,这个信号传导过程也受到供体变异性的影响。此外,Zheng等^[7]同样通过实验探究出压力驱动作用下的切应力大小影响细胞对刺激的反应时间,并发现切应力较大时有利于缩短该反应时间。液体压力梯度驱动流对细胞培养的影响的研究大多停留在实验层面,理论方面鲜有报道。当使用压力梯度驱动液体在较长的细胞培养腔流动时,力信号的传递从施力点（管口）到管段中央具有延迟性。Glawdel^[17]将电场驱动引入到细胞的体外培养,因为电场力属于体力,这样就可以避免上述压力梯度驱动力信号沿着腔体长度方向的耗散延迟,从而保证腔内各个位置细胞附近的流体切应力信号均匀。同细胞培养腔中的压力驱动技术类似,电场驱动或力-电协同驱动研究大都以实验为主,例如王淞等^[18]通过实验观察得出电场对表皮干细胞的增殖有影响,内源性生物电场具有引导其向阴极定向迁移的作用,Kumar等^[19]在压力驱动流中添加电场作用后细胞的成骨响应明显增强,而相关的力-电协同驱动下的细胞培养腔理论研究较少。

考虑到人体内复杂的生理载荷环境,压力和电场在组织内并存。因此,在体外进行细胞培养时应根据细胞的实际需要组合施加相应的物理场（压力梯度或电场）,使其所处的环境更接近于生理物理环境。以上研究均没有得到具体的外部物理驱动场与相应培养腔内液体流动行为的量化关系,而这是进一步研究流体刺激特别是力-电协同刺激细胞生长或分化机理的前提。为此,本文将分别建立压力梯度、电场和力-电协同三种驱动模式下的矩形细胞培养流动腔（平行板培养皿）理论模型,求解得到相应的流速、切应力和流率的解析解,揭示出培养腔中的液体流动规律,并将这些解答和实际的细胞培养等生理意义联系起来。该模型的建立有助于细胞微流控培养实验系统及细胞力学参数的优化设计。

1模型建立

假设培养液为不可压缩的牛顿流体,在培养腔内的流动可以用下面的式子（动量方程和连续性方程）来描述^[20,21]：

$\rho (\frac{\partial u}{\partial t}+u\cdot \nabla u)=-\nabla p+f+\mu \Delta u$ (1)

$\nabla \cdot u=0$ (2)

其中,u是流体速度矢量,p表示液体压强,f是流体所受的体力,ρ为细胞培养液的密度,μ是培养液的动力粘度。∂u/∂t和u·▽u分别是瞬态项和对流项,△为laplace算子,▽为梯度算子。对于生物反应器微流控芯片（如图1所示）上的矩形细胞培养流动腔（平行板培养皿）,忽略腔宽的影响,在直角坐标系（x, y）下,培养腔内的流体可以被简化为二维流动模型（上下腔壁间高为2h,腔长为l）,在固液界面上的边界条件为无滑移。

显示原图| 下载原图ZIP| 生成PPT

图1   细胞培养腔示意图

Fig. 1   Diagram of cell culture chamber

1.1压力梯度驱动流模型的建立

根据不可压缩条件,N-S动量方程(1)中的对流项将消失。由腔管的对称性,只需对矩形微腔管中心线的上半部分进行研究。由于微管通道宽度远大于高度,因此可以忽略矩形微通道的宽度对沿高度方向上的流速分布的影响。若微管高度y方向的速度为零, 只考虑x方向的流动,压强p仅是关于x和时间变量t的函数,忽略体力的影响,压力梯度驱动液流在xy截面上的N-S动量方程可以表示为：

$\rho \frac{\partial {{u}^{p}}(y,t)}{\partial t}=\mu \frac{\partial }{\partial y}(\frac{\partial {{u}^{p}}(y,t)}{\partial y})-\frac{\partial p}{\partial x}{{e}^{i{{\omega }^{p}}t}}$ (3)

其中,u^p(y, t)是压力梯度驱动时x方向的流速分量,∂p/∂x和ω^p分别是沿着x方向的压力梯度和振荡角频率（ω^p=2πf ^p）,其中-∂p/∂x=P₀/l,负号“-”表示压力沿流道方向下降,P₀是进出口压力差的幅值。上述方程满足边界条件：

${{u}^{p}}(y,t){{|}_{y=h}}=0\begin{matrix}{} & {} \\\end{matrix}\frac{\partial {{u}^{p}}(y,t)}{\partial y}{{|}_{y=0}}=0$ (4)

引入无量纲化参数：

$\begin{equation} \left. \begin{aligned} \overline{y}=\frac{y}{h}\\ \overline{t}=\frac{t}{1/{{\omega }^{p}}}\\ {{\overline{u}}^{p}}(\overline{y},\overline{t})=\mu {{u}^{p}}(y,t)\\ R{{e}^{p}}=\frac{\rho {{\omega }^{p}}{{h}^{2}}}{\mu } \end{aligned} \right\} \end{equation} $ (5)

控制方程(3)变为：

$R{{e}^{p}}\frac{\partial {{\overline{u}}^{p}}(\overline{y},\overline{t})}{\partial \overline{t}}=\frac{\partial }{\partial \overline{y}}(\frac{\partial {{\overline{u}}^{p}}(\overline{y},\overline{t})}{\partial \overline{y}})-{{h}^{2}}\frac{\partial p}{\partial x}{{e}^{i\overline{t}}}$ (6)

假设流速分布满足初始条件：

$\overline{u}(\overline{y},0)=0\begin{matrix} {} & {} & {} \\ \end{matrix}\frac{\partial \overline{u}(\overline{y},0)}{\partial \overline{t}}=0$ (7)

引入拉普拉斯（Laplace）变换：

$U(\overline{y},s)=\int_{0}^{+\infty }{\overline{u}(\overline{y},\overline{t}){{e}^{-s\overline{t}}}}d\overline{t}$ (8)

结合式（7）和式（8）对式（6）进行Laplace变换后有：

$\frac{\partial }{\partial \overline{y}}[\frac{\partial {{U}^{p}}(\overline{y},s)}{\partial \overline{y}}]-R{{e}^{p}}s{{U}^{p}}(\overline{y},s)-\frac{{{h}^{2}}}{(s-i)}\frac{\partial p}{\partial x}=0$ (9)

求解方程(9),利用海维赛(Heaviside)公式^[22]对其结果进行Laplace逆变换,并还原无量纲化参数得压力梯度驱动下的流速分布为：

${{u}^{p}}(y,t)=\frac{{{P}_{0}}{{h}^{2}}}{\mu lR{{e}^{p}}i}[1-\frac{\cosh (\sqrt{R{{e}^{p}}i}\frac{y}{h})}{\cosh (\sqrt{R{{e}^{p}}i})}]\cdot {{e}^{i{{\omega }^{p}}t}}$ (10)

切应力τ_yx分布和流率Q可以表示为：$\tau _{yx}^{{}}(y,t)=-\mu \frac{\partial u(y,t)}{\partial y}\begin{matrix} {} & {} \\ \end{matrix}Q(t)=2\int_{0}^{w}{dz}\int_{0}^{h}{u(y,t)dy}$ (11)

根据式（10）和式（11）可以得到压力梯度驱动作用下切应力和流率的分布规律为：

$\tau _{yx}^{p}(y,t)=\frac{{{P}_{0}}h}{l\sqrt{R{{e}^{p}}i}}\frac{\sinh (\sqrt{R{{e}^{p}}i}\frac{y}{h})}{\cosh (\sqrt{R{{e}^{p}}i})}\cdot {{e}^{i{{\omega }^{p}}t}}$ (12)

${{Q}^{p}}(t)=\frac{2w{{P}_{0}}{{h}^{3}}}{\mu lR{{e}^{p}}i}[1-\frac{\tanh (\sqrt{R{{e}^{p}}i})}{\sqrt{R{{e}^{p}}i}}]\cdot {{e}^{i{{\omega }^{p}}t}}$ (13)

以上各式中Re^p=ρω^ph²/ μ称为压力梯度驱动流的振荡雷诺数,w和l分别为矩形微管通道的宽度和长度。

1.2 电场驱动流模型的建立

同压力梯度驱动流的假设相同,不考虑x方向的压力梯度的影响,电场驱动下的方程（1）可以表示为：

$\rho \frac{\partial {{u}^{e}}(y,t)}{\partial t}=\mu \frac{\partial }{\partial y}(\frac{\partial {{u}^{e}}(y,t)}{\partial y})+{{\rho }^{e}}(y){{E}_{0}}{{e}^{i{{\omega }^{e}}t}}$ (14)

其中,u^e(y, t)是电场驱动时x方向的流速分量（只考虑x方向的流动）,ρ^e是液流中净电荷分布密度,E₀和ω^e分别是沿着x方向电场的幅值和振荡角频率（ω^e=2πf ^e）。上述方程满足边界条件：

${{u}^{e}}(y,t){{|}_{y=h}}=0\begin{matrix} {} & {} \\ \end{matrix}\frac{\partial {{u}^{e}}(y,t)}{\partial y}{{|}_{y=0}}=0$ (15)

由静电学理论可知单位体积内的净电荷密度ρ^e的值取决于双电层（EDL）电势,对于低的EDL电势ψ,若电势能远低于热能,运用Debye-Hückel的线性化近似后可得Poisson-Boltzmann方程^[23]：

$\frac{{{d}^{2}}\psi }{d{{y}^{2}}}=-\frac{{{\rho }^{e}}(y)}{\varepsilon }={{\kappa }^{2}}\psi $ (16)

且满足边界条件：

$\psi {{|}_{y=h}}={{\psi }_{0}}\begin{matrix} {} & {} \\ \end{matrix}\frac{d\psi }{dy}{{|}_{y=0}}=0$ (17)

求解得净电荷密度分布^[24]：

${{\rho }^{e}}(y)=-\varepsilon {{\kappa }^{2}}{{\psi }_{0}}\frac{\cosh (\kappa y)}{\cosh (\kappa h)}$ (18)

其中,ε=ε₀ε_r为介电常数,ψ₀为腔壁的Zeta电势,κ=$\sqrt{\text{2}{{n}_{0}}{{z}^{\text{2}}}{{e}^{\text{2}}}/{{\varepsilon }_{0}}{{\varepsilon }_{r}}{{k}_{b}}T}$为Debye-Hückel参数,n₀=1000N_AC为离子浓度,C为电解质的摩尔浓度,N_A为阿伏伽德罗常数,z为离子所带电荷数,e为单位电荷电量,${{k}_{b}}$为波尔兹曼常数,$T$为热力学温度（T=t₀+273.15）。引入无量纲化参数：

$\begin{equation} \left. \begin{aligned} \overline{y}=\frac{y}{h} \\ K=\kappa h \\ \overline{t}=\frac{t}{1/{\omega^{e}}} \\ R{e}^{e}=\frac{\rho {\omega ^e}{h^2}}{\mu } \\ \overline{u}^{e}(\overline{y},\overline{t})=\frac{{u^e}(y,t)}{U^{eo}} \\ U^{eo}=-\frac{\varepsilon {{\psi }_0}{E_0}}{\mu} \end{aligned} \right\} \end{equation}$ $ (19)

控制方程（14）变为：

$R{{e}^{e}}\frac{\partial {{\overline{u}}^{e}}(\overline{y},\overline{t})}{\partial \overline{t}}=\frac{\partial }{\partial \overline{y}}(\frac{\partial {{\overline{u}}^{e}}(\overline{y},\overline{t})}{\partial \overline{y}})+{{K}^{2}}\frac{\cosh (K\overline{y})}{\cosh (K)}{{e}^{i\overline{t}}}$ (20)

假设流速分布满足初始条件（7）,结合式（8）对式（20）进行Laplace变换后有：

$\frac{\partial }{\partial \overline{y}}[\frac{\partial {{U}^{e}}(\overline{y},s)}{\partial \overline{y}}]-R{{e}^{e}}s{{U}^{e}}(\overline{y},s)+\frac{{{K}^{2}}}{(s-i)}\frac{\cosh (K\overline{y})}{\cosh (K)}=0$ (21)

求解方程(21),利用Heaviside公式对其结果进行Laplace逆变换,并还原无量纲化参数得电场驱动下的流速分布为：

${{u}^{e}}(y,t)=-\frac{\varepsilon {{\psi }_{0}}{{E}_{0}}{{\kappa }^{2}}{{h}^{2}}}{\mu (R{{e}^{e}}i-{{\kappa }^{2}}{{h}^{2}})}\{\frac{\cosh (\kappa y)}{\cosh (\kappa h)}-\frac{\cosh [\sqrt{R{{e}^{e}}i}\frac{y}{h}]}{\cosh [\sqrt{R{{e}^{e}}i}]}\}\cdot {{e}^{i{{\omega }^{e}}t}}$(22)

根据式（11）和式（22）,最终可以得到电场驱动作用下切应力和流率的分布规律为：

$\tau _{yx}^{e}(y,t)=\frac{\varepsilon {{\psi }_{0}}{{E}_{0}}{{\kappa }^{2}}{{h}^{2}}}{R{{e}^{e}}i-{{\kappa }^{2}}{{h}^{2}}}\{\kappa \frac{\sinh (\kappa y)}{\cosh (\kappa h)}-\frac{\sqrt{R{{e}^{e}}i}}{h}\frac{\sinh [\sqrt{R{{e}^{e}}i}\frac{y}{h}]}{\cosh [\sqrt{R{{e}^{e}}i}]}\}\cdot {{e}^{i{{\omega }^{e}}t}}$ (23)

${{Q}^{e}}(t)=-\frac{2\varepsilon {{\psi }_{0}}{{E}_{0}}w{{\kappa }^{2}}{{h}^{3}}}{\mu (R{{e}^{e}}i-{{\kappa }^{2}}{{h}^{2}})}\{\frac{\tanh (\kappa h)}{\kappa h}-\frac{\tanh [\sqrt{R{{e}^{e}}i}]}{\sqrt{R{{e}^{e}}i}}\}\cdot {{e}^{i{{\omega }^{e}}t}}$ (24)

以上各式中Re^e=ρω^eh²/ μ称为电场驱动流的振荡雷诺数,w为矩形微管通道的宽度。

1.3 力-电协同驱动流模型的建立

考虑力-电协同作用时,方程（1）同时保留压力梯度项和体力项（电场力）。这时在xy截面方程变为：

$\rho \frac{\partial {{u}^{pe}}(y,t)}{\partial t}=\mu \frac{\partial }{\partial y}(\frac{\partial {{u}^{pe}}(y,t)}{\partial y})-\frac{\partial p}{\partial x}{{e}^{i{{\omega }^{p}}t}}+{{\rho }^{e}}{{E}_{0}}{{e}^{i{{\omega }^{e}}t}}$ (25)

其中ω^p为压力的振荡角频率（ω^p=2πf ^p）,ω^e为电场的振荡角频率（ω^e=2πf ^e）。方程（25）满足边界条件：

${{u}^{pe}}(y,t){{|}_{y=h}}=0\begin{matrix} {} & {} \\ \end{matrix}\frac{\partial {{u}^{pe}}(y,t)}{\partial y}{{|}_{y=0}}=0$ (26)

引用式（5）和式（19）的无量纲化参数,将式（25）及其边界条件（26）无量纲化后进行Laplace变换并求解可得：

${{U}^{pe}}(\overline{y},s)=-\frac{{{h}^{2}}}{R{{e}^{p}}s(s-i)}\frac{\partial p}{\partial x}[1-\frac{\cosh (\sqrt{R{{e}^{p}}s}\overline{y})}{\cosh (\sqrt{R{{e}^{p}}s})}]+\frac{{{K}^{2}}}{(s-i)(R{{e}^{e}}s-{{K}^{2}})}[\frac{\cosh (K\overline{y})}{\cosh (K)}-\frac{\cosh (\sqrt{R{{e}^{e}}s}\overline{y})}{\cosh (\sqrt{R{{e}^{e}}s})}]$ (27)

将式（27）进行Laplace逆变换并还原无量纲参数可得：力-电协同作用下的流速解答为压力梯度、电场单独作用时解答的代数和,即为：

${{u}^{pe}}(y,t)={{u}^{p}}(y,t)+{{u}^{e}}(y,t)$ (28)

切应力$\tau _{yx}^{pe}$分布和流率Q^pe可以表示为：

$\tau _{yx}^{pe}(y,t)=-\mu \frac{\partial {{u}^{pe}}(y,t)}{\partial y}\begin{matrix} {} & {} \\ \end{matrix}{{Q}^{pe}}(t)=2\int_{0}^{w}{dz}\int_{0}^{h}{{{u}^{pe}}(y,t)dy}$ (29)

将式（28）代入式（29）可得,力-电协同作用下切应力和流率的解析解也为压力梯度、电场单独作用所得结果的代数和：

$\tau _{yx}^{pe}(y,t)=\tau _{yx}^{p}(y,t)+\tau _{yx}^{e}(y,t)\begin{matrix} {} & {} \\ \end{matrix}{{Q}^{pe}}(t)={{Q}^{p}}(t)+{{Q}^{e}}(t)$ (30)

2 数值参数

细胞培养液中的电解质以NaCl为主,其浓度设定为10^-6mol/L,微通道腔壁上对应的Zeta电势为-0.2V^[25],其它参数见表1。温度的变化会影响培养液的粘度和密度,其对应关系见表2。

3 结果

3.1 时程曲线

图2绘制了压力梯度驱动频率为2Hz和电场驱动频率分别为1Hz和2Hz的时程曲线,各驱动力时程曲线均为余弦曲线,波峰代表沿着x轴正向的最大值,波谷则代表反向的最大值。从图3中可以看出力-电协同驱动结果是压力、电场驱动单独作用结果的代数和。图3（a）（c）（e）分别绘制了压力、电场及力-电协同三种模式在不同步驱动（ f ^p=2 Hz, f ^e=1 Hz）时腔心（y=0）流速、腔壁（y=h,细胞贴壁处）上切应力及流率的时程曲线。当压力梯度（t=0.25s）和电场（t=0.5s）载荷到达波谷时,各自的响应（流速、切应力和流率）同时到达波谷;压力梯度（t=0.5s）和电场（t=1s）载荷到达波峰时,各自的响应也同时到达波峰。力-电协同的响应值在部分时段相互削弱,幅值达到最小值（t=0.5s）;部分时段叠加增强,幅值可达到最大值（t=1s）。图3（b）（d）（f）则绘制了三种模式在同步驱动（ f ^p=2 Hz, f ^e=2 Hz）时的时程曲线,压力梯度和电场载荷到达波谷（t=0.25s）和波峰（t=0.5s）时,各自的响应值也同时到达波谷和波峰。力-电协同作用时,流速、切应力及流率同步叠加增强,其时程曲线与压力梯度和电场同时到达波峰和波谷,且幅值达到最大值。据此,本文主要研究增强效果比较明显的压力和电场在同步协同驱动的情况。

显示原图| 下载原图ZIP| 生成PPT

图2   压力和电场随时间的变化,Re( )表示复数取实部

Fig. 2   Time responses of pressure and electric field, the operator Re( ) gives the real part of the complex number

显示原图| 下载原图ZIP| 生成PPT

图3   流速、切应力和流率随时间的变化

Fig. 3   Time responses of FV, FSS and FFR

3.2 流速幅值分布

图4为压力梯度驱动下流速幅值随着压力（a）、腔高（b）、频率（c）、温度（d）的变化图,由图可知：从腔壁（y=±h）到腔心（y=0）,流速从零非线性增大到最大值,在培养腔内呈“拱形”分布。随着压力、腔高和温度的增大,腔内各处的流速均随之增大（见图4（a）（b）（d））,但频率增大时流速却减小(见图4（c）)。对比各处的流速变化可知,腔心处的流速变化最大,腔壁处的流速变化最小。

显示原图| 下载原图ZIP| 生成PPT

图4   压力梯度驱动下流速幅值在不同压力（a）、腔高（b）、频率（c）、温度（d）的分布

Fig. 4   Velocity amplitude profile under pressure gradient driven flow at different pressure (a),chamber height(b), frequency(c), temperature(d)

图5为电场驱动下流速幅值随着电场（a）、腔高（b）、频率（c）、温度（d）的变化图,流速在培养腔内呈“梯形”分布。在电场的作用下,培养腔内的流速在一定区域内（约y=-0.998h~0.998h）保持恒定。当电场强度和温度增大时,该恒定流速随之线性增大（见图5（a）（d））。腔高和频率增大时,相同区域内流速则不能保持恒定,腔心附近的流速会随之减小,但腔壁附近（约为y=±0.998h处）的流速却能保持不变（见图5（b）（c））。此外,由图5（d）可知温度对流速的影响较小,所以电场驱动下腔壁附近的流速主要受控于电场强度。

显示原图| 下载原图ZIP| 生成PPT

图5   电场驱动下流速幅值在不同电场（a）、腔高（b）、频率（c）、温度（d）的分布

Fig. 5   Velocity amplitude profile under electric field driven flow at different electric field (a),chamber height(b), frequency(c), temperature(d)

对比图4和图5可以发现,压力仅在腔心部分驱动液体产生流速,而电场驱动却能在整个腔体维持一定的流速。压力幅值（P₀=400~520Pa）和电场幅值（E₀=(1.2~2.0)×10⁵V/m）一定时,压力梯度驱动得到的最大流速幅值较电场驱动的大（约为2~6倍）。对比图4(a)(d)和图5(a)(d),场强和温度的增大均能使得两场作用下的流速增大,压力梯度驱动下增大趋势在腔心处明显,电场驱动下增大趋势则分布整个腔体。对比图4(b)和图5(b),腔高变化对腔心处流速的影响非常明显,随着腔高增大,压力梯度驱动下的腔心流速逐渐增大,而电场驱动下的腔心流速逐渐变小。当腔高变化量相同时,压力梯度驱动下腔心流速幅值的变化范围要远大于电场驱动下的变化（约为200倍）。由此可见:腔高变化对电场驱动液体流动行为的影响要远小于压力梯度驱动下的。对比图4(c)和图5(c),两场的驱动频率增大时,流速变化均出现在腔心附近且呈减小趋势。

图6绘制了力-电协同（同步）驱动下的流速幅值随着场强（a）、腔高（b）、频率（c）、温度（d）的变化图,在力-电协同作用下,腔内各处均能保持一定的流速。当压力/电场强度和温度的增大时,腔内各处的流速均随之增大,腔心处的流速变化要大于腔壁处的流速变化（见图6(a)(d)）。但腔高减小或频率增大时,腔心附近的流速减小,腔壁附近的流速则保持不变（见图6(b)(c)）。力-电协同作用时,压力（P₀=460Pa）和电场（E₀=2×10⁵V/m）一定时,腔壁附近（约为y=±0.95h处）的流速主要由电场驱动（约占80%）提供,腔心附近的流速主要由压力梯度驱动（约占70%）提供,可见电场的作用弥补了压力梯度驱动在腔壁附近流速小的劣势。对比图4、图5和图6可知：力-电协同作用时,腔心附近的流速可以通过压力、腔高、频率、温度四个参量进行调控,但腔壁附近流速通过电场调控较为明显。

显示原图| 下载原图ZIP| 生成PPT

图6   力-电协同驱动下流速幅值在不同场强（a）、腔高（b）、频率（c）、温度（d）的分布

Fig. 6   Velocity amplitude profile under pressure-electricity synergy driven flow at different field intensity(a),chamber height(b), frequency(c), temperature(d)

3.3 腔壁上的切应力幅值

图7绘制了压力、电场和力-电（f ^p =f ^e=1Hz）协同驱动三种驱动模式下,腔壁上（y=±h,细胞贴壁处）流体切应力幅值随压力（a）、电场（b）、腔高（c）、频率（d）、温度（e）的的变化规律。由图7(a)和(b)可知切应力幅值随着压力和电场的场强幅值增大呈线性增长关系,当切应力幅值增大量相同时,电场的变化值要远大于压力的变化值（约为50倍）,由此可见电场驱动时电场强度具有更大的调节空间。对比图7(a)和图7(b),压力（P₀=400~520Pa）驱动下的切应力仅为电场（E₀=(0.6~1)×10⁴V/m）驱动下切应力的1/4。当压力（P₀=460Pa）一定时,切应力幅值随着腔高增大线性增大(见图7(c)),随着频率的增大非线性减小(见图7(d)),随着温度的变化基本保持不变(见图7(e))。在电场（E₀=10⁴V/m）驱动下,腔高、频率和温度的改变时,切应力基本保持不变。力-电协同作用下切应力随各量的变化趋势和压力梯度驱动下的变化趋势相近,其切应力结果为力、电单独作用下切应力的代数和。对比三种驱动模式下切应力的变化可知：压力梯度驱动下,影响切应力变化的外因较多,可控性差;电场驱动下,切应力仅受控于电场强度,受其它条件约束少,可控性强。

显示原图| 下载原图ZIP| 生成PPT

图7   腔壁上切应力幅值随压力（a）、电场（b）、腔高（c）、频率（d）、温度（e）的变化

Fig. 7   Evolutions of FSSA with pressure (a), electrical field (b),chamber height(c), frequency(d), temperature(e) on the chamber wall

3.4 流率幅值

图8是三种驱动模式下流率幅值随压力（a）、电场（b）、腔高（c）、频率（d）、温度（e）的变化曲线。由图8(a) (b) (e)可知：当压力、电场和温度升高时,三种驱动模式下的流率都随之线性增大。当腔高增大时,压力梯度驱动及力-电协同驱动下的流率随之非线性增大,电场驱动下的流率保持不变（见图8(c)）。频率增大时,三种模式下的流率均随之非线性减小（见图8(d)）。与图7对比,在施加的压力和电场幅值一定（P₀=460Pa,E₀=10⁴V/m）的情况下,电场驱动得到的切应力比压力梯度驱动的大（约为4~9倍）,但电场驱动下的流率却远小于压力梯度驱动（约为1/30）。

显示原图| 下载原图ZIP| 生成PPT

图8   流率随压力（a）、电场（b）、腔高（c）、频率（d）、温度（e）的变化

Fig. 8   Evolutions of FFRA with pressure (a), electrical field (b),chamber height(c), frequency(d), temperature(e)

4 讨论

本文建立了基于微流控的矩形（腔宽远大于腔高）细胞培养流动腔的理论模型,该模型可以用来进一步量化其内部液体的压力、流速、腔壁（细胞贴壁处）的切应力和流率分布,并得到了压力梯度、电场、力-电协同驱动微液流动的解析解。根据Stavenschi等^[6]的实验参数将腔高设置为300μm、进出口压力差的幅值和频率分别为460Pa和1Hz时,模型计算出的流体切应力和流率分别为 0.92Pa 及 45.24 ml/min,理论计算值基本与实验结果吻合（如表3所示）。本文详细考察了外部物理场载荷（压力和电场驱动幅值、频率）和环境（腔高、温度）参数对细胞培养的微流动刺激的影响,同时也对比了各自物理场驱动的效果,为外加物理场刺激组织（细胞）生长的机理研究提供理论参考。

考虑到人体组织承受的生理载荷环境,本文采用周期性载荷的驱动方式（振荡流）,这样更接近细胞的体内生理环境^[28]。相关研究表明^[29]振荡流可以促进细胞在培养腔中更加均匀的分布和增殖。压力和电场单独驱动液体流动时,由时程曲线（图3）可知施加的驱动载荷和流体刺激结果（腔心流速、腔壁切应力、流率）均能同步达到峰值,这说明当培养液的粘性较小时,细胞感受力学信号的即时性非常好。当压力和电场协同作用时,结果呈线性叠加关系,可见压力和电场之间不存在相互干扰。在微流控培养过程中,细胞的正常生长/分化需要定量的流速、切应力和流率,三者缺一不可,但单独的压力和电场驱动却不能保证三者均在引起细胞响应的流体刺激范围（流速为400~800 μm/s^[30]、切应力为0.2~6 Pa^[31]）之内（见图4、5、7）。由图3~图8可知：两驱动场相比,电场驱动善于提供较大的腔壁流速和切应力,压力梯度驱动善于提供较大的流率。压力梯度驱动对细胞贴壁处流体的控制能力较弱,而电场恰好能弥补这一不足。因此,力-电协同驱动的效果（流速、切应力和流率均具有可调控性）将更加明显。从图3还可看出,当压力和电场驱动同步时（频率相同）,流体刺激结果同步线性叠加且相互不干扰,力-电协同作用下液流刺激效果增强最为明显。

流速的增大对细胞均匀分布生长有一定的促进作用,这可能是由于流速的改变会影响营养物质的运输速率,进而改变细胞的生长取向。速度梯度（y轴方向）的增大也会引起较强的剪切刺激,研究发现^[30]培养液流速会影响骨细胞形态和胞外基质的分泌和沉积,灌流速度为 400~800 μm/s 时成骨响应最佳。由于培养过程中细胞大多贴壁生长,因此控制腔壁附近的流速至关重要。在腔壁附近,压力梯度驱动流下的流速基本接近于零（见图4）,而电场驱动流却能保持相对稳定的流速（见图5）,这与song等^[32]所得流速图的结果大致吻合。由图6可知力-电协同作用时腔壁附近的流速大小主要由电场来调控：为了在腔壁附近获得相对稳定的流速,可通过改变电场强度来实现（见图5（a））。当电场强度控制为2×10⁵V/m时,距离腔壁约h/500处的流速大小保持为0.038m/s,由于流速正比于电场场强（与Bera等^[33]所得的流速与电场场强成线性增长关系的结果一致）,为了满足成骨的最佳响应值^[30]（400~800 μm/s）,电场强度应当减小50~100倍。如果流速代表物质的运输速率,流率则代表着单位时间内的物质运输量,保证腔壁处的流体刺激（流速及切应力（0.2~6 Pa^[31]））在细胞响应范围内的同时,为了满足细胞的新陈代谢——养料、废物的及时供应与排除,也需要保证一定的流率^[34]。流率因驱动方式的不同而不同,从图8中可以看出：在一定的驱动幅值内,电场驱动(E₀=10⁴V/m)下的流率远小于压力梯度驱动(P₀=460Pa)下的,当力-电协同作用时可获得较大的流率。压力梯度驱动下流率的波动范围受驱动载荷和环境参数的影响较大,可以根据实际需要选择力-电协同驱动模式弥补不足。

流体剪切应力是能刺激细胞响应的重要力学因素信号,在细胞微流控培养的生物反应器设计中是被着重考虑的流体刺激因子。一些研究^[35,36]根据流率估算切应力,且认为距离腔壁小于h/4的范围内均受切应力的影响,本文则根据牛顿内摩擦定律精确获得了腔壁处流体切应力的大小。腔壁处切应力随着压力和电场的增大均呈线性增大趋势,与压力梯度驱动相比腔高的变化对电场驱动下的切应力基本没有影响（图7）。根据腔壁处流速的取值（400~800 μm/s^[30]）可以确定电场的变化范围应为2000~4000V/m,相应得到的切应力范围约为1~2Pa。基因的表达决定着细胞的分化方向,总结Stavenschi等^[6]实验结果发现切应力为2Pa时对成骨的相关基因表达最有利,可见调节外部电场强度就能从基因水平上调控细胞的分化。由于电场驱动下的切应力受腔高影响较小,这样在使用电场驱动时勿需对腔高作过多要求,但可以适当调控腔高以保证合适的流率。压力梯度驱动下的切应力受腔高的影响较大（正比关系）,且压力梯度驱动幅值控制为460Pa以下时切应力幅值在1Pa以下,减小腔高甚至可以接近很低的切应力水平（见图7(b)）。Gao等^[5]发现极低切应力（10^-3Pa）能促使MSCs的迁移和成骨分化,这种极低切应力是非常值得关注的。例如在人静坐时,对骨组织的力刺激仅来源于肌肉的收缩,相比运动状态,这种微小的外力极易产生较低的切应力刺激。Zheng等^[7]发现MSCs在切应力的刺激下会产生收缩和再扩散现象,如实验中对细胞施加不同的切应力（1.3Pa、9.8Pa和14.7Pa）时会影响细胞开始收缩的时间,并发现切应力越大细胞从接种到开始收缩的时间越短。考虑到高强度电场取值的困难,压力梯度驱动所得的腔壁切应力较小等情况,可以采用力-电协同驱动的方式来达到细胞组织工程所需要的切应力刺激水平。

在研究载荷驱动频率对细胞生长的影响时,Tanaka等^[37]发现外部物理场驱动频率为2Hz时最有利于成骨细胞的矿化。结合图5（c）和图7（d）,频率的改变对腔壁处流速和切应力的影响均较小,因此频率对细胞的响应通路不是通过改变流速和切应力的大小来实现的。压力和电场驱动频率为2Hz时均能得到较大的流率,因而此频率可以为灌流培养过程提供设计参考。由此可见：细胞对载荷幅值的敏感性要比对载荷频率的敏感性强。本文选取的温度值在37℃左右小幅波动,对切应力的影响可以忽略（如图7（e））。Stolzing^[38]研究也发现在37℃的基础上适当降低温度不影响MSCs的活性,但分化程度增加、氧化损伤的积累降低、细胞凋亡减少。温度降低时,培养液粘度会增大,流速、流率随之减小,但减小的幅度不大（见图4（d）、图5（d）、图8（e））,所以在实际的灌流培养中可根据细胞的具体需要适当调控培养温度。另外,在微流控培养腔中,利用电场进行驱动时,通常是通过对浸入培养液中的电极施加电压信号实现的。而在这个过程中当电压较高时(>2.2V),电极附近就不可避免的会发生电化学反应,进而会导致附近培养液pH值的相应变化^[39,40]。对于非饱和的强电解质（NaCl）水溶液而言,电解主要以水解反应为主,当正负电极附近产生H⁺和OH^-离子后,在电场的作用下,带电离子向相反的电极运动,这些带电离子在各自相对的电极被完全中和。由此可见：电极的电解作用仅对电极附近的pH值有一定影响,对培养液的整体无影响。这与毕等^[41]的研究结果是相同的,他们设定的溶液初始pH为7.8,电场作用10min,溶液的导电性稳定后,阳极处的pH值变为7.0,阴极处的pH值变为8.5。杨等^[42]通过实验发现MSCs在pH值为7.2-7.4的培养液中生长活跃,pH值大于7.4时,逐渐停止增殖,pH值小于7.2时,增殖速度减慢。Lee等^[43]发现成骨细胞在pH=7.4的环境下生理现象一切正常,酸性条件下成骨细胞仍能够进行增殖,但随着培养基pH值的降低,成骨细胞的活性和增殖受到抑制,细胞凋亡增加。可见,细胞对碱性环境较为敏感,但对酸性却有一定的耐受性,且适应的酸性pH值有一定的可变空间。在进行细胞培养时,可将细胞接种于远离负电极的位置。由于电解速率不变,也可适当增大腔高使培养腔的体积增大,减缓培养液pH值的变化^[44]。

本文建立了力-电协同驱动下的细胞微流控培养腔模型,并分析了细胞灌流培养在压力梯度和电场驱动作用下流体刺激信号的传导规律,对比了两物理场单独和协同作用时的特点。许多实验研究者都尝试施加两个物理驱动场来探索最佳的流体刺激条件,但是通过实验来探索合适的刺激参数周期长、代价高且效果不明显。本文在外加物理场和细胞的响应流体切应力之间建立一个数学关系,根据细胞的切应力和流速的响应范围大致确定外部物理场的取值,这样可以为实验参数的初调提供参考。本模型的建立将细胞能感受到的流体刺激信号和外在的物理环境有机的联系在一起,是一个完整的信号传导系统,在特定驱动方式下,微流控细胞培养相关的研究内容能在这一模型下开展。模型仍存在着以下不足：⑴当忽略细胞高度时,细胞贴壁附近的流场会发生较大变化,对切应力结果有一定影响;⑵固液边界是理想无滑移的;⑶培养液假设为牛顿流体,非牛顿流体更接近真实体液的属性;⑷忽略了电粘性效应,培养液内的离子在压力梯度驱动流的作用下产生流动电流（场）,该流动电流（场）对外加电场及流体均有影响。

5 结论

人体通常处于复杂的受力（行走、静坐、跳跃）环境中,当这些外载荷力学信号传递到组织的微观结构尺度时以微液流动刺激细胞响应的方式来引起组织结构重建或者塑建,最终使得组织以最优结构适应外部力学环境。这是组织环境中的载荷信号传递及反馈机理,但这一机理的科学量化通路还不是很清楚,特别是细胞参与下的力学生物学机制还不是很明确。本文建立的力-电协同驱动矩形截面流的细胞培养腔理论模型,将影响细胞生长/分化的液体压力、流速、流率、切应力和外部物理环境（压力梯度、电场与温度等）联系在一起,可进一步为力-电耦合刺激细胞生长/分化机制研究及细胞组织工程参数优化提供参考。尽管模型还存在着不足,但初步可以得出下列参考结论：

1）力-电协同作用的流体刺激效果（解答）是单独压力梯度和电场驱动作用效果（解答）的代数叠加。

2）细胞培养腔内流体流速、切应力和流率均跟外部物理场（压力、电场）的场强呈线性正比关系。物理场的驱动频率增大时,压力梯度驱动下的腔壁切应力非线性递减而电场驱动下的腔壁流速和切应力则基本不变,但两场驱动下的流率均非线性减小。可见,细胞所感受的流体力学刺激信号（腔壁处的流速和切应力）对载荷幅值的变化相比频率更为敏感。

3）当培养腔的腔高增大时,压力梯度驱动下的腔壁切应力和流率分别呈线性和非线性增大,电场驱动下的腔壁流速、切应力及流率则均无变化。

4）培养液的温度（35℃~39℃）对压力、电场和力-电协同三种驱动模式下微液流动行为的影响均不明显。

5）在细胞响应范围内,电场驱动善于提供较大的流体切应力幅值而压力梯度驱动善于提供较大的流率幅值。力-电协同作用时,可以联合压力梯度和电场驱动的调控优势。

The authors have declared that no competing interests exist.

参考文献

文献选项

• 原文顺序
• 文献年度倒序
• 文中引用次数倒序
• 被引期刊影响因子